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

箫剑→残念

function parafit
%  k1->k-1，k2->k1，k3->k2,k4->k3,k5->k4
7 e. h2 d8 b6 }9 R8 s- F' s% k6->k6 k7->k7
% dGlcdt = k-1*C(Fru)-(k1+k2)*C(Glc);
% dFrudt = k1*C(Glc)-(k-1+k3+k4)C(Fru);
% dFadt = k(2)*C(Glc)+k4*C(Fru)+（k6+k7)*C(Hmf);
2 a* o5 G2 O4 T4 t0 [% dLadt = k(7)*C(Hmf);
%dHmfdt = k(3)*C(Fru)-（k6+k7)*C(Hmf);
4 v5 t. {8 P: q8 x4 F" ]4 Vclear all
5 \! F! [* m1 ^3 yclc
format long
0 w: y/ q! v3 g  G. n%        t/min   Glc    Fru        Fa   La   HMF/ mol/L
  Kinetics=[0    0.25    0           0    0       0
/ A% k: c$ G  C& T; R          15    0.2319    0.01257    0.0048    0    2.50E-04
6 j( D; ~6 D4 r) {9 x          30    0.19345    0.027    0.00868    0    7.00E-04
8 S# ]+ }( }! W3 G+ C; I          45    0.15105    0.06975    0.02473    0    0.0033
6 a$ \1 `- C( E4 a# p- Y          60    0.13763    0.07397    0.02615    0    0.00428
          90    0.08115    0.07877    0.07485    0    0.01405
          120    0.0656    0.07397    0.07885    0.00573    0.02143
          180    0.04488    0.0682    0.07135    0.0091    0.03623
8 h5 J: d. h# }' G          240    0.03653    0.06488    0.08945    0.01828    0.05452
0 p5 a- b5 y& \          300    0.02738    0.05448    0.09098    0.0227    0.0597
          360    0.01855    0.04125    0.09363    0.0239    0.06495];
) z, N* F+ p5 _- u- l% `0 K& Tk0 = [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];                  % 参数下限
# X' Z( U" j& S5 A$ `ub = [1  1  1  1  1  1  1  1  1  1];    % 参数上限
x0 = [0.25  0  0  0  0];
' }3 ?" E5 O0 v9 \: K9 ^yexp = Kinetics;                 % yexp: 实验数据[x1        x4        x5        x6]
5 T* P  [7 C% a9 n/ R% warning off
% 使用函数 ()进行参数估计
[k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
fprintf('\n使用函数fmincon()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
' A+ A; `( b+ e" S. [5 F$ i# [fprintf('\tk2 = %.11f\n',k(2))
9 n1 D2 Z3 A0 r# G1 w, Xfprintf('\tk3 = %.11f\n',k(3))
# j) D' ~' N5 S9 x- v( H3 Ofprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
! o' c  O7 `; T& T% J' S3 m9 gfprintf('\tk6 = %.11f\n',k(6))
- V: y$ V  X% D. `& lfprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
fprintf('\tk9 = %.11f\n',k(9))
6 w7 z, I( C1 L% x/ |" ~, g* e% p8 Pfprintf('\tk10 = %.11f\n',k(10))
fprintf('  The sum of the squares is: %.1e\n\n',fval)
) H4 d7 V2 |% wk_fm= k;
% warning off
% 使用函数lsqnonlin()进行参数估计
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
ci = nlparci(k,residual,jacobian);
! L8 J% @* x4 L: ^fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
$ z. A$ m& B% a. e; s- K  o5 dfprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
6 U7 D: q3 q0 t! p$ f0 A) ]0 ifprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
, O% C5 A9 a# [+ U4 c! kfprintf('\tk8 = %.11f\n',k(8))
/ S, a& W. n% @7 _) m8 A; v& cfprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
1 m+ `% i$ N2 z7 i0 a8 \fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
& ~8 ^$ G6 n' ]1 E" _8 z4 Wk_ls = k;
' N0 K# V) D3 [3 `7 m; xoutput
, k; N( g; c3 {6 I# w" O; rwarning off
  l) d' B- A4 r2 R% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
6 R1 z# o' ^0 b$ L0 x7 @0 rk0 = k_fm;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
* C) D! C$ T( u4 m2 t% Oci = nlparci(k,residual,jacobian);
$ B1 I% m9 J) p" G. P8 mfprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
$ d, @  z/ B/ i( x7 N9 Pfprintf('\tk4 = %.11f\n',k(4))
9 K4 |& x) b# |- {! V( {4 Mfprintf('\tk5 = %.11f\n',k(5))
5 X1 h& o: i9 T- W- q- K2 bfprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
. N8 Z" V: i  Ufprintf('\tk8 = %.11f\n',k(8))
fprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
4 [* e9 a' a/ @# B, C1 X2 j( _# Pfprintf('  The sum of the squares is: %.1e\n\n',resnorm)
7 O+ R% [/ |  }* Z: @k_fmls = k;
9 `; T0 x* w. b# q% X; B: C1 ]& q0 Noutput
% K- V9 i* M/ z1 htspan = [0 15 30 45 60 90 120 180 240 300 360];
9 k2 y; E- m# J" G[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
plot(t,x(:,1),t,yexp(:,2),'*');legend('Glc-pr','Glc-real')
" \4 }) B8 l! B6 {- z/ ~" c: x' [figure;plot(t,x(:,2:5));
$ X& e0 ?# W; P( o5 H/ Kp=x(:,1:5)
/ \5 v4 Z: {, U" Q. t! Bhold on
& t* d3 H7 N2 yplot(t,yexp(:,3:6),'o');legend('Fru-pr','Fa-pr','La-pr','HMF-pr','Fru-real','Fa-real','La-real','HMF-real')
* w( _& v& _* `$ G, ?; Q
# z0 a: H0 G. A+ l& }( ?; _9 {
0 q+ K8 x; \7 C
function f = ObjFunc7LNL(k,x0,yexp)
tspan = [0 15 30 45 60 90 120 180 240 300 360];
3 Z. Q: E  P" o8 z) e$ z  ?) |[t, x] = ode45(@KineticEqs,tspan,x0,[],k);   
y(:,2) = x(:,1);
) s& s! ]* X: `1 e% M; r  ^y(:,3:6) = x(:,2:5);
$ @5 F' V* k+ c6 A  n  j- Rf1 = y(:,2) - yexp(:,2);
f2 = y(:,3) - yexp(:,3);
$ I4 d! c6 D0 sf3 = y(:,4) - yexp(:,4);
" B* ^' V, p" K$ df4 = y(:,5) - yexp(:,5);
f5 = y(:,6) - yexp(:,6);
f = [f1; f2; f3; f4; f5];

# D; K% H, C' Y" b  z

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);   
3 G5 u4 g/ f0 ?y(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
/ c- f1 I6 y5 T! E9 y) [f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2)   ...
. n- R! V" y: a, n6 h" X! }/ g6 o    + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2)   ...
+ ~$ F6 J5 I) I    + sum((y(:,6)-yexp(:,6)).^2) ;
7 [. A9 g' V$ U, e% s8 K



) A; }5 ^3 E, z8 S5 j1 P+ [function dxdt = KineticEqs(t,x,k)
8 K+ \6 f% k3 r/ g8 d3 ]dGldt = k(1)*x(2)-(k(2)+k(3)+k(8))*x(1);
" v& a) q- F8 R$ I+ o6 X3 K9 cdFrdt = 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);
dHmdt = k(4)*x(2)-(k(6)+k(7)+k(10))*x(5);
0 ]% ~. ^$ y. u  v; u: }dxdt = [dGldt; dFrdt; dFadt; dLadt; dHmdt];