## estimation of parameters using lsqnonlin

### vinutha paravastu (view profile)

on 25 Aug 2019
Latest activity Commented on by Alan Weiss

### Alan Weiss (view profile)

on 26 Aug 2019
hello all, i am in the process of estimating parameters for my model
i wrote the following sample codes..pls tell me whether they are right. i ran the script i got the result as
Initial point is a local minimum.
Optimization completed because the size of the gradient at the initial point
is less than 1e-4 times the default value of the function tolerance.
Pls let me know whether my approach is right or any changes to be made....pls help me i am new to matlab...
1. My model :
function dc = kinetics(p,t,c)
k1=p(1);
k2=p(2);
k3=p(3);
k4=p(4);
k5=p(5);
k6=p(6);
dc =zeros(5,1);
dcdt(1) = (k1+k2+k3)*c(1);
dcdt(2) = (k1*c(1))-(k4*c(2));
dcdt(3) = (k2*c(1))-(k5*c(3));
dcdt(4) = (k3*c(1))-(k6*c(4));
dcdt(5)= k6*c(4);
end
2. my objective function to estimate parameters :
function obj= interest (p,t,c1meas,c2meas,c3meas,c4meas,c5meas)
global c1meas c2meas c3meas c4meas c5meas
cdata =[c1meas ,c2meas,c3meas,c4meas,c5meas]
%simulate model
c=simulate(p);
obj = c-cdata;
end
3. Simulate function containing ode45 to solve differential equations
function c = simulate(p)
global t c1meas c2meas c3meas c4meas c5meas
c=zeros(length(t),5);
c0=[1,0,0,0,0];
for i=1:length(t)-1
ts=[t(i),t(i+1)];
[tsol,csol]= ode45(@(t,c)kinetics2(p,t,c),ts,c0);
c(i+1,1)=c0(1);
c(i+1,2)=c0(2);
c(i+1,3)=c0(3);
c(i+1,4)=c0(4);
c(i+1,5)=c0(5);
end
end
4. my final script for calling all the above functions and using lsqnonlin
global t c1meas c2meas c3meas c4meas c5meas
t= [0.88 ; 0.96; 0.98;1.04 ; 1.05];
c1meas=[0.211 ;0.066 ;0.17 ; 0.455; 0.088];
c2meas =[0.666 ;0.165 ;0.083 ;0.047 ;0.009];
c3meas=[0.302 ;0.093; 0.155; 0.341 ;0.094];
c4meas=[0.237;0.084; 0.177; 0.404 ;0.082];
c5meas=[0.686 ;0.16 ;0.072; 0.041; 0.008];
%parameters initial guess
k1=0;
k2=0;
k3=1;
k4=1;
k5=1;
k6=1;
p0=[k1,k2,k3,k4,k5,k6];
p= lsqnonlin(@(p) interest(p),p0)
% show final objective
disp(['Final SSE Objective: ' num2str(interest(p))]);
% optimized parameter values
k1=p(1);
k2=p(2);
k3=p(3);
k4=p(4);
k5=p(5);
k6=p(6);
disp(['k1: ' num2str(k1)])
disp(['k2: ' num2str(k2)])
disp(['k3: ' num2str(k3)])
disp(['k4: ' num2str(k4)])
disp(['k5: ' num2str(k5)])
disp(['k6: ' num2str(k6)])
% calculate model with updated parameters
ci = simulate(p0);
cp = simulate(p);

Alan Weiss

### Alan Weiss (view profile)

on 26 Aug 2019
I suspect that you are not passing t properly; it is declared as global in some places, but you do not declare it global in the kinetics function or interest function.
Generally, you should avoid using global variables. Learn to parametrize functions.
Alan Weiss
MATLAB mathematical toolbox documentation