24小时热门版块排行榜    

查看: 2052  |  回复: 6

mljever

金虫 (著名写手)

[求助] 请教Matlab高手解多元非线性方程组中的参数

方程如下:
dXp/dt=μm*S*Xp/(Ks+S)-k1*Xp                   (1)
dXs/dt= k1*Xp-kd*Xs                                    (2)
X=Xp+Xs                                                       (3)
dX/dt=μm*S*Xp/(Ks+S)-kd*Xs                     (4)
dEt/dt=k2*Xs-k3*Et                                      (5)
-dS/dt=(k4*(1-x)^k5+k6)*0.9062*S*Ef/(1+0.9062*Ef)          (6)
其中式(4)是由前面(1)-(3)推导得到
方程组中X, Et, Ef, S, x为已知数据的变量(见下表),μm, kd, Ks, k1, k2, k3, k4, k5, k6为待拟合求解的参数。请高人们给帮帮忙吧,冰天雪地里跪谢啦!
t(h)        0        6        12        18        24        30        37        42        48        54        61        66        72        78        84        90        96        102        108        114        120        126        132
X        2.76        3.93        6.01        8.58        11.25        12.29        13.42        18.16        22.55        22.53        22.52        20.81        17.95        17.94        17.94        17.29        16.63        16.35        15.98        14.21        12.78        12.43        11.96
S        50        48.66        46.2        42.78        38.3        32.25        26.4        26.13        25.72        23.28        21.69        20.55        22.4        19.83        19.07        19.26        20.38        19.13        19.06        18.18        17.24        16.93        16.56
Ef        0        0.025        0.076        0.129        0.198        0.95        1.89        2.62        3.46        4.58        5.72        6.38        6.91        7.27        7.58        9.05        10.63        10.66        10.48        10.59        10.69        11.72        12.25
Et        0        0.046        0.133        0.215        0.31        1.237        2.21        2.973        3.835        4.94        6.069        6.717        7.281        7.601        7.9        9.38        10.985        10.993        10.811        10.906        10.99        12.017        12.542
x        0        0.0268        0.076        0.1444        0.234        0.355        0.472        0.4774        0.4856        0.5344        0.5662        0.589        0.552        0.6034        0.6186        0.6148        0.5924        0.6174        0.6188        0.6364        0.6552        0.6614        0.6688
回复此楼

» 收录本帖的淘帖专辑推荐

感兴趣的问题

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dbb627

荣誉版主 (著名写手)

【答案】应助回帖

mljever(金币+10): 谢谢版主 2011-11-23 10:00:30
Xp=dsolve('DXp=um*S*Xp/(Ks+S)-k1*Xp','t')

Xp =

C2/exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S))

X=dsolve('DX=um*S*(C2/exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S)))/(Ks+S)-kd*(X-(C2/exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S))))','X(0)=2.76','t') ;%替换Xp Xs

X=simplify(X)



X=
(exp((S*t*um)/(Ks + S) - (S*k1*t)/(Ks + S) - (Ks*k1*t)/(Ks + S))*(C2*Ks*kd + C2*S*kd + C2*S*um))/(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um) - ((C2*Ks*kd + C2*S*kd + C2*S*um)/(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um) - 69/25)/exp(kd*t)



Xs=simple(X-Xp)

Xs =
(C2*(Ks*kd + S*kd + S*um))/(exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S))*(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um)) - C2/exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S)) - ((C2*(Ks*kd + S*kd + S*um))/(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um) - 69/25)/exp(kd*t)

Et=dsolve('DEt=k2*((C2*(Ks*kd + S*kd + S*um))/(exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S))*(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um)) - C2/exp((t*(Ks*k1 + S*k1 - S*um))/(Ks + S)) - ((C2*(Ks*kd + S*kd + S*um))/(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um) - 69/25)/exp(kd*t))-k3*Et','Et(0)=0','t')

Et =simple(Et)
Et=
(k2*(69*Ks*k1 - 69*Ks*k3 + 69*S*k1 - 69*S*k3 - 69*S*um + 25*C2*Ks*k3 + 25*C2*S*k3 + 25*C2*S*um))/(25*exp(k3*t)*(k3 - kd)*(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um)) - ((k2*exp(t*(k3 - kd))*(69*Ks*k1 - 69*Ks*kd + 69*S*k1 - 69*S*kd - 69*S*um + 25*C2*Ks*kd + 25*C2*S*kd + 25*C2*S*um))/(k3 - kd) - (25*C2*k1*k2*exp((t*(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um))/(Ks + S))*(Ks + S)^2)/(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um))/(25*exp(k3*t)*(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um))


这里x,Ef应该是关于t的函数
你没给出关系是所以dS/dt不能积分求表达式,需要用数值微分求dS/dt 再根据x,Ef和S的数据用拟合的方法求k4,k5,k6

其他X Et中的参数μm, kd, Ks, k1, k2, k3的可直接用拟合方法求当然也可采用插值的方法
The more you learn, the more you know, the more you know, and the more you forget. The more you forget, the less you know. So why bother to learn.
2楼2011-11-04 22:05:29
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mljever

金虫 (著名写手)

谢谢版主,x,Ef是关于t的函数,但是不同时间点都有对应的数据,这个还需要公式是吗?
3楼2011-11-23 09:59:33
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dbb627

荣誉版主 (著名写手)

【答案】应助回帖

引用回帖:
3楼: Originally posted by mljever at 2011-11-23 09:59:33:
谢谢版主,x,Ef是关于t的函数,但是不同时间点都有对应的数据,这个还需要公式是吗?

不同时间点都有对应的数据,这只是个时间序列,没有确定的表达式啊,不能带入演算解析式,只能是算数值的
The more you learn, the more you know, the more you know, and the more you forget. The more you forget, the less you know. So why bother to learn.
4楼2011-11-23 11:02:28
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mljever

金虫 (著名写手)

引用回帖:
4楼: Originally posted by dbb627 at 2011-11-23 11:02:28:
不同时间点都有对应的数据,这只是个时间序列,没有确定的表达式啊,不能带入演算解析式,只能是算数值的

x=(S0-S)/S0  由于式中S是t的函数,故x为t的函数

Ef=Et-Ead, 其中Ead=(S/10^3)*Amax*Ef*Kad/(1+Ef*Kad)
                其中参数Amax, Kad由式A=Amax*Kad*E/(1+Kad*E)拟合得到
                数据:
         E   0.039   0.082   0.346   0.653   2.056   2.202   4.613   5.781
              A   0.536   0.352   6.347   7.686   9.352   10.94   15.337 16.329
5楼2011-11-23 13:34:19
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dbb627

荣誉版主 (著名写手)

【答案】应助回帖

mljever(金币+20): 版主太厉害啦,我再计算一下,如有问题再向您请教哈 2011-11-24 08:55:49
引用回帖:
5楼: Originally posted by mljever at 2011-11-23 13:34:19:
x=(S0-S)/S0  由于式中S是t的函数,故x为t的函数

Ef=Et-Ead, 其中Ead=(S/10^3)*Amax*Ef*Kad/(1+Ef*Kad)
                其中参数Amax, Kad由式A=Amax*Kad*E/(1+Kad*E)拟合得到
                数据:
   ...

General model:
     f(x) = Amax*Kad*x/(1+Kad*x)
Coefficients (with 95% confidence bounds):
       Amax =        18.1  (12.78, 23.41)
       Kad =      0.9062  (0.05316, 1.759)

Goodness of fit:
  SSE: 14.79
  R-square: 0.9408
  Adjusted R-square: 0.9309
  RMSE: 1.57

x=(S0-S)/S0  由于式中S是t的函数,故x为t的函数
Ef=Et-Ead, 其中Ead=(S/10^3)*18.1*Ef*0.9062/(1+Ef*0.9062)
Ef=Et-(S/10^3)*18.1*Ef*0.9062/(1+Ef*0.9062)--
-->解出Ef 用Et和S表达
Et=
(k2*(69*Ks*k1 - 69*Ks*k3 + 69*S*k1 - 69*S*k3 - 69*S*um + 25*C2*Ks*k3 + 25*C2*S*k3 + 25*C2*S*um))/(25*exp(k3*t)*(k3 - kd)*(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um)) - ((k2*exp(t*(k3 - kd))*(69*Ks*k1 - 69*Ks*kd + 69*S*k1 - 69*S*kd - 69*S*um + 25*C2*Ks*kd + 25*C2*S*kd + 25*C2*S*um))/(k3 - kd) - (25*C2*k1*k2*exp((t*(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um))/(Ks + S))*(Ks + S)^2)/(Ks*k3 - Ks*k1 - S*k1 + S*k3 + S*um))/(25*exp(k3*t)*(Ks*kd - Ks*k1 - S*k1 + S*kd + S*um))
x-->(S0-S)/S0  
替换下式中Ef和x 解出S的解析式 在拟合即可 只不过表达是太复杂了
-dS/dt=(k4*(1-x)^k5+k6)*0.9062*S*Ef/(1+0.9062*Ef)
The more you learn, the more you know, the more you know, and the more you forget. The more you forget, the less you know. So why bother to learn.
6楼2011-11-23 18:21:02
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

微分方程拟合可以用1stOpt来解决,很方便好用。
7楼2011-12-05 11:35:34
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 mljever 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 应该是93bebmhtak前后十一个字符比较关键 +16 Lanmanbaby 2026-08-09 25/1250 2026-08-11 15:51 by 医学老男孩
[论文投稿] 职称评审,求友友推荐见刊最快的期刊 +3 工厂打螺丝 2026-08-08 3/150 2026-08-11 13:42 by 会议编辑
[基金申请] 听说今天filecode变了 +26 布布和一二 2026-08-06 49/2450 2026-08-11 13:18 by WH3796
[基金申请] FileCode能看出啥? +7 要乐观耀哥 2026-08-10 17/850 2026-08-11 10:50 by 要乐观耀哥
[基金申请] filecode +9 documentary 2026-08-10 9/450 2026-08-11 10:08 by wk7465
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +8 家与远方 2026-08-10 13/650 2026-08-11 09:49 by 家与远方
[基金申请] 基金中了 +15 laoda193707 2026-08-06 15/750 2026-08-11 00:11 by jiafei2190
[基金申请] 国自然结果 +4 Vierhys 2026-08-10 8/400 2026-08-10 15:06 by Vierhys
[基金申请] 据悉今年马上要出结果了 +7 瞬息宇宙 2026-08-10 8/400 2026-08-10 12:42 by Vivilian
[基金申请] 2026国自然放榜时间 +9 布布和一二 2026-08-08 9/450 2026-08-10 11:22 by xxxx2020
[基金申请] 这样的filecode谁见过 +11 布布和一二 2026-08-08 22/1100 2026-08-10 11:10 by wmfsnow
[基金申请] filecode与中标关系的预测 +5 布布和一二 2026-08-07 5/250 2026-08-09 16:15 by 袁向阳007
[基金申请] fileCode有新解读? +10 Tide man 2026-08-08 18/900 2026-08-09 12:55 by 仁砚薪传
[基金申请] 关于filecode,很负责任的告诉大家 +6 爱看书的可乐 2026-08-08 7/350 2026-08-08 22:13 by a_niu
[基金申请] 关于filecode +4 布布和一二 2026-08-07 7/350 2026-08-07 22:55 by zhanghaozhu
[基金申请] 化学口download_prp&fileCode的固定段好像这几天一直没变,有变的大神么? +3 Tide man 2026-08-07 4/200 2026-08-07 22:39 by Tide man
[基金申请] 固定端突然变了,今天 +6 archvillain 2026-08-06 10/500 2026-08-07 16:03 by 医学老男孩
[基金申请] filecode变化情况 +6 布布和一二 2026-08-07 22/1100 2026-08-07 14:45 by 且听虎啸
[基金申请] 大家散了吧,后缀研究没有意义,别浪费时间了,过好目前的每一天,不要焦虑 +5 Tide man 2026-08-06 7/350 2026-08-07 13:11 by 医学老男孩
[基金申请] 求各位大神看下 100+6 hpkpkpkp 2026-08-05 33/1650 2026-08-06 14:49 by zhiyanjiang
信息提示
请填处理意见