fitting data to system of ode

조회 수: 1 (최근 30일)
Chinwendu Madubueze
Chinwendu Madubueze 2021년 10월 12일
답변: Bjorn Gustavsson 2021년 10월 13일
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)

채택된 답변

Bjorn Gustavsson
Bjorn Gustavsson 2021년 10월 13일
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개)

카테고리

Help CenterFile Exchange에서 Programming에 대해 자세히 알아보기

태그

Community Treasure Hunt

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

Start Hunting!

Translated by