24小时热门版块排行榜    

查看: 2082  |  回复: 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 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 EZ67DG6CvcKL 2026-08-31 5/250 2026-09-01 05:43 by o3EMxHUAQDtS
[公派出国] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 wsYGJnURYvm6 2026-08-31 3/150 2026-09-01 01:52 by o3EMxHUAQDtS
[公派出国] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 EZ67DG6CvcKL 2026-08-31 4/200 2026-09-01 01:30 by o3EMxHUAQDtS
[基金申请] 投票:  有多少人是今天查系统知道结果的? +17 爱看书的可乐 2026-08-26 19/950 2026-08-31 21:30 by xiangy672
[基金申请] 基金未中,这种答复是模板吗? +6 zhaosm1982 2026-08-27 7/350 2026-08-31 21:18 by qdxxmc
[基金申请] 哪位高人中了,把查询到的截图贴出来让我看看,让我长长见识 +6 yuleib84 2026-08-26 7/350 2026-08-31 19:46 by 鱼翔浅底1
[文学芳草园] 梦想 +5 myrtle 2026-08-26 5/250 2026-08-31 14:40 by hahaboy
[基金申请] 基金不中,共勉 +12 eulota 2026-08-26 12/600 2026-08-31 08:42 by ZJTJZ
[基金申请] 有没有仍没收到信息的 +7 德尚中行 2026-08-27 8/400 2026-08-30 20:52 by purplejack
[基金申请] 麻烦专家们看看评委们的意见(F口面上) +6 gdd2018 2026-08-28 11/550 2026-08-30 08:58 by 超级无敌华子
[基金申请] 我就是申请一个面上项目而已,这评审意见是按照杰青的条件评的吧? +6 gouxfjh 2026-08-28 11/550 2026-08-30 07:57 by gouxfjh
[基金申请] 国自然面上复盘~欢迎讨论 (金币+15) +15 晴天加油 2026-08-26 16/800 2026-08-29 18:28 by symmetry
[基金申请] 国自然评审意见 +13 wangmingqi 2026-08-28 19/950 2026-08-29 10:22 by Poppy1104
[基金申请] 基金系统什么内容也没有 30+4 winsaint 2026-08-27 9/450 2026-08-28 11:06 by maolC
[基金申请] 怎么查啊 +6 huang1991js 2026-08-26 6/300 2026-08-28 08:42 by winsaint
[基金申请] 面上合作单位盖章 +5 ssyjh 2026-08-27 5/250 2026-08-27 20:50 by gdfollow
[基金申请] 出来了 +9 trojank 2026-08-26 9/450 2026-08-26 14:25 by 宝贝虫子
[基金申请] 国合里面能看到了 +7 一怀馨秋 2026-08-26 7/350 2026-08-26 11:23 by zhaosm1982
[基金申请] 牛来!米来!面来! +8 beefly 2026-08-26 8/400 2026-08-26 08:37 by xuzhipiao
[基金申请] 某些机构,以效率低为荣,以效率低作为存在感 +9 yuleib84 2026-08-25 10/500 2026-08-25 17:14 by alexon
信息提示
请填处理意见