Kinetic parameter estimation and initial time setting

Hi: how can I set the initial time on this code to be at 0 min at an initial concentration of 100 mg/ml. This is a case study of a degradation reaction for which at 0 min, 100 mg/ml of the reactant was all present, but at 20 min of the reaction, the reactant started degrading into other products. As such, the X-axis should start at 0 min and not 20 min. Please I need assistance on this. The code is as below.

function HtwoT_How

% 2016 12 03

% NOTES:

%

% 1. The ‘theta’ (parameter) argument has to be first in your

% ‘kinetics’ funciton,

% 2. You need to return ALL the values from ‘DifEq’ since you are fitting

% all the values

function C=kinetics(theta,t)

c0=[100;0;0;0];

[T,Cv]=ode45(@DifEq,t,c0);

%

function dC=DifEq(t,c)

dcdt=zeros(4,1);

dcdt(1)=-theta(1).*c(1);

dcdt(2)= theta(1).*c(1)-theta(2).*c(2)-theta(3).*c(2);

dcdt(3)= theta(2).*c(2)-theta(4).*c(3);

dcdt(4)= theta(4).*c(3);

dC=dcdt;

end

C=Cv;

end

t=[25

35

45

55

65];

c=[20.76 5.93 2.77 69.54

16.37 5.72 0.6 77.31

19.88 3.9 0.19 76.03

8.01 2 0.18 89.81

17.73 0.15 0.16 81.96];

theta0=[1;1;1;1];

[theta,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat]=lsqcurvefit(@kinetics,theta0,t,c);

fprintf(1,'\tRate Constants:\n')

for k1 = 1:length(theta)

fprintf(1, '\t\tTheta(%d) = %8.5f\n', k1, theta(k1))

end

tv = linspace(min(t), max(t));

Cfit = kinetics(theta, tv);

figure(1)

plot(t, c, 'p')

hold on

hlp = plot(tv, Cfit);

hold off

grid

xlabel('Time (min)')

ylabel('Concentration (mg/ml)')

legend(hlp, 'C1', 'C2', 'C3', 'C4', 'Location','N')

end

Answer 1

[Emmanuel Nzediegwu]

_

Answer 2

[Alex Sha]

Hi, see the results below:

Root of Mean Square Error (RMSE): 4.00513564326846
Sum of Squared Residual: 320.822230419589
Correlation Coef. (R): 0.788386977278981
R-Square: 0.621554025943088

Parameter          	Best Estimate
--------------------	-------------
theta1	0.0529114059211692
theta2	0.337616166932669
theta3	0.0180685434553373
theta4	52.0236365257367

[Emmanuel Nzediegwu]

Hi Alex, Thank you, please can you share the code you used as I have other data to attempt accordingly