How should I approach solving differential equation as function of two variables r, k?
13 views (last 30 days)
I have a dataset containing two variables, time (hours) and colony size of yeast growth (10^6 cells/mL). I also have the differential equation dy/dt = r*y(t)*(1-y(t)/k) and have this so far:
function dydt = deq(t, y)
dydt = r * y * (1 - (y/k));
And my task is to solve the differential equation as a function of r, k and also compute the least squares error between the solution to the DEQ and dataset. For solving the DEQ, I have this so far:
t0 = 0; %initial time
tf = 40; %final time
y0 = 1; %first colony size at t = 0
h = 0.5; %step size?
k = 1; %initial guesses for constants?
r = 0.5;
tempFile = 'yeastGrowth_01.xlsx'; %dataset given
myData = xlsread(tempFile);
tVec = (t0 : h : tf);
[r, k] = ode45(deq, tVec, y0);
%plot(r, k, 'ko')
Alan Stevens on 6 May 2022
You need to call the ode with something like
[t, y] = ode45(@(t,y) deq(t,y,r,k),tVec,y0,r,k);
where you pass values of r and k to deq and get back the corresponding values of t and y.
Use fminsearch to find the values of r and k that minimize the resulting least squares error - see
John D'Errico on 6 May 2022
While you COULD use an ODE solver, this is something you would ONLY do if the ODE had no analytical solution. And this ODE has a simple solution. That is:
syms r k y(t)
ysol = dsolve(diff(y) == r * y * (1 - (y/k)), y(0) == 1)
ysol = simplify(ysol)
And now you have data that can be used to estimate r and k. As long as you have multiple data points, that will provide an estimate.
Now you could use a tool like the curve fitting toolbox, or you could use a tool from the optimization toolbox, so lsqcurvefit or lsqnonlin, or you could use nlinfit form the stats TB, or you could use a simple optimzer like fminsearch to minimize the sum of squares of errors of your data.
Only if you cannot solve the ODE as I did here would you use a solver and then apply it on top of an ODE solver like ODE45. Why is that something to be avoided if you can avoid it? There are multiple issues I can think of.
- The ODE solver itself is not exact. It has tolerances. So the numbers that come out ot ODE45 may only be accurate to within some tolerance. And that means the optimizer does not see the truly correct values, and that means the objective function it is trying to minimize is not exact. It need not even be a differentiable function And that can cause the optimizer to fail.
- An ODE solver can sometimes have convergence problems for some sets of parameters. And since an optimizer will send in all sorts of parameters into the ODE solver, you might see all sorts of garbage coming out. ODE45 is a classic example of a solver that will simply fail to converge if the problem is stiff. And some problems can be stiff for some sets of parameters, but not others.
- It is inefficient to solve an ODE many times, compared to solving it once, and then using a curve fitting tool on the analytical solution.
Since this is surely homework, and you may have been instructed to use ODE45 to solve the ODE, and then to fit the data to the solution as you want to do, you should do as you were told, IF that is the case. However, you should understand this is not remotely the solution method you SHOULD take for this problem.