Fitting data with ODE equation and multiple inputs

11 次查看(过去 30 天)
Hello,
I have an equation that needs to be fitted for multiple parameters but am genuinely stuck on how to do it. Please see the equation below:
It is important to note that z is a function of x and I have data for z(x) and dy/dx.
I am interested in fitting parameters for a, b, and n, but the problem does not seem entirely straightforward to me in MATLAB and I am hoping to get some guidance. I've looked at ways to try to fit the ODE using the guide written here and even tried looking at some analagous problems with Monod kinetics, but the way that the solutions have been posed for Monod kinetics is not entirely clear to me since they have data for all of the variables.
Below is some example data that can be used to try and fit.
x = [0 0.0258 0.0515 0.0773 0.1031 0.1289 0.1546]
dy/dx = [3.4301 5.5175 6.3542 5.9365 4.8554 3.6837 2.7136]
z = [943.0802 1136.4406 1308.0878 1439.5741 1533.7588 1602.7919 1648.6146]
Thank you in advance!
  2 个评论
Star Strider
Star Strider 2020-9-12
‘... z is a function of x ...’
Is it possible for us to know what that functional relationship is?
Dylan Kline
Dylan Kline 2020-9-12
Well, the data is experimental so x and z are measured independently, though an equation could be fit if that is what you are referring to?

请先登录,再进行评论。

采纳的回答

Star Strider
Star Strider 2020-9-12
It turns out that your differential equation has an analytical solution:
syms a b n x y(x) z y0
Eq = diff(y) == a*(1-y)^n * exp(-b/z);
ys = dsolve(Eq, y(0) == y0)
ys = simplify(ys, 'Steps',500)
yfcn = matlabFunction(ys, 'Vars',{[a,n,b,y0],x,z})
In the 'Vars' name-value pair, this [a,n,b,y0] concatenates those variables (parameters here) into a single row vector (in that order) called ‘in1’. That is the parameter vector I later refer to as ‘p’.
Using ‘yfcn’ as calculated:
x = [0 0.0258 0.0515 0.0773 0.1031 0.1289 0.1546];
dydx = [3.4301 5.5175 6.3542 5.9365 4.8554 3.6837 2.7136];
z = [943.0802 1136.4406 1308.0878 1439.5741 1533.7588 1602.7919 1648.6146];
xz = [x(:) z(:)];
yfcn = @(in1,x,z)-((in1(:,2)-1.0).*((-in1(:,4)+1.0).^(-in1(:,2)+1.0)./(in1(:,2)-1.0)+in1(:,1).*x.*exp(-in1(:,3)./z))).^(-1.0./(in1(:,2)-1.0))+1.0;
p0 = rand(1,4);
dydxfcn = @(p,xz) gradient(yfcn(p,xz(:,1),xz(:,2))./gradient(xz(:,1)));
B = lsqcurvefit(dydxfcn, p0, xz, dydx(:));
If you want to integrate ‘dydx’, the code changes to:
yv = cumtrapz(x, dydx);
objfcn = @(p,xz) yfcn(p,xz(:,1),xz(:,2));
B = lsqcurvefit(objfcn, p0, xz, yv(:));
with no other changes.
Both of these run without error. I leave it to you to determine if it/they produce the correct results. Note that nonlinear parameter estimation problems can be very sensitive to the initial parameter estimates, and may end up in a local minimum rather than the desired global minimum, with unacceptable parameters. There are ways to deal with this (the ga function being my favourite because I’ve had the most experience with it over the years), however even those do not always find the best parameter set.
  2 个评论
Dylan Kline
Dylan Kline 2020-9-12
I must say, this is an elegant solution and I genuinely appreciate your time in helping me with this. Also, thank you for being such an active member of this community... I saw your name on almost all of the posts with similar questions and you have truly been a help to us at large.
Star Strider
Star Strider 2020-9-12
As always, my pleasure!
I very much appreciate your compliment!

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by