fitting data to system of ode
2 views (last 30 days)
Show older comments
Chinwendu Madubueze
on 12 Oct 2021
Answered: Bjorn Gustavsson
on 13 Oct 2021
Hello everyone.
I am trying to fit data to a system of ODE using the following MatLab code.
clear all
realdata=[0 6; 1 12; 2 19; 3 25; 4 31; 5 38; 6 44; 7 60; 8 80;9 131;10 131;
11 259; 12 467; 13 688; 14 776; 15 1776; 16 1460; 17 1739; 18 1984; 19 2101; 20 2590;
21 2827; 22 3233; 23 3892; 24 3697; 25 3151; 26 3387; 27 2653; 28 2984; 29 2473; 30 2022;
31 1820; 32 1998; 33 1506; 34 1278; 35 2051; 36 1772; 37 1891; 38 399; 39 894; 40 397;
41 650; 42 415; 43 518; 44 412; 45 439; 46 441; 47 435; 48 579; 49 206; 50 130;
51 120; 52 143; 53 146; 54 102; 55 46; 56 45; 57 20; 58 31; 59 26; 60 11;
61 18; 62 27; 63 29; 64 39; 65 39];
aux=size(realdata);
%initial and end points in days
t0=0;
tend=realdata(aux(1),1);
time=t0:tend;
totalill=realdata(:,2);
%numerical constants of the model
LambdH=100;beta1=0.221;beta4=0.052;beta3=0.0206;beta2=0.0296;beta5=0.0637;omega=0.00578;d1=0.003;
muR=0.003;sigma=0.0528;muv=0.0485;gamma=0.0024;alpha1=0.4;alpha2=0.003;LambdR=0.172;
muH=0.005;
N=20000;
initialvalue=realdata(1,2);
e0=5;
i0=initialvalue-e0;
s0=N-i0;
r0=0;
sr0=2000;
Ir0=10;
v0=4;
%stepsize of the numerical method
stepsize=0.001;
deltaH = 0.005;
%differential system with the substitutions in system
%S=X(1) %E=X(2) %I=X(3) %P=X(4) %A=X(5) %H=X(6) %R=X(7) %D=X(8) dead people
system1=@(t,X)[
-beta1.*X(3).*X(1)-beta2.*X(6).*X(1)-beta3.*X(7).*X(1) - muH.*X(1)+LambdH + omega.*X(4);
beta1.*X(3).X(1)+ beta2.*X(6).*X(1) + beta3.*X(7).*X(1)-(sigma+muH).*X(2);
sigma.*X(2)-(muH+gamma).*X(3)-deltaH.*X(3);
gamma.*X(3)-(omega+muH).*X(4);
LambdR-beta4.*X(5).*X(6) - beta5.*X(5).*X(7)-(d1+muR).*X(5);
beta4.*X(5).*X(6)+ beta5.*X(5).*X(7)-(muR+d1).*X(6);
alpha1.*X(6)+alpha2.*X(3)-muv.*X(7);
];
[ts1,ys1]=fde12(1,system1,t0,tend,[s0;e0;i0;r0;sr0;Ir0;v0],stepsize);
% plot I+P+H
figure
% hold on
plot(time,totalill,'green-','LineWidth',2.5)
xlabel({'Time','(in days)'})
ylabel('Confirmed cases per day')
plot(ts1(1,1:end),ys1(3,:),'black','linewidth',2);
I keep getting the following error which i do not have an ideal on how to solve the problem. Please I need help. I am not perfect with MatLab. I am still learning. Thank you.
??? Improper index matrix reference.
Error in ==>
@(t,X)[-beta1.*X(3).*X(1)-beta2.*X(6).*X(1)-beta3.*X(7).*X(1)-muH.*X(1)+LambdH+omega.*X(4);beta1.*X(3).X(1)+beta2.*X(6).*X(1)+beta3.*X(7).*X(1)-(sigma+muH).*X(2);sigma.*X(2)-(muH+gamma).*X(3)-deltaH.*X(3);gamma.*X(3)-(omega+muH).*X(4);LambdR-beta4.*X(5).*X(6)-beta5.*X(5).*X(7)-(d1+muR).*X(5);beta4.*X(5).*X(6)+beta5.*X(5).*X(7)-(muR+d1).*X(6);alpha1.*X(6)+alpha2.*X(3)-muv.*X(7);]
Error in ==> fde12>f_vectorfield at 335
f = feval(Probl.fdefun,t,y) ;
Error in ==> fde12 at 121
f_temp = f_vectorfield(t0,y0(:,1),Probl) ;
Error in ==> datalassa at 52
[ts1,ys1]=fde12(1,system1,t0,tend,[s0;e0;i0;r0;sr0;Ir0;v0],stepsize)
0 Comments
Accepted Answer
Bjorn Gustavsson
on 13 Oct 2021
You will find the answer to this problem in a link to yet another answer (by Star Strider) here: fit-a-parameterised-ode-solution-to-experimental-data.
HTH
0 Comments
More Answers (0)
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!