Journal Pre-proofs

A mathematical COVID-19 model considering asymptomatic and symptomat-

ic classes with waning immunity

Nursanti Anggriani, Meksianis Z. Ndii, Rika Amelia, Wahyu Suryaningrat,

Mochammad Andhika Aji Pratama

PII: S1110-0168(21)00351-3
DOI: https://doi.org/10.1016/j.aej.2021.04.104
Reference: AEJ 2346

To appear in: Alexandria Engineering Journal

Received Date: 1 February 2021

Revised Date: 21 April 2021
Accepted Date: 27 April 2021

Please cite this article as: N. Anggriani, M.Z. Ndii, R. Amelia, W. Suryaningrat, M.A. Aji Pratama, A
mathematical COVID-19 model considering asymptomatic and symptomatic classes with waning immunity,
Alexandria Engineering Journal (2021), doi: https://doi.org/10.1016/j.aej.2021.04.104

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover
page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version
will undergo additional copyediting, typesetting and review before it is published in its final form, but we are
providing this version to give early visibility of the article. Please note that, during the production process, errors
may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2021 THE AUTHORS. Published by Elsevier BV on behalf of Faculty of Engineering, Alexandria University.
Manuscript Click here to view linked References

1
2
3
4
5
6
7
8
9 A mathematical COVID-19 model considering
10
11
asymptomatic and symptomatic classes with waning
12 immunity
13
14
15 Nursanti Anggriania , Meksianis Z. Ndiib , Rika Ameliaa , Wahyu
16 Suryaningrata , Mochammad Andhika Aji Pratamaa
17 a Department of Mathematics, Universitas Padjadjaran, Jln. Raya Bandung-Sumedang Km.
18 21 Jatinangor, Kab. Sumedang 45363 Jawa Barat, Indonesia
19 b Department of Mathematics, Faculty of Sciences and Engineering, The University of Nusa
20 Cendana, Kupang-NTT, Indonesia
21
22
23
24
25 Abstract
26
27 The spread of COVID-19 to more than 200 countries has shocked the public.
28
29 Therefore, understanding the dynamics of transmission is very important. In
30 this paper, the COVID-19 mathematical model has been formulated, analyzed,
31
32 and validated using incident data from West Java Province, Indonesia. The
33
model made considers the asymptomatic and symptomatic compartments and
34
35 decreased immunity. The model is formulated in the form of a system of dif-
36
37 ferential equations, where the population is divided into seven compartments,
38 namely Susceptible Population (S0 ), Exposed Population (E), Asymptomatic
39
40 Infection Population (IA ), Symptomatic Infection Population (YS ), Recovered
41 Population (Z), Susceptible Populations previously infected (Z0 ), and Quar-
42
43 antine population (Q). The results show that there has been an outbreak of
44
45 COVID-19 in West Java Province, Indonesia. This can be seen from the basic
46 reproduction number of this model, which is 3.180126127 (R0 > 1). Also, the
47
48 numerical simulation results show that waning immunity can increase the oc-
49 currence of outbreaks; and a period of isolation can slow down the process of
50
51 spreading COVID-19. So if a strict social distancing policy is enforced like a
52
53
∗ Corresponding author at Department of Mathematics, Universitas Padjadjaran, Jln. Raya
54
Bandung-Sumedang Km. 21 Jatinangor, Kab. Sumedang 45363 Jawa Barat, Indonesia
55
Email address: nursanti.anggriani@unpad.ac.id; nursanti.anggriani@gmail.com
56 (Nursanti Anggriani)
57
58
59 Preprint submitted to Alexandria Engineering Journal April 21, 2021
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 quarantine, the outbreak will lessen.
10
11 Key words: COVID-19, Basic Reproduction Ratio, waning immunity,
12 asymptomatic, previous infection, parameter estimation
13
14
15
1. Introduction
16
17
18 The spread of COVID-19 has shocked society and currently has transmitted
19
to more than 200 countries [1]. As of 04 April 2021, there are 130,998,190 con-
20
21 firmed cases, 2,853,280 death, and 105,447,782 recovered individuals [2]. It has
22
23 caused severe economic and social loss. The disease has been transmitted from
24 human to human via droplets [3]. Infected individuals may show symptomps
25
26 such as fever, cough, sore throat, rhinorrhea, myalgia or fatigue, phlegm, and
27 headache [3–5] with the body temperature of 39◦ C or above [5]. Individuals who
28
29 are infected by COVID-19 can show symptoms (symptomatic) or cannot show
30
symptoms (asymptomatic) but both types of individuals can transmit disease
31
32 [3]. The incubation period has been estimated between two two fourteen days
33
34 [6].
35 Research showed that there is possibility for infected individuals to be re-
36
37 infected by COVID-19. Currently it has been found that several recovered
38
individuals have been re-infected by COVID-19 and this can cause death from
39
40 fatal heart failure [7]. Of the 111 recovered patients, 5% of China and
41
42 10% of South Korea tested positive for COVID-19 [8]. This situation
43 contradicts the fact that after a person catches the virus and then recovers,
44
45 the individual will forman antibody that prevents the same virus from attack-
46 ing twice. Research showed that reinfected individuals have experienced viral
47
48 replication but did not neutralize antibodies, which implies that it is unlikely
49
that long-term protective immunity will occur in people with COVID-19 after
50
51 the first infection [9]. The virus’s immune response can be reduced within four
52
53 months to one year after infection [10]. The genetic basis of the innate immune
54 response affects the severity of COVID-19, it can also lead to more severe rein-
55
56 fection depending on antibodies generated against the bound virus but cannot
57
58
59 2
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 neutralize the same strain [10]. The reinfection COVID-19 case has a more se-
10
vere impact [11]. The reinfection occurs due to the decrease in the individual’s
11
12 immunity [12]. Understanding the effects of waning immunity is important.
13
14 Mathematical models can be used to understand the complex phenomena
15 such as population dynamics problem [13–16] and disease transmission dynam-
16
17 ics [17–21]. A compartment-based epidemic model in the form of system of
18
(fractional or integer) differential equations has been formulated to understand
19
20 disease transmission dynamics, where the human population is divided into dif-
21
22 ferent stages according to their status to the diseases [22–36]. A mathematical
23 SEIR model is mostly used as a basis for the model’s development for COVID-19
24
25 transmission [37, 38]. The SEIR model has been extended to include quarantine
26 compartment [39], to include symptomatic and asymptomatic classes [40]. The
27
28 models are mostly used to investigate the disease transmission dynamics in sev-
29
eral countries or provinces such as Indonesia [41], Hubei has been researched [42],
30
31 Pakistan [43]. In this paper, a modified SEIR model considering symptomatic
32
33 and asymptomatic cases from [44] has been formulated. The work focuses on
34 studying the effects of waning immunity. The model is validated against data of
35
36 COVID-19 incidence from West Java Province. The basic reproduction number
37
is calculated, and a global sensitivity analysis is performed. The model is then
38
39 used to determine he effects of waning immunity or reduced immunity to an
40
41 increase in the number of infected individuals.
42 The remainder of this paper is organized as follows. In Section 2, the con-
43
44 struction of the SEIR compartmental model. Next, the model’s mathematical
45 properties, such as the equilibrium points, Basic Reproduction Number, and the
46
47 existence of backward bifurcation, are detailed in Section 3. In Section 4, we
48
explain the real-world problem using the incidence data of West Java Province,
49
50 Indonesia. A discussion on the Basic Reproduction Number and the sensitiv-
51
52 ity analysis results are provided in Section 5. Finally, some conclusions are
53 presented in Section 6.
54
55
56
57
58
59 3
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 2. Model Formulation
10
11 We developed model of transmission of COVID-19 by considering asymp-
12
13 tomatic and symptomatic compartments and decreased immunity. The total
14
population is divided into Susceptible population (S0 ), Exposed population (E),
15
16 Asymptomatic infected population (IA ), Symptomatic infected population (YS ),
17
18 Recovered population (Z), Susceptible that previously infected (Z0 ), Quaran-
19 tine population (Q). The total number of population at time t is given by:
20
21
22 N (t) = S0 (t) + E(t) + IA (t) + YS (t) + Z(t) + Z0 (t) + Q(t).
23
24
The assumption used in the formulation of a mathematical model for the
25
26 spread of the COVID-19 disease is that individuals with symptoms will undergo
27
28 hospitalization or quarantine. Deaths experienced by latent, symptomatic,
29 asymptomatic, and quarantine individuals are caused by disease [45]. This
30
31 means that the death of the three individuals is a combination of natural death
32
and death due to disease. We assume that the µ1 parameter contained in com-
33
34 partments E, Ia , Ys is a death caused by COVID-19 plus a natural death factor.
35
36 People who have decreased immunity can catch COVID-19 again with a high
37 severity [11]. Hence, the second person infected will develop symptoms and be
38
39 hospitalized. The model is represented by the diagrams shown in Figure 1, with
40 the description of the parameters given in Table 1.
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59 4
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Table 1: Parameters Description
10
11
12 Parameters Descriptions
13
14 Λ The recruitment rate of susceptible population
15 β1 The probability of transmission from asymptomatic infected
16
17 people
18 β2 The probability rate of transmission from symptomatic in-
19
20 fected people
21
22 µ The natural mortality rate
23 µ1 Natural death rate plus COVID-19 death rate
24
25 α The probability of exposed people become infected
26 p The proportion of exposed people become infected
27
28 κ The rate of asymptomatic infected people become infected
29
symptomatic
30
31 q The rate of quarantine
32
33 γ1 The natural recovery rate of infected asymptomatic people
34 γ2 The natural recovery rate of infected symptomatic people
35
36 δ The recovery rate of quarantine people
37 ξ The probability rate of recovered people become susceptible
38
39 (waning immune)
40
41
42 So, based on the interaction diagram above, the COVID-19 spread mathe-
43
44 matics model constructed as follows:
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59 5
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
Figure 1: Interaction Diagram of Populations
25
26
27
28
29
30 dS0
31 = Λ − (β1 IA S0 + β2 YS S0 ) − µS0 (1)
32 dt
33 dE
= (β1 IA S0 + β2 YS S0 ) − αE − µ1 E (2)
34 dt
35 dIA
= pαE − κIA − γ1 IA − µ1 IA (3)
36 dt
37 dYS
= (1 − p)αE − qYS − γ2 YS + β1 IA Z0 + β2 Ys Z0
38 dt
39
+ κIA − µ1 YS (4)
40
41 dZ
42 = γ1 IA + γ2 YS + δQ − ξZ − µZ (5)
dt
43 dZ0
44 = ξZ − β1 IA Z0 − β2 Ys Z0 − µZ0 (6)
dt
45 dQ
46 = qYS − δQ − µ1 Q, (7)
dt
47
48 with S0 (0) ≥ 0, E(0) ≥ 0, IA (0) ≥ 0, YS (0) ≥ 0, Z(0) ≥ 0, Z0 (0) ≥ 0, Q(0) ≥
49
50 0 as the initial conditions.
51
52
53
54
55
56
57
58
59 6
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 3. Mathematical Analysis
10
11 Lemma 3.1. If the initial values S0 (0) > 0, E(0) > 0, IA (0) > 0, YS (0) >
12
13 0, Z(0) > 0, Z0 (0) > 0, and Q(0) > 0, the solution of
14
15 S0 (t), E(t), IA (t), YS (t), Z(t), Z0 (t), Q(t),
16
17 of system (1-7) are positif for all t > 0.
18
19 Proof. Assume that
20
21 X (t) = min{S0 (t), E(t), IA (t), YS (t), Z(t), Z0 (t), Q(t).}, ∀t > 0.
22
23 Clearly, X (0) > 0.
24
25
26 Assuming that there exist a t1 > 0 such that
27
28 X (t1 ) = 0 and X (t) > 0, for all t ∈ [0, t1 ),
29 If X (t1 ) = S0 (t1 ), then E(t) ≥ 0, IA (t) ≥ 0, YS (t) ≥ 0, Z(t) ≥ 0, Z0 (t) ≥ 0 for
30
31 all t ∈ [0, t1 ].
32
33
34 From the equation of model (1), we can obtain
35
36 dS0
≥ −β1 IA S0 − β2 YS S0 − µS0 , t ∈ [0, t1 ].
37 dt
38
39
40 Thus, we have
41  Z t1 h i 
42 S0 (t) ≥ S0 (0) exp − β1 IA S0 + β2 YS S0 + µS0 dt ,
43 0
44
45
46 which will be positive since exponential functions and initial solutions S0 (0)are
47 non-negative. Thus, S0 (t) > 0 for all t ≥ 0.
48
49
50
Similarly, we can also prove that
51
52
53
54 E(t) > 0, IA (t) > 0, YS (t) > 0, Z(t) > 0,
55
56 Z0 (t) > 0, Q(t) > 0.
57
58
59 7
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11 Lemma 3.2. All solution of system (1-7) are bounded for all t ∈ [0, t0 ]
12
13
Proof. Since N (t) = S0 (t) + E(t) + IA (t) + YS (t) + Z(t) + Z0 (t) + Q(t).
14
15 We get:
16 dN
17 = Λ − µ(S0 + Z + Z0 ) − µ1 (E + IA + YS + Q).
dt
18
19 Assume that µ = µ1 , to simplify the analysis process.
20 Then:
21 dN
22 = Λ − µN.
23 dt
24 Thus we have
25 Λ
26 0 ≤ lim sup N (t) ≤ ,
t→∞ µ
27
28 so all solutions of system (1-7) are ultimately bounded for all t ∈ [0, t0 ].
29
30 3.1. Non-endemic Equilibrium point
31
32 The non-endemic equilibrium point of the COVID-19 disease model is ob-
33
34 tained by setting IA = 0, E = 0, YS = 0, and substituting it into (1-7) to obtain:
35
36
37
38 P0 = (S0 0 , E 0 , IA 0 , YS 0 , Z 0 , Z0 0 , Q0 )
39 Λ 
40 = , 0, 0, 0, 0, 0, 0 , (8)
µ
41
42
43 3.2. Stability of Non-endemic Equilibrium Point
44 Theorem 3.3. The non-endemic equilibrium point of system (1-7) is locally
45
46 asymptotically stable whenever it exists.
47
48 Proof. By following Diekmann (2000) [46] substituting P0 from 8 into the Ja-
49
50 cobian matrix for the non-endemic equilibrium point is obtained:
51
52
53
54
55
56
57
58
59 8
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10  
−βΛ −βΛ
11 −µ 0 µ µ 0 0 0
12
 
βΛ βΛ
−α − µ
 
13  0 µ µ 0 0 0 
 
14 
 0 pα −κ − γ1 − µ1 0 0 0 0


15  
J(P0 ) =  0 (1 − p)α −q − γ2 − µ1 .
 
16 κ 0 0 0
 
17 
 0 0 γ1 γ2 −ξ − µ 0 δ


18  
−µ
 
19  0 0 0 0 ξ 0 
 
20
0 0 0 q 0 0 −δ − µ1
21
22
23 The characteristics of the polynomial is
24
1
25 P(λ) = ((λ + µ)P1 (λ)) = 0, (9)
26 µ
27 P1 (λ) = a0 λ6 + a1 λ5 + a2 λ4 + a3 λ3 + a4 λ2 + a5 λ + a6 .
28
29
30 From the polynomial (P(λ)) we get λ1 = −µ and for λi with i = 2, 3, . . . , 7
31 will be negative if aj > 0 where j = 0, 1, 2, . . . , 6, R0 < 1, a1 a2 > a0 a3 ,
32
33 a1 (a2 a3 + a0 a5 ) > a21 a4 + a0 a23 , and a1 a2 a4 > a0 (a1 a6 + a2 a5 ). Since the
34
coefficients in the characteristic equation P1 (λ) are complex, we proceed to
35
36 analyze the coefficient values numerically with β1 = β2 . The results of the
37
38 numerical analysis obtained (see Appendix A.), show that for λi with i =
39 2, 3, . . . , 7 negative. Because λj with j = 1, 2, 3, . . . , 7 is negative, it can be
40
41 concluded that the non-endemic equilibrium point of the system (1-7) is locally
42 stable, so Theorem 3.3 is proven.
43
44
45 3.3. Basic reproduction ratio
46
47 The Basic Reproduction Ratio (R0 ) is an important number in epidemiology,
48
which is defined as the number of secondary infections caused by one primary
49
50 infection in a population. We use the next-generation method to determine
51
52 R0 , the value of R0 an be obtained by finding the dominant eigenvalue F V −1 .
53 Where F and V are Jacobian matrices of f (newly infected matrices) and v
54
55 (exiting matrices) that are evaluated at the disease-free equilibrium point (P0 )
56
57
58
59 9
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 from 8. From the models (1-7) are obtained:
10    
11 0 Λβµ
1 Λβ2
µ −α − µ 0 0
12    
F = 0 , V =  αp −κ − γ1 − µ
   
13 0 0 0 
   
14 0 0 0 (1 − p)α κ −q − γ2 − µ1
15
16
and
17
18
 
− α (((p−1)βµ2 −pβ 1 )µ1 +(γ1 p−κ−γ1 )β2 −pβ1 (γ2 +q))Λ
(µ1 +α)(µ1 +κ+γ1 )(γ2 +q+µ1 )
Λ (β1 (γ2 +q+µ1 )+β2 κ)
µ (γ2 +q+µ1 )(µ1 +κ+γ1 )
Λ β2
µ (γ2 +q+µ1 ) 
19 
20 F V −1 = .
 
0 0 0
21  
22 0 0 0
23
24 The eigenvalues of (F V −1 ) are:
25
26
α (((p − 1) β2 − pβ1 ) µ1 + (γ1 p − κ − γ1 ) β2 − pβ1 (γ2 + q)) Λ
27 λ1 = − , λ2,3 = 0.
28 µ (µ1 + α) (µ1 + κ + γ1 ) (γ2 + q + µ1 )
29
30
31 Following the method described by Castilo Cavez et al. (2002) [47] the basic
32 reproduction number in the COVID-19 is:
33
34 α (((p − 1) β2 − pβ1 ) µ1 + (γ1 p − κ − γ1 ) β2 − pβ1 (γ2 + q)) Λ
35 R0 = ξ(F V −1 ) = − .
µ (µ1 + α) (µ1 + κ + γ1 ) (γ2 + q + µ1 )
36
37 Under certain conditions where the probability of transmission from infected
38
39 people same as from asymptomatic infected people hold β1 = β2 and the natural
40 recovery rate of infected people asymptomatic and symptomatic γ1 = γ2 = γ.
41
42 It obtained the reproduction number for this condition symbolized by R0β .
43
44 Where:
45
46 Λαβ(pq + µ1 + γ + κ)
47 R0β = .
µ(µ1 + γ + κ)(µ1 + α)(µ1 + γ + q)
48
49 3.4. Endemic Equilibrium Points
50
51 Theorem 3.4. An endemic equilibrium point of system
52
53 P1 = (S0 ∗ , E ∗ , IA ∗ , YS ∗ , Z ∗ , Z0 ∗ , Q∗ ) will exist if G > 0 and H > 0 or G < 0
54 and H < 0.
55
56
57
58
59 10
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Proof. The endemic point of this disease is endemic in certain areas for a certain
10
period, which releases the COVID-19 in the population. It is indicated by the
11
12 presence of compartments exposed to virus transmission E ∗ , Ia∗ , YS∗ at steady
13
14 state. By calculating model (1, 5-7) and setting the right hand side zero we
15 obtained:
16
17 Λ
18 S0∗ = ∗ ,
β(IA + YS∗ ) + µ
19
20 ∗ αEp
IA = ,
21 γ + κ + µ1
22 (1 − p)αE∗ + βI∗A Z∗0
YS∗ = ,
23 q + γ − βZ∗0
24 γ(I∗ + Y∗S ) + δQ∗
25 Z∗ = A ,
26 (µ + ξ)
27 ξZ∗
Z0∗ = ,
28 (β(IA + Y∗S ) + µ
∗

29 qY∗S
30 Q∗ = .
δ + µ1
31
32
33
34
35 By substituting S0∗ , IA
∗
, YS∗ , Z0∗ , Z0∗ , Q∗ to equation (2.2) and set the right
36 hand side equal to zero, obtained:
37
38
39
40 A2 E 2 + A1 E + A0 = 0, (10)
41
42
43 which this polynomial have to roots E = 0 or E = E ∗ which can be written
44 by
45 G
46 E∗ = ,
H
47
48 where G and H written on Appendix B. .
49
50 Because the denomerator of H always positif, the steady state E ∗ will exsist
51 if R0β > 1 and pq > 0, see attachment for the proof. The system of equations
52
53 (1-7) will have an endemic equilibrium point if G > 0 and H > 0 or G < 0 and
54 H < 0. This condition indicates that the system of equations (1-7) has a unique
55
56 endemic equilibrium point.
57
58
59 11
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 3.5. Stability of Endemic Equilibrium Point
10
11 Theorem 3.5. The endemic equilibrium point of the system (P1 ) is locally
12 asymptotically stable whenever it exists.
13
14
15 Proof. Following method on the proof of Theorem 2, Based on the method of
16 proof of Theorem 2, by substituting P1 is obtained characteistic polynomial
17
18
19 Q(λ) =a0 λ7 + a1 λ6 + a2 λ5 + a3 λ4 + a4 λ3 +
20
21 a5 λ2 + a6 λ + a7 . (11)
22
23 From the polynomial (Q(λ)) we get λi with i = 1, 2, 3, . . . , 7 will be negative
24
25 if aj > 0 where j = 0, 1, 2, . . . , 7, a1 a2 > a0 a3 , a1 (a2 a3 + a0 a5 ) > a21 a4 + a0 a23 ,
26 a1 a2 (a3 a4 + a0 a7 ) > a0 a3 (a1 a6 + a2 a5 ), and a2 a5 > a0 a7 . Since the coefficients
27
28 in the characteristic equation Q(λ) are complex, we proceed to analyze the
29 coefficient values numerically. The results of the numerical analysis obtained
30
31 can be seen in the Appendix C. It satisfies the Routh-Hurwitz’s criteria so
32
that the endemic point is locally asymptotically stable whenever it exists. These
33
34 results will remain consistent using the parameter values in Table 2.
35
36
37 3.6. Global Stability of the equilibria
38
39 Theorem 3.6. The non-endemic equilibrium point (P0 ) is globally asymptoti-
40 cally stable if
41 β1 Λ β2 Λ
42 < µ1 and < µ1 .
µ µ
43
44 Proof. Refer to global proving by Tewa et al. (2009)[48], let
45  
46 P0 = (S0 0 , E 0 , IA 0 , YS 0 , Z 0 , Z0 0 , Q0 ) = Λµ , 0, 0, 0, 0, 0, 0 is the non-endemic
47 equilibrium point of system (1-7).
48
49
50 Define the Lyapunov function
51
52  ln S0 
53 V (t) = S0 − S0∗ − S0∗ ∗ + E + IA + YS + Z + Z0 + Q.
S0
54
55 Differentiating with respect to time yields
56
57
58
59 12
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10 dS0
11 dV dE dIA dYS Z Z0 Q
=(S0 − S0∗ ) dt∗ + + + + + +
12 dt S0 dt dt dt dt dt dt
13 dV (S0 − S0∗ )(Λ − µS0 )
14 = − (S0 − S0∗ )(β1 IA + β2 YS ) + β1 IA S0 + β2 YS S0
dt S0∗
15
16 − αE − µ1 E + pαE − κIA − γIA − µ1 IA + (1 − p)αE − qYS − γYS
17
18 + β1 IA Z0 + β2 YS Z0 + κIA − µ1 YS + γIA + γYS + δQ − ξZ − µZ
19
20 + ξZ − β1 IA Z0 − β2 YS Z0 − µZ0 + qYS − δQ − µ1 Q
21 dV (S0 − S0∗ )(Λ − µS0 )
22 = − (S0 − S0∗ )(β1 IA + β2 YS ) + β1 IA S0 + β2 YS S0
dt S0
23
24 − µ1 (E + IA + YS + Q) − µ(Z + Z0 )
25
dV 1 1
26 =Λ(S0 − S0∗ )( − ∗ ) + (β1 S0∗ − µ1 )IA + (β2 S0∗ − µ1 )YS − µ1 (E + Q)
27 dt S0 S0
28 − µ(Z + Z0 )
29
30 dV 1 1 β1 Λ β2 Λ
=Λ(S0 − S0∗ )( − ∗) + ( − µ1 )IA + ( − µ1 )YS − µ1 (E + Q)
31 dt S0 S0 µ µ
32
33 − µ(Z + Z0 ).
34
35
36
dV
37 The value of dt will be negative if
38
39 β1 Λ β2 Λ
< µ1 and < µ1 .
40 µ µ
41
42 By following LaSalle’s extension on Lyapunov’s method [49], disease-free
43 equilibrium P0 is globally asymptotically stable.
44
45 This concludes the proof.
46
47
48 4. Sensitivity analysis
49
50 This section presents a global sensitivity analysis of the model. We use the
51
52 combination of Latin Hypercube Sampling (LHS) and Partial Rank Correlation
53 Coefficient (PRCC) to determine the most influential parameters of the model
54
55 [50]. LHS is stratified sampling without replacement. The parameter distribu-
56
tion is divided into equation probability intervals and then is sampled. Each
57
58
59 13
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
1
11
12
13 0.8
14
15 0.6
16
17 0.4
18 1
19
20 0.2 2
PRCC

21 q
22 0 p
23 k
24 -0.2
1
25
26 -0.4
27
28
29 -0.6
30
31 -0.8
32
33 -1
34 0 100 200 300 400 500 600 700
35
Time (days)
36
37
38 Figure 2: PRCC over time when we measure against the increasing number of infected indi-
39 viduals.
40
41
42 parameter interval is sampled once and the entire range of each parameter is
43
explored. A matrix is then generated which consists of N rows for the number
44
45 of samples and k columns for the number of varied parameters. The model
46
47 solution is then generated using the combination of parameters (each row). The
48 outcome of interest is the increasing number of infected individuals. The result
49
50 of sensitivity analysis is given in Figure 2.
51 It showed that the waning immunity (ξ) is one of the influential parameters.
52
53 When the value of waning immunity increases, the number of infected individ-
54
55 uals also increases. This means that waning immunity would contribute to the
56 increasing number of infected individuals. Therefore, an analysis of effects of
57
58
59 14
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 waning immunity is of importance. The parameters k, p are also influential and
10
has positive relationship. This means that the rate of asymptomatic become
11
12 infected and the proportion of exposed individuals become infected contributes
13
14 to an increasing number of infected individuals. When these parameter values
15 increases, the number of infected individuals decreases.
16
17
18
5. A Case Study
19
20
21 In this section, we estimated the parameters β1 ,β2 , and γ against data
22 of West Java, Indonesia. The data are obtained from the website https:
23
24 //pikobar.jabarprov.go.id/table-case/ . We estimate the parameter val-
25
26 ues by minimizing the sum of squared error. The parameters b and a are
27 estimated against the data for the first 30 days. It is sufficient since the aim is
28
29 to obtain the general insights of the values of parameters β1 , β2 and q in the
30 early period of the outbreak. The other parameter values are obtained from
31
32 literature and are given in Table 2.
33
34
35
36
37
38
39
40
41
42
43
44
45
46 (a) (b)
47
48 Figure 3: Fitting Parameter from Confirmed Cases (a) and Cumulative Death (b)
49
50
51 The lsqnonlin built-in function in MATLAB is used for the parameter esti-
52 mation.
53
54
55
56
57
58
59 15
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 We minimize the sum of squared error as
10
11 n
X 2
12 SE = (Qt − gt (x)) . (12)
13 t=1
14
where Qt is the number of active cases of Q up to day t, respectively, while
15
16 gt (x) is the number of active cases for Q up to day t from the model’s solution,
17
18 respectively. The transmission rate, β0 and β1 , the quarantine rate q are then
19 estimated using the “lsqnonlin” built-in function in MATLAB. The case fatality
20
21 rate is estimated using the linear regression method.
22 The initial conditions used Table 3. The initial conditions for susceptible
23
24 individuals are an approximate total population in West Java. The fitted values
25
of β1 , β2 and q. The values are then used in the numerical simulation. The plot
26
27 of model’s solution and data is given in Figure 3. With these parameter values,
28
29 the reproduction number for West Java R0 = 3.180126127. This means that
30 an outbreak happens and the control needs to be implemented to minimize the
31
32 risk of infections.
33
34
35 6. Numerical Simulation
36
37 This numerical simulation is designed to support the results of the analy-
38
39 sis discussed in the previous section. We set the parameter by curve fitting
40
41 from actual case of COVID-19 in West Java Province, Indonesia. We applied
42 Runge-Kutta-Fehlberg (RKF) method in MAPLE software, to solve the ordi-
43
44 nary differential equations of model 1-7 using the parameters in Table 2 and
45 Table 3. RKF method is one of the most popular numerical approach because
46
47 it is quite accurate, stable, to program [51].
48
Figure 4 show the endemic incidence where the susceptible population (S0 )
49
50 decreases as a result of transmission from the symptomatic and asymptomatic
51
52 infected population. Hereafter, this increases the latent population (E), the
53 asymptomatic infected population (IA ), the symptomatic infected population
54
55 (YS ), the recovered population (Z), the susceptibility to previously infected
56 populations (Z0 ), and the quarantine population (Q).
57
58
59 16
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Table 2: Parameters Values
10
11
12 Parameter Value Unit Source
13 107
14 Λ 365×65 people × day −1 Estimated
15 −7 −1
β1 1.727 ×10 (people × day) Fitting
16 −1
17 β2 7.478 ×10−8 (people × day) Fitting
18 µ 1
day −1 Estimated
19 365×65

20 µ1 0.082 day −1 Fitting

21
22 α 1
5.2 day −1 [52]
23 p 0.2 N/A Assumed
24
25 κ 0.19 day −1 [53]
26 q 1.026 ×10 −6
day −1
Fitting
27
28 γ1 1
10 day −1 [52]
29 1 −1
γ2 14 day Assumed
30
31 δ 0.1 day −1 [54]
32 −1
33 ξ 0.02 day Assumed
34
35 Table 3: Initial Values of each Compartments
36
37
38 Compartment S(0) E(0) Ia(0) Y s(0) Z(0) Z0 (0) Q(0)
39
40 Initial Values 107 100 100 100 100 100 5
41
42
43 However, after the 20th day, the latent (E) and asymptomatic infected pop-
44
45 ulation (IA ) decreased, this is because the latent population and the asymp-
46
tomatic infected population became the removed population. While symp-
47
48 tomatic human populations (YS ) have declined due to an increase in quarantined
49
50 populations (Q), the recovered (Z) and susceptible that previously infected pop-
51 ulations (Z0 ) have consequently increased.
52
53 Figure 5 shows the number of quarantine population Q and the cumulative
54 population of quarantine. The number of quarantine compartment populations
55
56 increases as the asymptomatic infection increases. The peak occurs at 30 days
57
58
59 17
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
(a)
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36 (b)
37
38 Figure 4: Dynamical Population of each Compartment: (a) Population of S, Z0, Z (b) Pop-
39 ulation of E, Y s, & Ia
40
41
42 where the number reaches 70 and after that decreases. At the end of the 400th
43
44 day, the number of cumulative quarantine reaches 3200.
45
46
47 6.1. The Effect of Waning Immunity
48 In this section, we simulate the sensitivity analysis for the effect of param-
49
50 eter ξ, related to waning immunity issue, which describes the probability rate
51
52 of recovered people become susceptible, and the probability rate of suscepti-
53 ble people that previously infected become asymptomatic infected, respectively.
54
55 Using the parameters and initial values in Table 2 and Table 3, except for ξ, we
56 choose ξ = 0.001, 0.01, 0.1, 1.
57
58
59 18
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21 (a) Population of Active Quarantine
22
23
24
25
26
27
28
29
30
31
32
33
(b) Population of Cumulative Quarantine
34
35
Figure 5: Dynamical Population of (a) Active Quarantine (b) Cumulative Quarantine
36
37
38
39 Figure 6 and Figure 7 show the effect of increasing the value of ξ and τ . In
40 these simulations the peak time of disease spread do not change, but at the time
41
42 after the peak has been passed, the more value of ξ and τ multiply the number
43 of Asymptomatic infected population (YS ).
44
45 The effect of changes in the value of the probability rate of recovered people
46
become susceptible (ξ) on the E, IA , YS , and Z compartments is shown in Fig-
47
48 ure 7. The changes value of the parameter ξ did not have a significant impact on
49
50 the number of compartments E and IA . The number of populations E(t) and
51 IA (t) is relatively unchanged for every ξ ∈ [0.001, 0.1]. This means that changes
52
53 in the reinfected parameter value do not really affect the number population
54 Exposed (E) and asymptomatic infected population (IA ).
55
56 However, the higher the value of the parameter ξ, the higher the population
57
58
59 19
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
(a) Population of Active Quarantine
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39 (b) Population of Cumulative Quarantine

40
41 Figure 6: Simulation of The effect of Waning Immunity (ξ) on (a)Symptomatic Infected
42 Population (YS ) and (b) Quarantine Population (Q).
43
44
45 of YS after passing the peak of the spread, and the lower the population Z.
46
47 When this parameter is greater, the recovered population (Z) decreases due to
48 the loss of immunity and returns to the susceptible population that previously
49
50 infected (Z0 ). Where population Z0 can be re-infected to become asymptomatic
51 infected population (YS ).
52
53
54
55
56
57
58
59 20
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
(a) (b)
22
23
24
25
26
27
28
29
30
31
32
33
34
35 (c) (d)
36
37 Figure 7: Simulation of The effect of Waning Immunity ξ with respect to time (t) for each
38 Compartments (a) E, (b) IA , (c) YS , and (d) Z
39
40
41 6.2. The Effect of Quarantine
42
43 Figure 8 show that with increase the value of quarantine parameter (q), the
44
45 peak size and the final size of symptomatic infected population (YS ) is decrease.
46 This show that Increasing the intensity of quarantine policy may press the spread
47
48 of COVID-19. Figure 9 show that changes in the value of quarantine parameters
49
50 to the population of each compartments E, IA , YS , and Z are presented in 3-
51 dimensional changes in time. When the value of the quarantine parameter q is
52
53 increased, within the range [0.01, 1], the number of infected populations can be
54 reduced. This is indicated by the reduction in the peak value of Exposed (E),
55
56 Asymptomatic Infected (IA ), and Symptomatic Infected (YS ), in the change in
57
58
59 21
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 the value of q. Meanwhile the population in the quarantine compartment (Q)
10
is increasing from population of Symptomatic Infected which did Quarantine.
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25 Figure 8: Dynamical population of Symptomatic Infected Popu-
26
lation (YS ) in changes of quarantine parameter (q).
27
28
29
30
31
32
33
34
35
36
37
38
39 (a) (b)
40
41
42
43
44
45
46
47
48
49
50
51 (c) (d)
52
53
Figure 9: Simulation of The effect of quarantine parameter (q) with respect to time (t) for
54
55 each Compartments (a) E, (b) IA , (c) YS , and (d) Q.
56
57
58
59 22
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 7. Discussion and Conclusion
10
11 We have formulated a mathematical model of COVID-19 transmission by
12
13 considering infected individuals with symptoms and asymptomatic, as well as
14
decreased immunity, validated with data from West Java Province, Indonesia.
15
16 The compartment-based model is formulated as a system of differential equa-
17
18 tions, where the population is divided into Susceptible Populations (S0 ), Ex-
19 posed Populations (E), Asymptomatic Infection Populations (IA ), Symptomatic
20
21 Infection Populations (YS ), Recovered Populations (Z), Susceptible Populations
22 previously infected (Z0 ), and Quarantine Population (Q). Then the model is
23
24 analyzed mathematically, the results show that there are two equilibrium points,
25
26 namely a disease-free equilibrium point and an endemic equilibrium point. Be-
27 sides, with the next-generation matrix method, the Basic Reproduction Number
28
29 (R0 ) for West Java Province is obtained of 3.180126127. This means that West
30 Java Province is affected by the COVID-19 outbreak and controls are needed to
31
32 minimize the risk of transmitting COVID-19. Stability and sensitivity are ana-
33
lyzed to determine the parameters that influence the spread of COVID-19. The
34
35 simulation results show that the factor of decreasing immunity can affect the
36
37 spread of COVID-19. This is because when the increase in immunity decreases,
38 the infected population increases. Meanwhile, reinfection has no significant ef-
39
40 fect on the number of exposed and infected asymptomatic populations and the
41 isolation period can slow the spread of COVID-19 in West Java Province, In-
42
43 donesia. The results obtained can be used as a reference for the early prevention
44
45 of the spread of COVID-19 in West Java.
46
47
48 Acknowledgments
49
50 The work was supported by Universitas Padjadjaran, with contract number
51
52 1735/UN6.3.1/LT/2020 through Hibah Riset Data Pustaka dan Daring.
53
54
55
56
57
58
59 23
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Authors’s Contributions
10
11 N Anggriani, R Amelia, & W Suryaningrat contributed to the study de-
12
13 sign, model formulation, model analysis, and numerical simulation. MZ Ndii
14
designed sensitivity analysis and performed the case study including parameter
15
16 estimation. MAA Pratama complete and verify the analysis. All authors have
17
18 read and agreed to the published version of the manuscript.
19
20
21 Conflict of interest
22
23 The authors declare that they have no conflict of interest.
24
25
26 References
27
28
29 [1] Q. Li, X. Guan, P. Wu, X. Wang, L. Zhou, Y. Tong, ... & Z. Feng. Early
30 transmission dynamics in Wuhan, China, of novel coronavirus–infected
31
32 pneumonia, New England journal of medicine., 382 (2020), 1199-1207.
33
34 [2] Worldometer, COVID-19 coronavirus pandemic, (2020). Avaliable from:
35
36 https://www.worldometers.info/coronavirus/
37
38 [3] J. F. W. Chan, S. Yuan, K. H. Kok, K. K. W. Yo, H. Chu et al., A familial
39
40 cluster of pneumonia associated with the 2019 novel coronavirus indicating
41 person-to-person transmission: a study of a family cluster, Lancet, 395
42
43 (2020), 514–523.
44
45 [4] N. Chen, M. Zhou, X. Dong, J. Qu, F. Gong et al., Epidemiological and
46
47 clinical characteristics of 99 cases of 2019 novel coronavirus pneumonia in
48 Wuhan, China: a descriptive study, Lancet, 395 (2020), 507–513.
49
50
51 [5] C. Huang, Y. Wang, X. Li, L. Ren, J. Zhao et al., Clinical features of
52 patients infected with 2019 novel coronavirus in Wuhan, China, Lancet,
53
54 395 (2020), 497–506.
55
56
57
58
59 24
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 [6] Centers for Disease Control and Prevention, Situation
10
Summary COVID-19, 2020. Retrieved March 23, 2020
11
12 (https://www.cdc.gov/coronavirus/2019-nCoV/index.html)
13
14 [7] B. Shayak, M. M. Sharma, M. Gaur, & A. K. Mishra. Impact of reproduc-
15
16 tion number on multiwave spreading dynamics of COVID-19 with tempo-
17
18 rary immunity: a mathematical model, International Journal of Infectious
19 Diseases(2021). Preprint at doi: 10.1016/j.ijid.2021.01.018
20
21
[8] S. K. Law, A. W. N. Leung, & C. Xu., Is reinfection possible after recovery
22
23 from COVID-19?, Hong Kong Medical Journal, 26 (2020), 264–265.
24
25 [9] P. Selhorst, S. V. Ierssel, J. Michiels, J. Mariën, K. Bartholomeeusen, E.
26
27 Dirinck, S. Vandamme, H. Jansens,& K. K. Ariën.Symptomatic SARS-
28
CoV-2 reinfection of a health care worker in a Belgian nosocomial out-
29
30 break despite primary neutralizing antibody response, Clinical Infectious
31
32 Diseases, (2021). Preprint at doi.org/10.1093/cid/ciaa1850
33
34 [10] M. Galanti & J. Shaman. Direct Observation of Repeated Infections with
35
Endemic Coronaviruses, The Journal of Infectious Diseases , 3 (2020), 409-
36
37 415.
38
39 [11] H. Ledford. Coronavirus Re-infections: Three Questions Scientists Are Ask-
40
41 ing, Nature Research Journals, 7824 (2020), 168-169.
42
43 [12] C. Giannitsarou, S. Kissler, & F. Toxvaerd. Waning Immunity and the Sec-
44
45 ond Wave: Some Projections for SARS-CoV-2, CEPR Centre for Economic
46 Policy Research, (2020).
47
48
49 [13] S. Kumar, R. Kumar, C. Cattani, & B. Samet, Chaotic behaviour of frac-
50 tional predator-prey dynamical system, Chaos, Solitons & Fractals, 135
51
52 (2020), 109811.
53
54 [14] S. Kumar, S. Ghosh, B. Samet, & E. F. D. Goufo, An analysis for heat
55
56 equations arises in diffusion process using new Yang-Abdel-Aty-Cattani
57
58
59 25
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 fractional operator, Mathematical Methods in the Applied Sciences, (2020),
10
1-19.
11
12
13 [15] S. Kumar, R. Kumar, R. P. Agarwal, & B. Samet, A study of frac-
14 tional Lotka-Volterra population model using Haar wavelet and Adams-
15
16 Bashforth-Moulton methods, Mathematical Methods in the Applied Sci-
17
18 ences, (2020), 1-15.
19
20 [16] S. Kumar, A. Kumar, B. Samet, & H. Dutta, A study on fractional
21
host–parasitoid population dynamical model to describe insect species, Nu-
22
23 merical Methods for Partial Differential Equations, (2020), 1-20.
24
25 [17] N. Anggriani, H. Tasman, M.Z. Ndii, A.K. Supriatna, E. Soewono, E Sire-
26
27 gar, The effect of reinfection with the same serotype on dengue transmission
28
dynamics, Applied Mathematics and Computation, 349 (2019), 62-80.
29
30
31 [18] M. Z. Ndii, Z. Amarti, E. D. Wiraningsih, A. K. Supriatna, Rabies epidemic
32 model with uncertainty in parameters: crisp and fuzzy approaches, IOP
33
34 Conference Series: Materials Science and Engineering, 332 (2018), 012031.
35
36 [19] S. Kumar, A. Kumar, B. Samet, J. F. Gmez-Aguilar,& M. S. Osman, A
37
38 chaos study of tumor and effector cells in fractional tumor-immune model
39 for cancer treatment, Chaos, Solitons & Fractals, 141 (2020), 110321.
40
41
42 [20] B. Ghanbari, S. Kumar, & R. Kumar, A study of behaviour for immune and
43 tumor cells in immunogenetic tumour model with non-singular fractional
44
45 derivative, Chaos, Solitons & Fractals,133 (2020), 109619.
46
47 [21] M. Z Ndii, A. K. Supriatna, Stochastic mathematical models in epidemiol-
48
49 ogy, Information, 20 (2017) 6185-6196.
50
51 [22] H. Alrabaiah, Mo. A. Safi, M. H. DarAssi, B. Al-Hdaibat, S. Ullah, M. A.
52
53 Khan, S. A. A. Shah, Optimal control analysis of hepatitis B virus with
54 treatment and vaccination, Results in Physics, 19 (2020), 103599
55
56
57
58
59 26
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 [23] M. A. A. Oud, A. Ali, H. Alrabaiah, S. Ullah, M. A. Khan, & S. Islam, A
10
fractional order mathematical model forCOVID-19 dynamics with quaran-
11
12 tine,isolation, and environmental viral load, Advances in Difference Equa-
13
14 tions, 106 (2021), 1-19.
15
16 [24] Z. Liu, P. Magal, O. Seydi, & G. Webb., Predicting the Cumulative Number
17
of Cases for the COVID-19 Epidemic in China from Early Data. Preprint
18
19 at doi: 10.20944/preprints202002.0365.v1 (2020a).
20
21 [25] Z. Liu, P. Magal, O. Seydi, & G. Webb., Understanding Unreported Cases
22
23 in the COVID-19 Epidemic Outbreak in Wuhan, China, and the Importance
24 of Major Public Health Interventions, Biology 9(2020b).
25
26
27 [26] M. Pierre and G. Webb, The Parameter Identification Problem for SIR
28 Epidemic Models: Identifying Unreported Cases, Journal of Mathematical
29
30 Biology 77 (2018).
31
32 [27] K. Biswas, A. Khaleque, & P. Sen, COVID-19 spread: Reproduction of
33
data and prediction using a SIR model on Euclidean network, pp. 4–7,
34
35 2020, [Online]. Available: http://arxiv.org/abs/2003.07063.
36
37 [28] K. Biswas & P. Sen, Space-time dependence of corona virus
38
39 (COVID-19) outbreak, no. 1, pp. 5–7, 2020, [Online]. Available:
40 https://arxiv.org/abs/2003.03149.
41
42
43 [29] M. A. Khan & A. Atangana, Modeling the dynamics of novel
44 coronavirus(2019-nCov) with fractional derivative,Alexandria Engineering
45
46 Journal, 59(2020),2379-2389.
47
48 [30] M. A. Khan, A. Atangana, & E. Alzahrani, The dynamics of COVID-19
49
with quarantined and isolation, Advances in Difference Equations, 1 (2020),
50
51 1-22.
52
53 [31] M. S. Alqarni, M. Alghamdi, T. Muhammad, A. S. Alshomrani, & M.
54
55 A. Khan, Mathematical modeling for novel coronavirus (COVID-19) and
56
control, Numerical Methods for Partial Differential Equations, (2020).
57
58
59 27
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 [32] Y. Chu, A. Ali, M. A. Khan, S. Islam, & S.Ullah, Dynamics of fractional or-
10
der COVID-19 model with a case study of Saudi Arabia,Results in Physics,
11
12 21 (2021), 103787.
13
14 [33] S. Kumar, R. Kumar, S. Momani,& Samir Hadid, A study on fractional
15
16 COVID-19 disease model by using Hermite wavelets, Mathematical Methods
17
18 in the Applied Sciences, (2021), 1-17.
19
20 [34] M. A. Khan, S. Ullah, & S. Kumar, A robust study on 2019-nCOV out-
21
breaks through non-singular derivative. The European Physical Journal
22
23 Plus, 136 (2021), 1-20.
24
25 [35] K. M. Safare, V. S. Betageri, D. G. Prakasha, P. Veeresha, & S. Kumar,
26
27 A mathematical analysis of ongoing outbreak COVID-19 in India through
28
nonsingular derivative, Numerical Methods for Partial Differential Equa-
29
30 tions, 37 (2021), 1282-1298.
31
32 [36] S. Kumar, R. P. Chauhan, S. Momani, & Samir Hadid, Numerical inves-
33
34 tigations on COVID-19 modelthrough singular and non-singularfractional
35
operators, Numerical Methods for Partial Differential Equations, (2020),
36
37 1-27.
38
39 [37] Q. Lin, S. Zhao, G. Daozhou, Y. Lou, S. Yang, et al., A conceptual model
40
41 for the coronavirus disease 2019 (COVID-19) outbreak in Wuhan, China
42
43 with individual reaction and governmental action, International Journal of
44 Infectious Diseases., 93 (2020), 211–216.
45
46 [38] B. Ivorra, M. R. Ferrandez, M. V. Perez, & A. M. Ramos., Mathematical
47
48 modeling of the spread of the coronavirus disease 2019 (COVID-19) taking
49
50 into account the undetected infections. The case of China, Communications
51 in Nonlinear Science and Numerical Simulation. 88 (2020).
52
53
[39] B. Tang, F. Xia, S. Tang, N. L. Bragazzi, Q. Li, et al., The effectiveness of
54
55 quarantine and isolation determine the trend of the COVID-19 epidemic in
56
57
58
59 28
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 the final phase of the current outbreak in China, International Journal of
10
Infectious Diseases, 96 (2020), 288-293.
11
12
13 [40] J. P. Arcede, R. L. C. Anan, C. Q. Mentuda, & Y. Mammer., Accounting
14 for Symptomatic and Asymptomatic in a SEIR-type model of COVID-19.
15
16 Preprint at https://arxiv.org/abs/2004.01805 (2020).
17
18 [41] M. Z. Ndii, A.R. Mage, J.J. Messakh, & B.S. Djahi. Optimal vaccina-
19
20 tion strategy for dengue transmission in Kupang city, Indonesia. Heliyon,
21
6(2020),e05345.
22
23
24 [42] H. Sun, Y. Qiu, H. Yan, Y. Huang, Y. Zhu, et al. ., Tracking
25 and Predicting COVID-19 Epidemic in China Mainland. Preprint at
26
27 https://www.medrxiv.org/content/10.1101/2020.02.17.20024257v1 (2020).
28
29 [43] S. Ullah & M. A. Khan. Modeling the impact of non-pharmaceutical inter-
30
31 ventions on the dynamics of novel coronavirus with optimal control analysis
32 with a case study. Chaos, Solitons & Fractals, 139 (2020), 110075.
33
34
35 [44] N. Anggriani, A.K. Supriatna, & E. Soewono. A Critical Protection Level
36 Derived from Dengue Infection Mathematical Model Considering Asymp-
37
38 tomatic and Symptomatic Classes. Journal of Physics Conference Series.
39 423 (2013), 1-10.
40
41
42 [45] V. Kontis, J. E. Bennett, T. Rashid, R. M. Parks, J. Pearson-Stuttard,
43 M. Guillot, P. Asaria, B. Zhou, M. Battaglini, G. Corsetti, M. McKee, M.
44
45 Di Cesare, C. D. Mathers and M. Ezzati, Magnitude, demographics and
46 dynamics of the effect of the first wave of the COVID-19 pandemic on all-
47
48 cause mortality in 21 industrialized countries, Nature Medicine, 26 (2020),
49
50 1919–1928.
51
52 [46] O. Diekmann & J. Heesterbeek., Mathematical Epidemiology of Infectious
53
Diseases, Wiley Series in Mathematical and Computational Biology, (2000).
54
55
56
57
58
59 29
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 [47] C. Chavez, C., Z. Feng, W. Huang., On the computation of R0 and its role
10
in global stability. In: Math ApproachesEmerg Reemerg Infect Dis Intro.
11
12 IMA, 125 (2002), 229–250.
13
14 [48] J. J. Tewa, J. L. Dimi, and S. Bowong,. Lyapunov functions for a dengue
15
16 disease transmission model, Chaos, Solitons, & Fractals, 39 (2009), 936-
17
18 941.
19
20 [49] J.P. LaSalle., The stability of dynamical systems, Philadelphia: SIAM,
21
(1976).
22
23
24 [50] S. Marino, I. B. Hogue, C. J. Ray, D. E.Kirschner, A methodology for
25 performing global uncertainty and sensitivity analysis in systems biology,
26
27 Journal of Theoretical Biology, 254 (2008) 178-198.
28
29 [51] A.M. Islam, A comparative study on numerical solutions of initial value
30
31 problems (IVP) for ordinary differential equations (ODE) with Euler and
32 Runge Kutta Methods, American Journal of Computational and Applied
33
34 Mathematics,5 (2015) 393–404.
35
36 [52] A. J. Kucharski, T. W. Russell, C. Diamond, Y. Liu, J. Edmunds,S. Funk,
37
38 R. M. Eggo, F. Sun, M. Jit, J. D. Munday, N. Davies, A. Gimma, K. van
39 Zandvoort, H. Gibbs, J. Hellewell, C. I. Jarvis, S. Clifford, B. J. Quilty, N.
40
41 I. Bosse, S. Abbott, P. Klepac, S. Flasche, Early dynamics of transmission
42
43 and control of COVID-19: a mathematical modelling study, The Lancet
44 Infectious Diseases, 20 (2020), 553-558.
45
46 [53] W. J. Guan , Z. Y. Ni, Y. Hu, W. H. Liang, C. Q. Ou, J. X. He, et al.,
47
48 Clinical Characteristics of Coronavirus Disease 2019 in China, New Engl.
49
50 J. Med., 382 (2020), 1708-1720.
51
52 [54] D. Aldila, S. H.A. Khoshnaw, E. Safitri, Y.R. Anwar, A.R.Q. Bakry, B. M.
53
Samiadji, D. A. Anugerah, M. F. Farhan GH, I.D. Ayulani, S. N. Salim, A
54
55 mathematical study on the spread of COVID-19 considering social distanc-
56
57
58
59 30
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 ing and rapid assessment: The case of Jakarta, Indonesia, Chaos, Solitons
10
and Fractals, 139 (2020), 110042.
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 31
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Appendix A. Proof of numerical analysis of Theorem 3.3
10
11 Coefficient polynomial (P1 (λ)) :
12
13
14 a0 =1,
15
16 a1 =1.030393017,
17 a2 =0.3892774217,
18
19 a3 =0.06516694579,
20
21 a4 =0.004442996019,
22
23 a5 =0.00006597971860, &
24
25 a6 =2.773132296 × 10−9
26
27
28
29 Value of:
30
31
a1 a2 = 0.4011087370,
32
33 a0 a3 = 0.06516694579,
34
35 a1 (a2 a3 + a0 a5 ) = 0.4683242878,
36
37 a21 a4 + a0 a23 = 0.008963903101,
38
39 a1 a2 a4 = 0.001782124518,
40
41 a0 (a1 a6 + a2 a5 ) = 0.00002568727211.
42
43
44
45
46 Appendix B. Characteristic polynomial of Exposed Compartment
47
48
49
50 √
G = K + (((−µ1 2 γ + ((−q − γ)δ − κ γ − γ 2 )µ1 + (−γ 2 + (−κ − q)γ − κ q)δ)µ
51
52 + β Λ µ1 2 + (δ β Λ + Λ β γ + (κ + qp)Λ β)µ1 + (Λ β γ + (2 qp + κ)Λ β)δ)ξ
53
54 + (β Λ µ1 2 + (δ β Λ + Λ β γ + (κ + qp)Λ β)µ1 + (Λ β γ + (κ + qp)Λ β)δ)µ)α
55
56 + (−µ1 3 γ + ((−q − γ)δ − κ γ − γ 2 )µ1 2 + (−γ 2 + (−κ − q)γ − κ q)δ µ1 )µ ξ,
57
58
59 32
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 where
10
11
12
13 K =((µ1 γ + (q + γ)δ)2 (κ + γ + µ1 )2 (α + µ1 )2 ξ 2 + 2 (κ + γ + µ1 )(δ + µ1 )β Λ α (α + µ1 )
14
15 (−µ1 2 γ + ((−γ + (2 p − 1)q)δ − γ (κ + γ + qp))µ1 + δ (q + γ)(qp − κ − γ))ξ
16
17 + β 2 α2 Λ2 (δ + µ1 )2 (κ + γ + qp + µ1 )2 )µ2 + 2 ξ β Λ α ((κ + γ + µ1 )(−µ1 3 γ+
18
((−2 γ + (2 p − 1)q)δ − γ (κ + γ + qp))µ1 2 + ((−γ + (2 p − 1)q)δ − 2 γ 2
19
20 + ((−p − 1)q − 2 κ)γ + q(qp − κ))δ µ1 − δ 2 (κ + γ)(q + γ))(α + µ1 )ξ
21
22 + (µ1 2 + (κ + γ + δ + qp)µ1 + δ (κ + γ))(δ + µ1 )β Λ α (κ + γ + qp + µ1 ))µ
23
24 + ξ 2 (µ1 2 + (κ + γ + δ + qp)µ1 + δ (κ + γ))2 β 2 Λ2 α2 ,
25
26
and
27
28
29 H ={2 α β δ pqξ (µ1 + α)}.
30
31
.
32
33 By substituting the parameter value from Table 2 we have endemic equilib-
34
35 rium point:
36
37 E = 1536.477836,
38
39 IA = 158.8583370,
40
41 Q = 0.01936422900,
42
43 S0 = 679.0681919,
44
YS = 3434.980193,
45
46 Z = 17931.49917,
47
48 Z0 = 577.7849236.
49
50
51
52
So, the non-endemic equilibrium point is stable if it exists, because it satisfies
53
54 Routh-Hurwit’s criterion. These results remain consistent using the parameter
55
56 values in the Table 2, except β1 = β2 = 0.7477942169036 × 10−10 .
57
58
59 33
60
61
62
63
64
65
 1
2
3
4
5
6
7
8
9 Appendix C. Proof of numerical analysis of endemic equilibrum point
10
11 Coefficient polynomial (Q(λ)) :
12
13
14 a0 = 1.000000000,
15
16 a1 = 2.171964310,
17
a2 = 1.885854850,
18
19 a3 = 0.84368388,
20
21 a4 = 0.20819336,
22
23 a5 = 0.0279316,
24
25 a6 = 0.001468969,
26
27 a7 = 0.00001178844740.
28
29 .
30
31 Value of:
32
33
a1 a2 = 4.096009428,
34
35 a0 a3 = 0.8436838800,
36
37 a1 (a2 a3 + a0 a5 ) = 3.516403565,
38
39 a1 2 a4 + a0 a23 = 1.693939876,
40
41 a1 a2 (a3 a4 + a0 a7 ) = 0.7195098093,
42
43 a0 a3 (a1 a6 + a2 a5 ) = 0.04713281469,
44
45 a2 a5 = 0.05267494333,
46 a0 a7 = 0.00001178844740.
47
48
49
50
51
52
53
54
55
56
57
58
59 34
60
61
62
63
64
65
Conflict of Interest

Conflict of Interest and Authorship Conformation Form

Please check the following as appropriate:

All authors have participated in (a) conception and design, or analysis and
✓ interpretation of the data; (b) drafting the article or revising it critically for
important intellectual content; and (c) approval of the final version.
This manuscript has not been submitted to, nor is under review at, another journal
✓ or other publishing venue.

The authors have no affiliation with any organization with a direct or indirect
✓ financial interest in the subject matter discussed in the manuscript

The following authors have affiliations with organizations with direct or indirect
✓
financial interest in the subject matter discussed in the manuscript:

Author’s Name Affiliation

Nursanti Anggriani Department of Mathematics, Universitas Padjadjaran,
Indonesia
Meksianis Z. Ndii Department of Mathematics, The University of Nusa
Cendana, Indonesia
Rika Amelia Department of Mathematics, Universitas Padjadjaran,
Indonesia
Wahyu Suryaningrat Departemen of Mathematics, Universitas Padjadjaran,
Indonesia
M Andhika Aji Pratama Departemen of Mathematics, Universitas Padjadjaran,
Indonesia