fitting the solution of odes to data using fmincon with nonlinear restraints

1 view (last 30 days)
I have some data (abswl2) that changes with time. I can fit the data to a set of ODEs that represent a two-step reaction with both steps reversible (four rate constants). The data is readily fit using lsqcurvefit; however because of the nature of the system of ODEs, a unique solution does not exist. To get around this, I would like to perform the fit using fmincon enforcing a nonlinear restraint to guide the minimization to physically meaningful results. I have tried several variations of the script below with the error that there are insufficient input arguments. I would greatly appreciate any advice from the community.
The most recent version script is as follows:
%%Extracts Data
data=dlmread('02pAromASD.o3a');
wl1=416; wl2=392;
time=data(1,2:end-5); time=time';
wl=data(2:end,1);
idx1=find(wl>(wl1-0.5)&wl<(wl1+0.5),1,'first');
idx2=find(wl>(wl2-0.5)&wl<(wl2+0.5),1,'first');
abs=data(2:end-5,2:end-5);
abswl1=abs(idx1,:);abswl1=abswl1';
abswl2=abs(idx2,:);abswl2=abswl2';abswl2=abswl2-min(abswl2);
abs=[abswl1,abswl2];
%%Initial constants, upper & lower bounds
K0=rand(1,4);
lb=zeros(1,4);
ub=Inf*ones(1,4);
%%Unconstrained Minimization with Constraints
[K1,fval,exitflag,output]=fmincon(@ssq,K0,[],[],[],[],lb,ub,@dissoc);
%%Objective Functions
function err=ssq(K,time,abswl2)
err=sum((twostepreverse(K,time)-abswl2).^2)
end
%%Solution to differential equations
function C=twostepreverse(K,time)
Y0=[1;0;0];
[tSol,YSol]=ode113(@diffeq,time,Y0);
function dYdt = diffeq(time,Y)
a=Y(1);
b=Y(2);
c=Y(3);
k1=K(1);
km1=K(2);
k2=K(3);
km2=K(4);
L=4;
dadt=-k1*L*a+km1*b;
dbdt=k1*L*a-(km1+k2)*b+km2*c;
dcdt=k2*b-km2*c;
dYdt=[dadt;dbdt;dcdt];
end
C=YSol(:,3);
end
%%Nonlinear Constraints
function [c,ceq]=dissoc(K)
k1=K(1);
km1=K(2);
k2=K(3);
km2=K(4);
%Nonlinear inequality constraints
c=[];
%Nonlinear equality constraints
ceq=((km2*km1)/(k1*k2))-1;
end
After running the script, this error is returned:
Not enough input arguments.
Error in test2>ssq (line 25)
err=sum((twostepreverse(K,time)-abswl2).^2)
Error in fmincon (line 536)
initVals.f = feval(funfcn{3},X,varargin{:});
Error in test2 (line 21)
[K1,fval,exitflag,output]=fmincon(@ssq,K0,[],[],[],[],lb,ub,@dissoc);
Caused by:
Failure in initial objective function evaluation. FMINCON cannot continue.

Accepted Answer

Walter Roberson
Walter Roberson on 26 Feb 2018
You are using
[K1,fval,exitflag,output]=fmincon(@ssq,K0,[],[],[],[],lb,ub,@dissoc);
That is going to call ssq with a single parameter that is a vector. Your ssq function expects three inputs and will fail if it attempts to access the second or third inputs as those are missing.
  1 Comment
John Hackett
John Hackett on 27 Feb 2018
Edited: John Hackett on 27 Feb 2018
That works! Thanks. I ended up defining the least squares as an anonymous function:
[K1,fval,exitflag,output]=fmincon(@(K)sum((onesteprev(K,time)-abswl2).^2),K0,[],[],[],[],lb,ub,@dissoc);

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!