帮忙做下统计显著性检验和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);
% dFrudt = k1*C(Glc)-(k-1+k3+k4)C(Fru);
+ C5 O. a: _( f  T7 U% dFadt = k(2)*C(Glc)+k4*C(Fru)+（k6+k7)*C(Hmf);
! C% t; n$ k: E, N& C% dLadt = k(7)*C(Hmf);
& a3 d" B2 q: v3 H' s* ]& u%dHmfdt = k(3)*C(Fru)-（k6+k7)*C(Hmf);
clear all
clc
format long
# _. u; e( W: q) v& L4 F3 r& |% M4 y%        t/min   Glc    Fru        Fa   La   HMF/ mol/L
  Kinetics=[0    0.25    0           0    0       0
* g$ i( ?5 ]7 i- n8 h          15    0.2319    0.01257    0.0048    0    2.50E-04
          30    0.19345    0.027    0.00868    0    7.00E-04
          45    0.15105    0.06975    0.02473    0    0.0033
          60    0.13763    0.07397    0.02615    0    0.00428
' w# I# L; f+ s/ a. m8 Z1 ]          90    0.08115    0.07877    0.07485    0    0.01405
5 n" M  \1 X. E6 b. k          120    0.0656    0.07397    0.07885    0.00573    0.02143
- F: k5 {2 [1 L; p          180    0.04488    0.0682    0.07135    0.0091    0.03623
          240    0.03653    0.06488    0.08945    0.01828    0.05452
          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];                  % 参数下限
; g1 Y* g2 G, \3 j6 R# Zub = [1  1  1  1  1  1  1  1  1  1];    % 参数上限
9 k2 Q% k. l' }+ i( C3 o) R9 v' K0 e  Yx0 = [0.25  0  0  0  0];
yexp = Kinetics;                 % yexp: 实验数据[x1        x4        x5        x6]
% warning off
: A; F" y* N0 J" q$ U4 h( {' k% 使用函数 ()进行参数估计
' n9 [7 c! W9 J$ R, @[k,fval,flag] = fmincon(@ObjFunc7Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
, M+ X3 e6 o5 H+ l4 L+ S1 L8 ?fprintf('\n使用函数fmincon()估计得到的参数值为:\n')
) Y/ t2 E( z0 E4 ^: B5 l  xfprintf('\tk1 = %.11f\n',k(1))
fprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
9 e6 K6 {9 w3 g9 I( F; b& Q# z2 f( v( @; Hfprintf('\tk5 = %.11f\n',k(5))
fprintf('\tk6 = %.11f\n',k(6))
. O* t8 j3 N1 N- Tfprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
! u! H3 b% Y( a/ r3 I2 |+ b% D, v! Ffprintf('\tk9 = %.11f\n',k(9))
fprintf('\tk10 = %.11f\n',k(10))
fprintf('  The sum of the squares is: %.1e\n\n',fval)
. T4 R4 |  j# Z; ^- ]1 O5 [k_fm= k;
& C% S2 A  a8 x4 p1 H% warning off
% 使用函数lsqnonlin()进行参数估计
+ p: _1 ~% K1 \6 g. t; y5 B0 \+ s[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
! c% L* M# x+ g3 ^* V. B- f5 S* [3 J    lsqnonlin(@ObjFunc7LNL,k0,lb,ub,[],x0,yexp);      
1 u2 V% g. i4 o* |ci = nlparci(k,residual,jacobian);
1 v8 Y6 e& r9 i( E" s8 Lfprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\tk1 = %.11f\n',k(1))
, ^+ G; _* f9 i) v# `' F, Cfprintf('\tk2 = %.11f\n',k(2))
# c+ D* d8 {# t" yfprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
' Z9 o0 G; |( f, d( X6 h+ gfprintf('\tk5 = %.11f\n',k(5))
' R* o, w" g# E  H6 v/ t1 @( F! ifprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
fprintf('\tk8 = %.11f\n',k(8))
/ n& ~* ?% H( z$ Y9 N6 F8 ]7 ufprintf('\tk9 = %.11f\n',k(9))
7 ]" G$ J0 ]+ N2 s% N" K7 D' tfprintf('\tk10 = %.11f\n',k(10))
fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
k_ls = k;
output
warning off
6 b9 H# ^% Q  u* e6 r% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
$ R' Z) R% K: ik0 = k_fm;
[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')
- F; t# e' Y+ q8 ?, n+ Yfprintf('\tk1 = %.11f\n',k(1))
, n% K0 z) K6 A2 @3 Afprintf('\tk2 = %.11f\n',k(2))
fprintf('\tk3 = %.11f\n',k(3))
fprintf('\tk4 = %.11f\n',k(4))
fprintf('\tk5 = %.11f\n',k(5))
& k3 P# w6 N7 D4 ]: Qfprintf('\tk6 = %.11f\n',k(6))
fprintf('\tk7 = %.11f\n',k(7))
4 n' d5 A$ S1 e, kfprintf('\tk8 = %.11f\n',k(8))
fprintf('\tk9 = %.11f\n',k(9))
/ `$ B0 `$ b6 Y& lfprintf('\tk10 = %.11f\n',k(10))
9 z! c, }" e" |fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
. S6 b, A3 f3 Yk_fmls = k;
output
tspan = [0 15 30 45 60 90 120 180 240 300 360];
" q! D$ h3 H  J) B3 i[t x] = ode45(@KineticEqs,tspan,x0,[],k_fmls);
figure;
7 y4 r5 T/ u+ S6 Q- s/ lplot(t,x(:,1),t,yexp(:,2),'*');legend('Glc-pr','Glc-real')
4 E: s& a! ?: q4 w+ q; B0 nfigure;plot(t,x(:,2:5));
p=x(:,1:5)
, `2 _1 Q5 b# K/ v7 l3 L9 {hold on
plot(t,yexp(:,3:6),'o');legend('Fru-pr','Fa-pr','La-pr','HMF-pr','Fru-real','Fa-real','La-real','HMF-real')



function f = ObjFunc7LNL(k,x0,yexp)
& U' E  s1 r, ^2 g) p2 y  N9 Xtspan = [0 15 30 45 60 90 120 180 240 300 360];
[t, x] = ode45(@KineticEqs,tspan,x0,[],k);   
y(:,2) = x(:,1);
4 H7 _9 f! O7 y2 D, Gy(:,3:6) = x(:,2:5);
f1 = y(:,2) - yexp(:,2);
  {( S: V2 Y2 }* [& K8 Bf2 = y(:,3) - yexp(:,3);
f3 = y(:,4) - yexp(:,4);
f4 = y(:,5) - yexp(:,5);
f5 = y(:,6) - yexp(:,6);
1 o- Y, U% {/ ^f = [f1; f2; f3; f4; f5];



function f = ObjFunc7Fmincon(k,x0,yexp)
4 G" h! L- l6 k; Ltspan = [0 15 30 45 60 90 120 180 240 300 360];
[t x] = ode45(@KineticEqs,tspan,x0,[],k);   
  a2 ~+ d3 F# my(:,2) = x(:,1);
$ f) \5 H% K% }) O8 O$ [0 by(:,3:6) = x(:,2:5);
f =  sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2)   ...
5 w8 z* |' H5 h, G4 [    + sum((y(:,4)-yexp(:,4)).^2) + sum((y(:,5)-yexp(:,5)).^2)   ...
8 Q" d' j1 R: D$ s* H9 ?    + sum((y(:,6)-yexp(:,6)).^2) ;
! q8 Q' ]' v% B  n6 @4 |( B
- {; z8 A2 R3 H; a6 Y7 \- q


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);
/ q. g" f1 I9 R5 O# j/ \dFadt = k(3)*x(1)+k(5)*x(2)+(k(6)+k(7))*x(5);
dLadt = k(7)*x(5);
# C' G# D' l7 j6 b( N  k# EdHmdt = k(4)*x(2)-(k(6)+k(7)+k(10))*x(5);
dxdt = [dGldt; dFrdt; dFadt; dLadt; dHmdt];
  r, L( `$ S# y" z' ~( d

董事长之友

蒙的一比 大哥