帮忙做下统计显著性检验和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);
" D7 c4 @2 \1 j0 b* ?+ l1 ?% dFrudt = k1*C(Glc)-(k-1+k3+k4)C(Fru);
# p. V3 x. w! R( V2 Z% dFadt = k(2)*C(Glc)+k4*C(Fru)+（k6+k7)*C(Hmf);
% dLadt = k(7)*C(Hmf);
%dHmfdt = k(3)*C(Fru)-（k6+k7)*C(Hmf);
clear all
2 O5 `) K& u$ Z# d. Jclc
! ^8 Y( c( d- K0 Y: @& @format long
, L+ Q+ u- `2 [4 ?, N%        t/min   Glc    Fru        Fa   La   HMF/ mol/L
# t7 y1 x' D  g6 t: w; w  Kinetics=[0    0.25    0           0    0       0
4 ?5 B: G4 R4 \' P/ \          15    0.2319    0.01257    0.0048    0    2.50E-04
          30    0.19345    0.027    0.00868    0    7.00E-04
' H1 j; h8 E1 O- Y% J' I- h' _          45    0.15105    0.06975    0.02473    0    0.0033
          60    0.13763    0.07397    0.02615    0    0.00428
# C, u& H2 H9 x" p          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
3 M. {0 k* T9 E' _          240    0.03653    0.06488    0.08945    0.01828    0.05452
/ B1 \5 Z3 }; L$ S8 M- L7 C2 j, p, g          300    0.02738    0.05448    0.09098    0.0227    0.0597
$ K* x& K; F3 F! l; H! {9 f& W3 t          360    0.01855    0.04125    0.09363    0.0239    0.06495];
& a# V$ a, k2 bk0 = [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];                  % 参数下限
ub = [1  1  1  1  1  1  1  1  1  1];    % 参数上限
/ O: A! n7 q9 H/ w1 p- Y+ F0 wx0 = [0.25  0  0  0  0];
yexp = Kinetics;                 % yexp: 实验数据[x1        x4        x5        x6]
% warning off
5 m' P- f6 C$ H9 O2 L/ e) H% 使用函数 ()进行参数估计
! H* z; N. ?/ X[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))
/ B5 Y6 v( L8 N3 Z) kfprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
2 v( P1 _$ K5 v5 ?# v9 sfprintf('\tk6 = %.11f\n',k(6))
$ s* x5 X" r* w' Y1 h# rfprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
+ ~# r9 \2 W# q8 }3 n2 pfprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
; G* G# S( S9 H9 Y. g$ n4 lfprintf('  The sum of the squares is: %.1e\n\n',fval)
8 k+ A& I9 V* X1 x2 Hk_fm= k;
% warning off
7 J* x1 c0 P7 S1 V5 C! i% 使用函数lsqnonlin()进行参数估计
  B- F) f  `! |7 |% S[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
7 ^. I! ]0 |" mci = nlparci(k,residual,jacobian);
. D* b8 q  [. z$ [% ?* O% s4 F6 k7 ofprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
- S& v% R& q! l) c4 g! ]fprintf('\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;
& L4 T0 i1 |8 R+ a: H* youtput
0 w8 B: d  a9 d4 w9 gwarning off
% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
+ A# b; C8 j- C& f. a9 uk0 = k_fm;
# U: Z( P$ v* ~$ e& |[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
ci = nlparci(k,residual,jacobian);
fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
# a7 q; }: y: |& o( V/ Nfprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
, y  F) y6 `. V9 V& R5 `: zfprintf('\tk3 = %.11f\n',k(3))
/ `: v  s8 }, n1 U' t9 ^: ~6 N; zfprintf('\tk4 = %.11f\n',k(4))
8 N! a  `0 U* M& H0 F# }fprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
7 x: ?( H) f; D; \% Ofprintf('\tk8 = %.11f\n',k(8))
' u5 P6 c& g6 s6 H) ufprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
' B; F+ H" t% u" e' O  g3 [fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
( ?' a& G8 T8 S* Z2 I( S5 ck_fmls = k;
0 J$ T( x% T2 z. V+ Moutput
, ]" J( {( M7 x- j+ xtspan = [0 15 30 45 60 90 120 180 240 300 360];
/ Q4 ?  N3 X/ F+ y. l; I: U9 W[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
  J9 K9 }" J; E3 _! dplot(t,x(:,1),t,yexp(:,2),'*');legend('Glc-pr','Glc-real')
( @  d$ B1 K1 Z. ]1 hfigure;plot(t,x(:,2:5));
p=x(:,1:5)
hold on
1 j2 z4 R5 b$ s; R% h/ Oplot(t,yexp(:,3:6),'o');legend('Fru-pr','Fa-pr','La-pr','HMF-pr','Fru-real','Fa-real','La-real','HMF-real')



, u0 T3 s2 `0 A/ M8 [5 F  \function f = ObjFunc7LNL(k,x0,yexp)
tspan = [0 15 30 45 60 90 120 180 240 300 360];
[t, x] = ode45(@KineticEqs,tspan,x0,[],k);   
; P9 B+ D2 M3 x* |y(:,2) = x(:,1);
y(:,3:6) = x(:,2:5);
- E& ?( K: ~) p! v( u$ C; Vf1 = y(:,2) - yexp(:,2);
. A% r! S1 E, x1 i$ c# j+ ^- ~. tf2 = y(:,3) - yexp(:,3);
f3 = y(:,4) - yexp(:,4);
( k( u) a1 K6 [' Sf4 = y(:,5) - yexp(:,5);
f5 = y(:,6) - yexp(:,6);
f = [f1; f2; f3; f4; f5];
' ^0 a$ b9 w( e* |! o

% K3 Q* z! y( Y" r* g
6 G+ H! }; i/ efunction 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);
: |# Q; N7 }5 \3 |y(:,3:6) = x(:,2:5);
f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2)   ...
; ?9 F% }. l8 C" j4 z$ k. Y- _* o    + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2)   ...
    + sum((y(:,6)-yexp(:,6)).^2) ;




function dxdt = KineticEqs(t,x,k)
dGldt = k(1)*x(2)-(k(2)+k(3)+k(8))*x(1);
# a7 N4 m6 c  z# zdFrdt = k(2)*x(1)-(k(1)+k(4)+k(5)+k(9))*x(2);
' W$ d5 t% G9 k+ {5 v- P- e; Q3 EdFadt = k(3)*x(1)+k(5)*x(2)+(k(6)+k(7))*x(5);
1 d, V+ v2 T& _& C' T$ j/ S* jdLadt = k(7)*x(5);
dHmdt = k(4)*x(2)-(k(6)+k(7)+k(10))*x(5);
dxdt = [dGldt; dFrdt; dFadt; dLadt; dHmdt];

4 u0 C$ [9 ]- T2 k( ~