# How do I fit coupled differential equations to experimental data?

9 views (last 30 days)
CP on 13 Nov 2015
Commented: Star Strider on 18 Nov 2015
I have a set of differential equations that I would like to fit to experimental data to get chemical kinetics parameters out of. The equations are:
d(Mzb(t))/dt = -[Mzb(t)-Mzb0]/T1b - kb*Mzb(t) + kf*Mzf(t);
d(Mzf(t))/dt = -[Mzf(t)-Mzf0]/T1f - kf*Mzf(t) + kb*Mzb(t);
The data I have is of Mzb(t) and Mzf(t) at different time points. From that data, I can glean Mzb0 and Mzf0 (the equilibrium values Mzb(t) and Mzf(t), which will be constants), and from other data I know the constants T1b and T1f. I also know that this will be a least-squares optimization to get the best fit for the kb and kf parameters. The part I don't understand is how to code all of this-- I'm a novice with the syntax and would appreciate any and all help about it.

Star Strider on 13 Nov 2015
Two possibilities: Monod Kinetics and Curve Fitting and Optimizing a Simulation or Ordinary Differential Equation. Those should give you the information you need.
##### 2 CommentsShowHide 1 older comment
Star Strider on 18 Nov 2015
My pleasure.
First, there are several ways to be sure your ODE function returns a column vector. Probably the easiest is to add a second dimension reference to your equations:
xdot(1,:) = -(x(1)-B(1))/T1b - B(3).*x(1) + B(4).*x(2);
xdot(2,:) = -(x(2)-B(2))/T1f - B(4).*x(2) + B(3).*x(1);
Those will return a column vector. Another way is to have as the first assignment statement in your function:
xdot = zeros(2,1);
choose the one you like best!
Note that you can pass several parameters to your ODE function, including the initial conditions, that will produce the best fit to your data. You just have to be sure to parse them out appropriately within your objective function (here ‘Mzbt’).
Second, I discourage your using global variables. It is best to pass them as parameters in your function:
function S = Mzbt(B,t,T1b,T1f,Mzb0,Mzf0)
Then in your script file that calls lsqnonlin or nlinfit, use an anonymous function to be sure the fitting function only ‘sees’ the arguments of the function you want it to:
objfcn = @(B,t) Mzbt(B,t,T1b,T1f,Mzb0,Mzf0);
You can use that as your objective function, or include it as the appropraite argument in, for example, lsqcurvefit:
B = lsqcurvefit(@(B,t) Mzbt(B,t,T1b,T1f,Mzb0,Mzf0), ... );
The function will pick up the other arguments from your workspace, since you have already defined them.
This is not a straightforward problem, as you have discovered. I will be here to help you as much as you need.