MATLAB code for solving system of ODE with unknown parameters
9 views (last 30 days)
Show older comments
Hello, I was wondering if someone can help me write a code in MATLAB for solving system of nonlinear ODE with unknown parameters that need to fit the experimental data. I tried doing it by using < this > , but it didn't work. Or maybe I didn't do it properly because I didn't really understand what was what and how the variables were defined. My problem include this system:
- dP/dt = (b - r0(N))P + ri(N)Q ;
- dQ/dt = r0(N)P - (ri(N) + μ)Q ;
- dD/dt = μQ - dD;
where
- P,Q,D are function in t (time) ;
- N = P + Q + D ;
- r0(N) = kN/(aN + 1);
- ri(N) = r/(N + m);
and a, b, d, k, m, μ and r are unknown parameters (that can be estimated within an interval, for every parameter; b, μ, d are positive). Initial conditions are:
- P(0) = 0.1 (or 0.01);
- Q(0) = D(0) = 0;
My problem includes finding approximative solution (or just in a form that can be graphed) and fitting parameters to experimental data that includes values for N(t) in some period of time (and then making a plot of P, Q, D, N functions and data on the same graph), for example:
t N(t)
4,6 0,0016308;
5,7 0,0032148;
6,7 0,0056140;
7,9 0,0118598;
8,8 0,0201500;
9,8 0,0275380;
11,0 0,0345460;
12,0 0,0660800;
12,9 0,0789320;
15,2 0,1550000;
Can someone explain me how to work this out in MATLAB? Please start from scratch, write me everything I need to do (like how to define P, Q, D and parameters, what goes in what function, how to use toolboxes and such). I am sorry if this is too much from me to ask, but I struggle with this for weeks now.
Thank you very much.
0 Comments
Answers (1)
Alan Weiss
on 28 Oct 2014
I am not going to give you a step-by-step tutorial on how to formulate ODE for solution in MATLAB. I will tell you that you have an optimization problem.
Set the vector
x = [a, b, d, k, m, μ, r]
I mean, a = x(1), b = x(2), etc. Then you can take a guess x(0) as initial parameters, and try to solve the optimization problem of minimizing the sum of squares of differences of your ODE solution from the data. The optimizer will adjust x to try to achieve the minimal difference.
For examples, see Curve Fitting via Optimization, lsqnonlin with a Simulink model, and nonlinear curve fitting with lsqcurvefit.
For solving an ODE given an initial point, see Ordinary Differential Equations (a very long section, but one that has a bunch of examples).
Good luck,
Alan Weiss
MATLAB mathematical toolbox documentation
3 Comments
Torsten
on 5 Nov 2014
In the sequel, exp_t and exp_y are your measurement data.
a = 1.012;
b = 0.35;
d = 0.02085;
k = 1.25;
m = 3.42;
mi = 0.058;
r = 1.2025;
p0(1)=a;
p0(2)=b;
p0(3)=d;
p0(4)=k;
p0(5)=m;
p0(6)=mi;
p0(7)=r;
p_estimate = lsqnonlin(@(p)odefit(exp_t,exp_y,p),p0);
function err = odefit(exp_t,exp_y,p)
a=p(1);
b=p(2);
d=p(3);
k=p(4);
m=p(5);
mi=p(6);
r=p(7);
y0 = [1.0e-003;0;0];
[t,y] = ode45(@(t,y)f(t,y,a,b,d,k,m,mi,r),exp_t,y0)
err = y-exp_y; % compute error between experimental y and fitted y
end
function dydt = f(t,y,a,b,d,k,m,mi,r)
dydt = [ (b - k*(y(1) + y(2) + y(3))/(a*(y(1) + y(2) + y(3)) + 1))*y(1) + r*y(2)/(y(1) + y(2) + y(3) + m) k*(y(1) + y(2) + y(3))*y(1)/(a*(y(1) + y(2) + y(3)) + 1) - (r/(y(1) + y(2) + y(3) + m) + mi)*y(2) mi*y(2) - d*y(3)]; end
Torsten
on 5 Nov 2014
err = y-exp_y
has to be modified since in your model, there are three instead of one y(i)'s.
I leave it to you to make the necessary changes.
Best wishes
Torsten.
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!