Fitting an ODE set with time dependent parameters to data

7 次查看(过去 30 天)
I have this set of ODE that I would like to fit to measured data of two time-dependent parameters Pa and Pe:
Both equations decrease exponentially as a function of time. To solve these equations I understood this needs to be programmed differently by integrating the exponential terms on forehand, while to fitting ODE to data is pretty straight forward.
Is it possible to do both at the same time? I have been breaking my head over this for a few days and I hope anyone here could help me out.

采纳的回答

Star Strider
Star Strider 2020-9-30
You could certainly estimate it by integrating it with ode45, however it has an analytic solution (kindly provided by the Symbolic Math Toolbox):
syms Pe(t) Pa(t) Pa0 Pe0
a = sym('a', [1 6]);
DEq1 = diff(Pe) == a(1)*exp(a(2)*t)-a(3) * (Pe-Pa+1)
DEq2 = diff(Pa) == a(4)*exp(a(5)*t)-a(6) * (Pe-Pa+1)
Eqs = dsolve(DEq1,DEq2, Pa(0)==Pa0, Pe(0)==Pe0)
Pa = simplify(Eqs.Pa, 'Steps',500)
Pe = simplify(Eqs.Pe, 'Steps',500)
Eqs = matlabFunction([Pa;Pe], 'Vars',{[a Pa0, Pe0],t})
that after a bit of editing resolves to:
Eqs = @(in1,t)[-exp(-t.*(in1(:,3)-in1(:,6))).*((in1(:,6).*exp(t.*(in1(:,3)-in1(:,6)))-(in1(:,1).*in1(:,6).*exp(t.*(in1(:,2)+in1(:,3)-in1(:,6))))./(in1(:,2)+in1(:,3)-in1(:,6))+(in1(:,4).*in1(:,6).*exp(t.*(in1(:,3)+in1(:,5)-in1(:,6))))./(in1(:,3)+in1(:,5)-in1(:,6)))./(in1(:,3)-in1(:,6))-(in1(:,6).*(-in1(:,1).*in1(:,3)+in1(:,2).*in1(:,3)-in1(:,1).*in1(:,5)+in1(:,2).*in1(:,4)+in1(:,1).*in1(:,6)+in1(:,2).*in1(:,5)+in1(:,3).*in1(:,4)-in1(:,2).*in1(:,6)+in1(:,3).*in1(:,5)-in1(:,3).*in1(:,6).*2.0-in1(:,4).*in1(:,6)-in1(:,5).*in1(:,6)-in1(:,7).*in1(:,3).^2-in1(:,7).*in1(:,6).^2+in1(:,8).*in1(:,3).^2+in1(:,8).*in1(:,6).^2+in1(:,3).^2+in1(:,6).^2-in1(:,7).*in1(:,2).*in1(:,3)-in1(:,7).*in1(:,2).*in1(:,5)+in1(:,7).*in1(:,2).*in1(:,6)-in1(:,7).*in1(:,3).*in1(:,5)+in1(:,7).*in1(:,3).*in1(:,6).*2.0+in1(:,7).*in1(:,5).*in1(:,6)+in1(:,8).*in1(:,2).*in1(:,3)+in1(:,8).*in1(:,2).*in1(:,5)-in1(:,8).*in1(:,2).*in1(:,6)+in1(:,8).*in1(:,3).*in1(:,5)-in1(:,8).*in1(:,3).*in1(:,6).*2.0-in1(:,8).*in1(:,5).*in1(:,6)))./((in1(:,3)-in1(:,6)).*(in1(:,2)+in1(:,3)-in1(:,6)).*(in1(:,3)+in1(:,5)-in1(:,6))))-(in1(:,1).*in1(:,5).*in1(:,6).*exp(in1(:,2).*t)-in1(:,2).*in1(:,3).*in1(:,4).*exp(in1(:,5).*t))./(in1(:,2).*in1(:,5).*(in1(:,3)-in1(:,6)))-(in1(:,2).*in1(:,3).*in1(:,4)-in1(:,1).*in1(:,5).*in1(:,6)-in1(:,7).*in1(:,2).*in1(:,3).*in1(:,5)+in1(:,8).*in1(:,2).*in1(:,5).*in1(:,6))./(in1(:,2).*in1(:,5).*(in1(:,3)-in1(:,6)));
-(in1(:,1).*in1(:,5).*in1(:,6).*exp(in1(:,2).*t)-in1(:,2).*in1(:,3).*in1(:,4).*exp(in1(:,5).*t))./(in1(:,2).*in1(:,5).*(in1(:,3)-in1(:,6)))-(in1(:,3).*exp(-t.*(in1(:,3)-in1(:,6))).*((in1(:,6).*exp(t.*(in1(:,3)-in1(:,6)))-(in1(:,1).*in1(:,6).*exp(t.*(in1(:,2)+in1(:,3)-in1(:,6))))./(in1(:,2)+in1(:,3)-in1(:,6))+(in1(:,4).*in1(:,6).*exp(t.*(in1(:,3)+in1(:,5)-in1(:,6))))./(in1(:,3)+in1(:,5)-in1(:,6)))./(in1(:,3)-in1(:,6))-(in1(:,6).*(-in1(:,1).*in1(:,3)+in1(:,2).*in1(:,3)-in1(:,1).*in1(:,5)+in1(:,2).*in1(:,4)+in1(:,1).*in1(:,6)+in1(:,2).*in1(:,5)+in1(:,3).*in1(:,4)-in1(:,2).*in1(:,6)+in1(:,3).*in1(:,5)-in1(:,3).*in1(:,6).*2.0-in1(:,4).*in1(:,6)-in1(:,5).*in1(:,6)-in1(:,7).*in1(:,3).^2-in1(:,7).*in1(:,6).^2+in1(:,8).*in1(:,3).^2+in1(:,8).*in1(:,6).^2+in1(:,3).^2+in1(:,6).^2-in1(:,7).*in1(:,2).*in1(:,3)-in1(:,7).*in1(:,2).*in1(:,5)+in1(:,7).*in1(:,2).*in1(:,6)-in1(:,7).*in1(:,3).*in1(:,5)+in1(:,7).*in1(:,3).*in1(:,6).*2.0+in1(:,7).*in1(:,5).*in1(:,6)+in1(:,8).*in1(:,2).*in1(:,3)+in1(:,8).*in1(:,2).*in1(:,5)-in1(:,8).*in1(:,2).*in1(:,6)+in1(:,8).*in1(:,3).*in1(:,5)-in1(:,8).*in1(:,3).*in1(:,6).*2.0-in1(:,8).*in1(:,5).*in1(:,6)))./((in1(:,3)-in1(:,6)).*(in1(:,2)+in1(:,3)-in1(:,6)).*(in1(:,3)+in1(:,5)-in1(:,6)))))./in1(:,6)-(in1(:,2).*in1(:,3).*in1(:,4)-in1(:,1).*in1(:,5).*in1(:,6)-in1(:,7).*in1(:,2).*in1(:,3).*in1(:,5)+in1(:,8).*in1(:,2).*in1(:,5).*in1(:,6))./(in1(:,2).*in1(:,5).*(in1(:,3)-in1(:,6)))];
where ‘in1’ are the parameters ‘a(1)’ to ‘a(6)’ in order, and the intital conditions ‘Pa(0)’ and ‘Pe(0)’ (also in that order) as a row vector, so the initial parameter estimates must be a row vector or the function will throw an error. (This is a pecularity of the functions created by the matlabFunction function.) The data you are fitting must be in row vectors as well, since they will be fitting [Pa; Pe] as row vectors, with ‘Pa’ the first row and ‘Pe’ the second row. I would use lsqcurvefit with these.
To clarify, unless I am missing something, these are time dependent equations, and do not appear to have time-dependent parameters, since the parameters do not appear to vary with time.
  4 个评论
Jos Rozema
Jos Rozema 2020-10-1
编辑:Jos Rozema 2020-10-1
Thank you for the suggestions. That will be a useful starting point. :)

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by