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

箫剑→残念

function parafit
' h8 A3 K$ }* M. v: i+ g4 |- G%  k1->k-1，k2->k1，k3->k2,k4->k3,k5->k4
7 u7 x" f4 m2 F2 D8 f% k6->k6 k7->k7
; ?$ n; M$ P: c4 ^/ _+ t% G4 Y% dGlcdt = k-1*C(Fru)-(k1+k2)*C(Glc);
% dFrudt = k1*C(Glc)-(k-1+k3+k4)C(Fru);
7 r* I- X% c; Q% dFadt = k(2)*C(Glc)+k4*C(Fru)+（k6+k7)*C(Hmf);
- Q, p3 _+ K) C+ p) d7 D* \" [& b% dLadt = k(7)*C(Hmf);
%dHmfdt = k(3)*C(Fru)-（k6+k7)*C(Hmf);
& ^, m3 \6 g0 ^: Qclear all
1 O3 h: c& s; h9 G1 Z; B  F: yclc
format long
%        t/min   Glc    Fru        Fa   La   HMF/ mol/L
0 O# L" R/ U, x  s( F( n% F) S1 a, _  Kinetics=[0    0.25    0           0    0       0
          15    0.2319    0.01257    0.0048    0    2.50E-04
          30    0.19345    0.027    0.00868    0    7.00E-04
6 p2 J, ~1 f" ^( \          45    0.15105    0.06975    0.02473    0    0.0033
          60    0.13763    0.07397    0.02615    0    0.00428
$ b, _: a+ L0 U1 j          90    0.08115    0.07877    0.07485    0    0.01405
          120    0.0656    0.07397    0.07885    0.00573    0.02143
2 n% z. [+ j- q" R" O          180    0.04488    0.0682    0.07135    0.0091    0.03623
          240    0.03653    0.06488    0.08945    0.01828    0.05452
4 p2 P+ ?  m) x2 D4 s% o9 V, Q          300    0.02738    0.05448    0.09098    0.0227    0.0597
          360    0.01855    0.04125    0.09363    0.0239    0.06495];
k0 = [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];                  % 参数下限
$ a: O* Z7 T" S: T& z8 Y& c* sub = [1  1  1  1  1  1  1  1  1  1];    % 参数上限
x0 = [0.25  0  0  0  0];
* L! \& \( _$ G) z  u  }# Q2 ?yexp = Kinetics;                 % yexp: 实验数据[x1        x4        x5        x6]
% warning off
0 K* [" ]* H5 _$ s, O% 使用函数 ()进行参数估计
' h7 U" U1 S# @7 t, e# d% w/ ~[k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
2 ?. b& z- U- r7 xfprintf('\n使用函数fmincon()估计得到的参数值为:\n')
7 y9 J7 w: x2 V6 e& yfprintf('\tk1 = %.11f\n',k(1))
6 h1 u* ~/ X: [6 @+ M4 S3 lfprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
% t$ q  w! q6 A' y, ^; O7 Afprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
fprintf('\tk9 = %.11f\n',k(9))
7 d4 [: G' [1 t+ F  q5 S" X, K$ [fprintf('\tk10 = %.11f\n',k(10))
5 h. w: ]+ B7 B$ sfprintf('  The sum of the squares is: %.1e\n\n',fval)
k_fm= k;
% warning off
% 使用函数lsqnonlin()进行参数估计
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
0 g7 @# Z, C9 R6 v' B    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
. ^+ z7 ]; P+ _" o% [2 Cci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
6 x  G; h4 v! o) z5 d* s$ Q; |2 pfprintf('\tk1 = %.11f\n',k(1))
$ t3 g3 z  V1 pfprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
; J& l3 f) B! K1 H. G& w: Jfprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
4 p$ j% q1 e& [- h+ Efprintf('\tk6 = %.11f\n',k(6))
$ y7 j( [- b1 E  ?! ~fprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
2 k$ v8 J, U  N/ [" f- T# u1 Z8 U: sfprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
4 f" i& G+ n2 ~fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
9 h% }3 s7 T9 L2 ^8 r2 Tk_ls = k;
output
' }. c$ ~6 n. G1 f6 @warning off
  I/ g3 i4 k; r9 ~" \  e- ~0 U: p% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
8 E3 r! N4 ]* z* c' ~6 Ok0 = k_fm;
/ y, t1 F5 U4 h8 C3 N6 v' k[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
+ M. w6 L+ n& _' h% u' z    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
ci = nlparci(k,residual,jacobian);
fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
# Z3 ^8 _2 C( z. h9 Z- Rfprintf('\tk1 = %.11f\n',k(1))
& [; E) N; a7 h0 e3 x# Efprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
1 \) \5 }1 W, ?0 _( I2 kfprintf('\tk4 = %.11f\n',k(4))
  p& q, k: G- yfprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
' w" `5 ~) U2 K. v2 U- @* Z  hfprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
2 w& N. }$ Z; k. z4 y$ hfprintf('\tk9 = %.11f\n',k(9))
& G" e3 p- ?9 J- Ifprintf('\tk10 = %.11f\n',k(10))
: m0 K1 u$ e/ e5 `+ G+ G. K. {fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_fmls = k;
' H9 J( Z& g8 U* }! }. k( Voutput
0 q8 n3 Q  e* c' i1 l9 s# r8 V# ktspan = [0 15 30 45 60 90 120 180 240 300 360];
[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
7 B8 p$ N2 s; f4 t1 {3 splot(t,x(:,1),t,yexp(:,2),'*');legend('Glc-pr','Glc-real')
figure;plot(t,x(:,2:5));
9 h0 ?4 t! k1 N3 k" c+ Tp=x(:,1:5)
hold on
  ~7 c  O) Y$ O* Y/ }. e+ Lplot(t,yexp(:,3:6),'o');legend('Fru-pr','Fa-pr','La-pr','HMF-pr','Fru-real','Fa-real','La-real','HMF-real')



6 Y' c7 E+ U6 qfunction f = ObjFunc7LNL(k,x0,yexp)
tspan = [0 15 30 45 60 90 120 180 240 300 360];
5 x& E8 u% g  _' i( D0 l[t, x] = ode45(@KineticEqs,tspan,x0,[],k);   
y(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
4 N6 V: N# m5 x1 H3 w! @f1 = y(:,2) - yexp(:,2);
2 i7 ?) H8 w4 F0 X8 P. [" Kf2 = y(:,3) - yexp(:,3);
f3 = y(:,4) - yexp(:,4);
+ F+ z) G- h* A: B0 L! I, Af4 = y(:,5) - yexp(:,5);
: [) w6 |6 \) j$ L; Q! I. ]( sf5 = y(:,6) - yexp(:,6);
- ~6 }" h( P3 u# u* y: ~f = [f1; f2; f3; f4; f5];

( a% }+ }( `+ X1 `

function f = ObjFunc7Fmincon(k,x0,yexp)
  I  J8 b* I& S& ^% etspan = [0 15 30 45 60 90 120 180 240 300 360];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);   
  ]. i( a. y) @3 Oy(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2)   ...
# r: d; P( ~4 o" s6 D: d- Y    + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2)   ...
8 M4 \" q6 |  p. Z    + sum((y(:,6)-yexp(:,6)).^2) ;

3 M# K0 x  w% v
+ v3 s6 d- L$ j7 P& L" L& C6 m

function dxdt = KineticEqs(t,x,k)
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);
0 k; i( S/ s3 {dFadt = k(3)*x(1)+k(5)*x(2)+(k(6)+k(7))*x(5);
" Y/ l. t4 f6 {* z6 ~3 w- z/ MdLadt = k(7)*x(5);
dHmdt = k(4)*x(2)-(k(6)+k(7)+k(10))*x(5);
dxdt = [dGldt; dFrdt; dFadt; dLadt; dHmdt];
4 L+ Y) {1 ^- i) d5 A! }