ode_covid19_v3.m
2.42 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
%ODE model to describe the spread and clinical progression of COVID-19
function dx=ode_covid19_v3(t,x,p,Tlock,s00,s11,s22)%,Tlock2,s02,s12)
be = p(1);
b0=p(2);
b1=p(3);
b2=p(4);
b3=p(5);
N=10^5;
be=be/(N);
b0=b0/(N);
b1=b1/(N);
b2=b2/(N);
b3=b3/(N);
FracSevere = p(6);
FracCritical = p(7);
FracMild = 1-FracSevere-FracCritical;
FracAsym = p(8);
IncubPeriod = p(9);
DurMildInf = p(10);
DurAsym = p(11);
DurHosp = p(12);
TimeICUDeath = p(13);
ProbDeath=p(14);
PresymPeriod=p(15)*IncubPeriod;
CFR=ProbDeath*FracCritical/100;
a1=1/PresymPeriod; %presymptomatic period of transmission
a0=1/(IncubPeriod-PresymPeriod); %true latent period
f=FracAsym;
g0=1/DurAsym;
g1=(1/DurMildInf)*FracMild;
p1=(1/DurMildInf)-g1;
p2=(1/DurHosp)*(FracCritical/(FracSevere+FracCritical));
g2=(1/DurHosp)-p2;
u=(1/TimeICUDeath)*(CFR/FracCritical);
g3=(1/TimeICUDeath)-u;
s0 = (t<Tlock(2)) + s00(1) * (t>=Tlock(2)) * (t<Tlock(3)) + s00(2)*(t>=Tlock(3)) * (t<Tlock(4)) + s00(3)*(t>=Tlock(4));
s1 = (t<Tlock(2)) + s11(1) * (t>=Tlock(2)) * (t<Tlock(3)) + s11(2)*(t>=Tlock(3))* (t<Tlock(4)) + s00(3)*(t>=Tlock(4));
s2 = (t<Tlock(1)) + s22 * (t>=Tlock(1));
%s3 = (t<Tlock(1)) + s33 * (t>=Tlock(1));
%susceptibles individuals (not infected)
S = x(1);
%exposed individuals, infected but no symptoms and no transmission
E0 = x(2);
%exposed individuals, infected but no symptoms and can transmit
E1 = x(3);
%infected individuals but asymptomatic infection
I0 = x(4);
%infected individuals, mild infection
I1 = x(5);
%infected individuals, severe infection
I2 = x(6);
%infected individuals, critical infection
I3 = x(7);
%recovered
R = x(8);
%dead
D = x(9);
%dS = - (be *(1/(1+s0^n))* E1 + b0*(1/(1+s0^n))*I0 + b1*(1/(1+s1^n))*I1 + b2*(1/(1+s2^n))*I2 + b3*(1/(1+s2^n))*I3)*S; % + k_e0s * E0;
%dE0 = (be *(1/(1+s0^n))* E1 + b0*(1/(1+s0^n))*I0 + b1*(1/(1+s1^n))*I1 + b2*(1/(1+s2^n))*I2 + b3*(1/(1+s2^n))*I3)*S -a0*E0; % - k_e0s * E0;
dS = - ((be * s0)* E1 + (b0*s0)*I0 + (b1*s1)*I1 + (b2*s2)*I2 + (b3*s2)*I3)*S; % + k_e0s * E0;
dE0 = ((be *s0)* E1 + (b0*s0)*I0 + (b1*s1)*I1 + (b2*s2)*I2 + (b3*s2)*I3)*S -a0*E0; % - k_e0s * E0;
dE1 = a0 * E0 - a1*E1; % -g4*E1;
dI0 = f * a1 * E1 - g0*I0;
dI1 = (1-f) * a1 * E1 - g1*I1 - p1*I1;
dI2 = p1 * I1 - g2*I2 - p2*I2; % -u2*I2;
dI3 = p2 * I2 - g3 *I3 - u * I3;
dR = g0*I0 + g1*I1 + g2*I2 + g3*I3; % +g4 * E1;
dD = u * I3; % + u2 * I2;
dx=[dS;dE0;dE1;dI0;dI1;dI2;dI3;dR;dD];