帮忙做下统计显著性检验和K值的误差以及灵敏度分析

箫剑→残念

function parafit
%  k1->k-1，k2->k1，k3->k2,k4->k3,k5->k4
% k6->k6 k7->k7
% dGlcdt = k-1*C(Fru)-(k1+k2)*C(Glc);
. E% V% G5 @$ y  z2 A% dFrudt = k1*C(Glc)-(k-1+k3+k4)C(Fru);
% dFadt = k(2)*C(Glc)+k4*C(Fru)+（k6+k7)*C(Hmf);
% dLadt = k(7)*C(Hmf);
) F- p: I* J- c$ r% s%dHmfdt = k(3)*C(Fru)-（k6+k7)*C(Hmf);
$ Y& Z3 i$ T& U& ]" L; g* gclear all
clc
& L! `+ J8 ^2 p! z) Zformat long
%        t/min   Glc    Fru        Fa   La   HMF/ mol/L
) K& l1 @. J, {% ?5 X# {  Kinetics=[0    0.25    0           0    0       0
& p$ q3 |  K, r% A; Z# r          15    0.2319    0.01257    0.0048    0    2.50E-04
          30    0.19345    0.027    0.00868    0    7.00E-04
. D' ^0 e3 x4 P6 q, L          45    0.15105    0.06975    0.02473    0    0.0033
, s' v5 `* r3 }3 K. l8 P          60    0.13763    0.07397    0.02615    0    0.00428
( X. o0 O: ^1 v( k6 ~2 ?0 T8 L- n          90    0.08115    0.07877    0.07485    0    0.01405
3 J4 s; U4 R- {' j  _8 F          120    0.0656    0.07397    0.07885    0.00573    0.02143
+ A. v  q6 U" i8 J: ~' f; q          180    0.04488    0.0682    0.07135    0.0091    0.03623
  h7 ~  {/ h# M+ t( d" C          240    0.03653    0.06488    0.08945    0.01828    0.05452
" M3 A! `( y6 ?; d, i4 z/ f          300    0.02738    0.05448    0.09098    0.0227    0.0597
6 T$ v0 \( }& n5 B0 W- Y6 H$ i          360    0.01855    0.04125    0.09363    0.0239    0.06495];
6 `: a8 `( H# W) l5 z3 Ck0 = [0.0000000005  0.0000000005  0.0000000005  0.00000000005  0.00005  0.0134  0.00564  0.00001  0.00001  0.00001];        % 参数初值
lb = [0  0  0  0  0  0  0  0  0  0];                  % 参数下限
  d3 z1 t5 H: t8 R2 ]8 ~) cub = [1  1  1  1  1  1  1  1  1  1];    % 参数上限
3 M" {# p9 s2 f, \3 d8 {* rx0 = [0.25  0  0  0  0];
yexp = Kinetics;                 % yexp: 实验数据[x1        x4        x5        x6]
% warning off
# k3 {0 B- x2 ?, f% 使用函数 ()进行参数估计
[k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
fprintf('\n使用函数fmincon()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
( _9 D. B4 U! l% l2 Wfprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
' X) V$ M* ^) Y) u# Q) lfprintf('\tk6 = %.11f\n',k(6))
5 `* _9 l: ]6 W% Vfprintf('\tk7 = %.11f\n',k(7))
, ]7 t( d) U0 w* sfprintf('\tk8 = %.11f\n',k(8))
fprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
fprintf('  The sum of the squares is: %.1e\n\n',fval)
k_fm= k;
4 e) n# B0 I+ C2 C2 J! c$ h7 T  q+ H% Z% warning off
8 Y5 f4 R, Z! q$ T, d% 使用函数lsqnonlin()进行参数估计
7 U# C9 D& u' h[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
+ G2 D6 ]) {6 |3 m+ \( V    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
5 j# A- {$ ^6 m1 L; C+ pfprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
( {) Q4 E2 B; Z7 Nfprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
2 N* T8 h* }' {1 v3 d7 r# P2 ?, Sfprintf('\tk8 = %.11f\n',k(8))
  `' T) _) M* y3 e% C$ Ofprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_ls = k;
6 i" y% n1 ^' ^& D7 F# boutput
warning off
9 w# r3 Y5 A' O( m! H4 C% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
& A, k8 u- F* lk0 = k_fm;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
ci = nlparci(k,residual,jacobian);
* M( U" V% Y% b) S! @fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
) G+ I0 h( Z, U6 f; |( hfprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
0 k0 H$ U6 a2 ?+ `4 w) [fprintf('\tk6 = %.11f\n',k(6))
! L: ]4 z: W0 e5 Y  y" l( \fprintf('\tk7 = %.11f\n',k(7))
4 m3 e% w; D2 X1 S; N+ qfprintf('\tk8 = %.11f\n',k(8))
4 j( u' c1 ~$ M0 ^* Z3 Xfprintf('\tk9 = %.11f\n',k(9))
4 i- C7 p4 |  Vfprintf('\tk10 = %.11f\n',k(10))
5 x3 p1 z6 T! X- F! K6 ffprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_fmls = k;
8 L) V' D7 T: X* q; M  [: Woutput
tspan = [0 15 30 45 60 90 120 180 240 300 360];
9 _- o% `9 ^$ @: G2 n" L[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
plot(t,x(:,1),t,yexp(:,2),'*');legend('Glc-pr','Glc-real')
! ?. Q& @$ }8 Ifigure;plot(t,x(:,2:5));
p=x(:,1:5)
! y, Q1 J2 P- Q2 h) d+ T/ Shold on
plot(t,yexp(:,3:6),'o');legend('Fru-pr','Fa-pr','La-pr','HMF-pr','Fru-real','Fa-real','La-real','HMF-real')
. U' R4 m' D1 G3 x, [3 q
6 T: U# Z2 V% E2 e- i2 P. k+ l9 g

function f = ObjFunc7LNL(k,x0,yexp)
tspan = [0 15 30 45 60 90 120 180 240 300 360];
2 s0 P) Q% S2 f, s[t, x] = ode45(@KineticEqs,tspan,x0,[],k);   
' r; C. G; D. k- a) T5 D9 V4 vy(:,2) = x(:,1);
* b5 _5 U) E6 d: q& K& ^% ?. P& fy(:,3:6) = x(:,2:5);
f1 = y(:,2) - yexp(:,2);
f2 = y(:,3) - yexp(:,3);
f3 = y(:,4) - yexp(:,4);
. o5 C# A  p+ Jf4 = y(:,5) - yexp(:,5);
0 H7 H4 q7 N- n. Q* U! I% Bf5 = y(:,6) - yexp(:,6);
f = [f1; f2; f3; f4; f5];



function f = ObjFunc7Fmincon(k,x0,yexp)
tspan = [0 15 30 45 60 90 120 180 240 300 360];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);   
y(:,2) = x(:,1);
7 P" b: ^+ W& U- a5 o+ W: j6 Uy(:,3:6) = x(:,2:5);
f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2)   ...
- t" c- G, X+ m7 n/ A3 T    + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2)   ...
    + sum((y(:,6)-yexp(:,6)).^2) ;
1 m2 ^: |1 Q1 w4 ]
  b1 z2 b  e) E: M, u. R4 q
( a2 _" n) Q3 ^7 ?4 V$ C

function dxdt = KineticEqs(t,x,k)
* {  c$ I! z# ?; A5 S) i: @dGldt = k(1)*x(2)-(k(2)+k(3)+k(8))*x(1);
dFrdt = k(2)*x(1)-(k(1)+k(4)+k(5)+k(9))*x(2);
dFadt = k(3)*x(1)+k(5)*x(2)+(k(6)+k(7))*x(5);
dLadt = k(7)*x(5);
* i3 m" ?/ ^& v8 adHmdt = k(4)*x(2)-(k(6)+k(7)+k(10))*x(5);
dxdt = [dGldt; dFrdt; dFadt; dLadt; dHmdt];

: U6 \3 A* V5 _1 ]0 j  S% a
+ A# ^3 J4 @+ A

董事长之友

蒙的一比 大哥