Fitting data with ODE equation and multiple inputs
显示 更早的评论
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
2020-9-12
‘... z is a function of x ...’
Is it possible for us to know what that functional relationship is?
Dylan Kline
2020-9-12
采纳的回答
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
2020-9-12
Star Strider
2020-9-12
As always, my pleasure!
I very much appreciate your compliment!
更多回答（0 个）
类别
在 帮助中心 和 File 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 
