CC BY
11
Simulation of the First and the Second Waves of COVID-19 Spreading in Russian Federation Regions Using an Agent-Based Model
Mikhail Kirillin1*, Aleksandr Khilov1+, Valeriya Perekatova1, Ekaterina Sergeeva1, Daria Kurakina1, Ilya Fiks1, Nikolay Saperkin1,2, Ming Tang3, Yong Zou3, Elbert Macau4, and Efim Pelinovsky1,5
1 Institute of Applied Physics RAS, 46 Ul'yanov str., Nizhny Novgorod 603950, Russia
2 Privolzhsky Research Medical University, 10/1, Minin and Pozharsky Sq., Nizhny Novgorod 603950, Russia
3 School of Physics and Electronic Science, East China Normal University, 3663 N. Zhongshan Rd., Shanghai 200062, China
4 Instituto de Ciencias e Tecnologia, Universidade Federal de Sao Paulo, R. Sena Madureira, 1500 - Vila Clementino, Sao Paulo 04021-001, Brasil
5 National Research University - Higher School of Economics, 25/12 Bolshaya Pecherskaya str., Nizhny Novgorod 603155, Russia
e-mail: *[email protected], +[email protected]
Abstract. The COVID-19 pandemics remains one of the largest global challenges. Necessity of effective systemic aids for the minimization of losses leads to the requirement of adequate models allowing to predict the impact of different factors on the spread of the disease. Agent-based simulation models provide a suitable solution with the possibility to accurately account for such factors as age structure of a population, characteristics of isolation, self-isolation strategies and testing strategies, presence of super-spreaders etc. In this paper we report on the results of simulating the spread of COVID-19 in several representative regions of Russia using an agent-based model with a general pool combined with the simulation of population testing strategy. The model accounts for the following key epidemiologic characteristics: population age distribution, reproducibility rate, distributions of infectivity period, a period of clinical manifestation, and age-dependent probability of critical disease. It is demonstrated that the daily epidemiologic curves can be predicted well for different territories with the same model parameters, except for the initial number of infected agents and region-dependent testing as well as isolation strategies, which are considered to be tuning parameters of the model. The developed approach can be further expanded to other regions of different countries, while the determined model parameters could be used as starting values for such simulations. © 2023 Journal of Biomedical Photonics & Engineering.
Keywords: dynamics and control of epidemics; COVID-19; agent-based modeling.
Paper #3568 received 02 Dec 2022; revised manuscript received 21 Dec 2022; accepted for publication 21 Dec 2022; published online 13 Feb 2023. doi: 10.18287/JBPE23.09.010302.
1 Introduction
Prediction of further development of the COVID-19 pandemic and timely introduction of preventive measures require reliable tools for epidemic spread simulations. Several classes of models have been currently employed
for forecasting the spread of infections. Even the simplest logistic model based on the first order ordinary differential equation has provided surprisingly reasonable accuracy in describing the first wave of COVID-19 in different countries [1-4]. Moreover,
employment of its predecessor, the Gompertz model, provides even better agreement of the predicted values with the official statistical data [5]. Regression models are known to provide us with rapid estimations of the spread of infections [6, 7]. Non-adaptive models allow obtaining the prognosis for any chosen period of time, however, they ignore local perturbations of epidemiologic characteristics, therefore they are not suitable for short-term prognosis. On the contrary, adaptive models are predominantly employed for short-term prognosis only. Autoregressive moving average models are applicable for short-term prognosis [8], while autoregressive integrated moving average model is applicable both for short-term and long-term prognosis [9]. Dynamic Bayesian networks are employed only for short-term prognosis and predominantly in the form of Markov models [10]. Neural networks and other machine learning based methods can be applied only for short-term prognosis. Furthermore, only feedforward neural networks and backpropagation algorithm can address the prognosis of infections spread [11].
Dynamic systems based on differential equations, which belong to a class of compartmental models, are employed for long-term prognosis. The pioneering work by Kermack and McKendrick [12] proposed SIR model, in which the population is divided into three groups -susceptible (S), infected (I), recovered (R) - and their interaction is described with non-linear differential equations. Modern compartmental models employ larger number of groups, accounting for exposed (E), hospitalized (H), critical (C), dead (D), and those in quarantine (Q) or isolation (J). Compartmental models are widely employed for the modeling of the spread of coronavirus infections (SARS, MERS, COVID-19). For example, SEIR (susceptible, exposed, infective, recovered) model has been applied for modeling of COVID-19 spread in African countries [13], Iran [14], Indonesia [15], Spain and Italy [16], India [17, 18]. The early COVID-19 spread and the efficacy of governmental measures are discussed in paper [19] also with regard to a SEIR model. Modified SEIR models may also account for vaccination process [18, 20]. Compartmental models are easy to construct, however, they do not capture random factors of the infections spreading nor individual characteristics of population members.
Individually-oriented models include so-called agent-based models, which can be applied both for short-term and long-term prognosis of the spread of infection. Every member of a population (referred to as an agent) is described by certain constant and variable characteristics, and the rules of interactions between the agents are determined. Agent-based models were proved to be effective in the description of the propagation of infections, such as Ebola virus disease [21] and influenza [22], regarding different sizes of population.
Agent-based models were also applied for the modeling of COVID-19 spread, for example, the modeling of the development and the regress of the infection in the various cities such as Helsinki [23], New York [24], and Singapore [25]. Certain agent-based
models were based on previously developed models for the prediction of influenza pandemics, for example, the NotreDame-FRED model [26], the model of Ferguson's research group from Imperial College London [27], and the model for simulation of the COVID-19 pandemic in Australia [28].
A number of studies utilized agent-based models for assessing the impact of universal face mask wearing [29], digital contact tracking [30-32] and social distancing [32-34] as well as analyzing various intervention scenarios [24, 27, 32]. Other agent-based models were developed for modeling COVID-19 transmission in small communities such as a university [35], a supermarket [36], and a small town [37].
It is worth noting, that the validation strategy for a COVID-19 spread model should be based on the agreement with the official statistical data on daily new cases and daily deaths. These "measurement results" may have certain errors connected with the imperfection of the employed "measurement instrument", namely, a testing system. Building a model of a testing system is a particularly important task, since the simulated "measurement results" taken in consideration are actually the combined results of both the model of the disease spread and the model of testing system. Improper model of the testing system may result in incorrect results of the entire model, although the disease spread part can be valid. Moreover, this is the feedback system, since the test positive result affect the isolation strategy and, hence, further disease spread.
The aim of this paper is the development of an agent-based model capable of simulating the progress of the COVID-19 burst in different regions of the Russian Federation. Another important issue to be resolved in this study is the determination of the key model parameters that can provide the agreement of the simulated dynamics and actual statistics on daily new infection cases and deaths associated with COVID-19. The developed approach includes a testing strategy model based on official data concerning daily tests performed.
2 Materials and Methods
2.1 Agent-Based COVID-19 Spread Model
2.1.1 Single Agent States
The agent-based model operates on a set of elements ^ where each element an e ^ is characterized by Boolean variables representing binary states OnJ, where n is the element number, and j is the binary state. The elements are traditionally called agents and each of them represents a person in the considered population. The number of elements in ^ usually attributed as cardinality N = is chosen in accordance with the population of the considered region. The agent-based model employed in this paper is a general pool model, which is schematically shown in Fig. 1. In the general pool model,
Fig. 1 Basic schematic of the employed agent-based model with a general pool.
/', days
infected
contagious
manifestation of symptoms
V
critical state
t„5 end of critical period
recovery start
1 1 fn <Pi ■ _ ^ death L > Pc (fin) in >Pm I >
Fig. 2 Schematic of a general disease scenario.
all the agents that leave their homesteads on a particular day, interact with other agents, and have the equal probability to be infected. The latter is determined by the number of infected agents, that are not isolated and also go out to the "pool".
Each agent has basic properties (age h„) and the following Boolean states: "infected" (an1), "contagious" (an2), "symptomatic" (an3), "in critical state" (an4), "dead" (a„5), "recovered" (an6), "isolated" (an7), "positive test" (a„8) that are set as 'false' (0) at the initialization of the simulation and are further governed by the Monte Carlo principle depending on a scenario. The basic properties of each agent are predefined randomly in accordance with the official statistical data on the population age structure.
If for an agent n the state "infected" becomes 'true' (an1 = 1) on day i = 41, a specific random scenario for each infected person is generated based on the predefined probability density distributions f (t).
The schematic of a general disease scenario is shown in Fig. 2. The probabilities of symptoms manifestation pm and critical period pc(hn) are defined separately, and the latter depends on an individual's age h„. Any agent with the status "in critical state" (an4 = 1) has a chance to die on any day of critical period with the given death probability pd. The following scenario time points (tnk)
are determined in accordance with predefined probability density distributions fk(x) based on available statistical data on disease progression: day when person becomes contagious (tn2), symptoms manifestation day (tn3), critical state day (tn4), critical state overcome day (tn5), and recovery start day (tn6),
tnk = tn1 + Fk-1(?„), k = 2, 3; (1)
tnk = tnk-1 + Fk-1(?„), k = 4, 5; (2)
t6 I fn + y L> Pc(hn)
. fn + td + tn - fn 1 In < Pc(hn)
where tdis the typical disease duration, , ^ are random numbers uniformly distributed on the interval [0, 1], Fk(t) is the cumulative distribution function of corresponding parameter:
Fk(t) = jA/fc(T)dT,
(4)
and Fk_1(^n) represents the inverse function. In turn, these parameters determine the corresponding time points of the status change:
a„2 = 1 | i = tn2, (5)
a„3 = 1 | i = tn3, \n < Pm, (6)
4 _ J 1 1 ^ ^ni^n < Pm'^in < Pc (hn)
*-nJ
0 | i = pm < pc(hn)
(7)
where , are random numbers uniformly distributed on the interval [0, 1], in here a single call for each random number value is supposed.
2.1.2 Implementation of the General Pool Concept
In each day the set of daily recovered agents is determined among those agents who reached the recovery start day in accordance with the rule:
9Î: { an | i > tn6, an5 = 0, < pr},
an6 = 1, an2 = 0, an3 = 0, an7 = 0,
an8 = 0 | an£ 5R.
(8)
(9)
The daily death probability pd during the critical period is predefined but can vary with time reflecting the current situation with medical care: an overflow of the region hospitalization capacity may decrease in the quality of provided treatment therefore resulting in an increase in the probability of a death. The set of daily deaths D is determined in accordance with the following rule:
£>: { an | an4 = 1, < pd}, an5 = 1, | an £ D.
(10) (11)
where pr is the daily recovery probability constructed in a similar way to the death probability. In contrast to SIR-like models [12-20], where the total numbers of susceptible, infected and recovered persons are governed by differential equations, in the agent-based model, we consider the evolutions of the sets of susceptible agents (S), infected agents (3), recovered agents Cft), and dead agents (D) together with the set of contagious agents (£) that can transmit infection to other agents (£ £ 3) as follows:
The isolation of agents in the developed model rests on two mechanisms: the first is the self-isolation, which means that the agent does not go out to the general pool on a particular day basing either on official governmental recommendations, or on early onset of a disease; the second mechanism implies forcing an agent to stay at home or hospitalization in the case of COVID-positive test result. As so, depending on the introduced anti-epidemic measures accounted in the simulation, agents also may get a binary status "isolated" (an7 = 1) if a positive COVID-19 test happens or due to self-isolation, if they fulfill the restrictive measures on a current day. Isolated agents do not interact with the general pool and cannot transmit infection. The developed simulation model accounts for the efficiency of the following restrictive measures with the employment of what is known as self-isolation index. This is an empirical value similar to one firstly introduced by Yandex (Russia) during the COVID-19 outbreak, which represents a cumulative parameter reflecting population social activity based on both traffic information and activities in different internet services. In the simulations, the self-isolation index (Is) varies between 0 and 5. Day-scale time resolution allows us to account for the self-isolation index variations associated with vacations, New Year holidays etc. We introduce, and the quarantine quality k which is assumed to be proportional to Is and is, described by the empirical law:
0.5 + 0.7/,
k = -
5
(17)
This linear relation was constructed empirically in order to fit the official data. The set © of isolated agents that obey restrictive rules and do not interact with the general pool on a current day is determined based on current k value. Day-scale temporal resolution allows to account for self-isolation index variations associated with vacations, New Year holidays etc.
({an | an3 = 0 , Çn < k } U
{an I an3 = 1, Çn <Pe} U
{an | an4 = 1} U {an | an8 = 1}) \ £>,
' = 1 | an £ ©.
(18)
(19)
{an | an5 = 1}, (12)
S: {an | an1 = 0}\ (13)
3: {an | an1 = 1}\ ©, (14)
{an | an6 = 1}\ ©, (15)
£: {an | an2 = 1}\ ©. (16)
The disease transmission rate is governed by the parameter r that is defined as an average number of agents infected by a spreader in the population during one week given that no restrictions or isolation aids are imposed. Given that the probabilities to be infected from different agents are independent, the daily infection probability pis for a given susceptible agent interacting with the general pool in a particular day i is calculated as:
pis = (r |Œ n(CT\Œ)|)/(7N).
(20)
n
a
n
The probability of being infected p^ is equal for all agents who visit the pool in the current day (those who are not isolated or self-isolated, the probability is determined in accordance with the current self-isolation index). Hence, the set of agents infected on a particular day i is determined as:
3: {an | an £ S\®, <ps},
an1 = 1 | an £ 3,
(21) (22)
where %n is a random number uniformly distributed on the interval [0, 1].
2.1.3 Connection with SIRD-Like Models
The proposed approach can be directly related to the typical SIRD-like models. The variables considered in SIRD-like models could be derived as the cardinalities of the considered agent sets:
S = |©|, (23)
I = |3|, (24)
R = |K|, (25)
D = |©|, (26)
which dynamics can be evaluated from basic equations resembling the discrete form of the classical SIRD model:
Ir = Ii-1 + AI, - AR, - AD„ (27)
S, = Si-1 - A/, + 1, (28)
Ri = Ri-1 + ARi, (29)
D, = Di-1 + AD„ (30)
where i denotes current day number, while
A/, = |3|, (31)
AR, = |W|, (32)
AD, = |©|, (33)
are the daily numbers of infected, recovered and dead agents, respectively.
2.1.4 Testing System Model
An important part of the simulation is the testing system, since the daily statistics on the number of incident cases and deaths of numerical simulations have to be compared to the real data. Obviously, these data depend on the region-specific testing strategies, which may
significantly vary between different regions. In addition, they could serve as tuning parameters of the entire model. In the simulations the tested group X is determined each day and consists of three parts: agents with symptoms agents infected by the infected agents U revealed on the previous day, and random agents 3 from the rest of population:
X = M u U u 3, M: Äp(|M|,{aJ «n = 1}), U: Äp(|U|,{aJ = 1}),
(34)
(35)
(36)
3: flp(|3l, (SU3U 5R)\{oJ < = 1}), (37)
where (/, is a random permutation of j elements from set The numbers of agents in each of subgroups U, 3 are defined in accordance with the following rules:
|OT| = d |{an | an3 = 1, an8 = 0}|,
|U| = Ncont |S3i-il,
131 = - |M| - |U|,
(38)
(39)
(40)
where parameter d characterizes the fraction of the tested agents among all symptomatic agents that have no positive test results, and Ncont describes the average number of infected agents tested as contacts of the agents with positive tests Q i-1 revealed on the previous day i-1. These parameters may vary depending on testing capacity and psychological status of population. The total number of tests N/est taken at day i is a predefined parameter, which can be obtained from the official statistics data. The parameters d and Ncont are empirical and are chosen to satisfy the condition:
|OT| + |U| < N/est.
(41)
The set of daily detected cases Q is determined based on the predefined test accuracy ac:
Q: {an | an £ X, < ac }, an8 = 1 | an £ Q.
(42)
(43)
It is also assumed in the model that all the agents with positive test results are subject to isolation until their recovery:
an1 = 1 | an £ Q.
(44)
The number of daily detected infected agents is the number of elements in Q:
Aßi= |Q|.
(45)
The dynamics AlO'), AD(i), AQi(i) are the major dependencies derived in simulations: AA(i), AQ,{i) could be compared to the official statistical data, while AIi(i) reveals the real progression of the epidemic.
2.2 Model Parameters Used in the Simulations
The probability distributions for different parameters of a disease scenario (Fig. 3) are based on the available statistical data on COVID-19 [38-41] and determine the probabilities and durations of each period of the disease: period before being contagious, incubation period, disease development period, critical period, and the probability of critical state occurrence depending on the agent's age. The simulations reported in this study consider only the first two waves of the pandemic (period from February 26th, 2020 to April 06th, 2021), at which
0.30
0.25
0.20
0.10 0.05 0.00
1 2 3 4 5 Period before contagious status, x = t^ - days
0.30
0.250.20-
^ 0.15
0.100.05-
0.00
Disease development period, t = t4 - t3 days
0.1,
0.01,
0.001,
the effect of vaccination omitted in the described model could be considered to be negligible. Moreover, it is assumed that the abovementioned parameters do not vary with time, since both these waves are caused by the COVID strains with similar contamination abilities, except death occurrence probability during the critical stage, since it also depends on the availability and quality of medical care.
These distributions were already used in our previous study [42] aiming at predicting the progression of the first pandemic wave. In this study, the simulations of the epidemic progress were performed for Moscow (N = 11.4106), Nizhny Novgorod region (N = 3.3-106) and Novosibirsk region (N = 2.7-106).
6 8 10 12 14
0 - t1
Incubation period,t3 - days
5 6 7 8 9 10 11 12 13 Critical state, t = t5 - t^, days
0-10 11-17 18-34 35-50 51-65 66+ Age group hn
Fig. 3 Probability density distributions of the disease scenario parameters: period from infection to the contagious status x = tn2-tnl (a), incubation period x = tn3-tnl (b), period from the disease manifestation to the critical state x = t„4-t„3 (c), critical state period x = 45-44 (d), and the age-dependent probability of critical state occurrence pc(hn) (e) employed in simulations.
4
4
5
6
7
8
Table 1 Values of main parameters employed in simulations.
Parameter
Abbreviation Value
Probability of symptoms manifestation pm 0.3
Daily death probability in critical state Pd 0.04
Average weekly number of agents infested by a spreader r 4.6
Probability of isolation decision for agents with symptoms manifestation
0.5
The reproducibility rate value was determined in the course of fitting the official statistics for the first wave only [42], resulting in the value of r = 4.6 for all the considered regions. The number of initial infected agents for all regions was the same as previously reported in Ref. [42].
It was assumed that both the first and the second wave were caused by the COVID strains with similar contamination abilities, so this parameter was kept constant in all the simulations. The same assumption concerns the probability of symptoms manifestation pm which was predefined 0.3, which corresponds to 70% of silent (symptomless) cases typically reported in the beginning of the pandemic. Death probability pd was considered to be time-independent for all the regions, but its value was different for different regions, which can be associated with difference of treatment quality. It is assumed that the statistics of new deaths strongly correlates with the statistics of new cases, and the pd coefficient was determined by finding the best fits time series after the fit of the series for daily new cases for each region.
In order to fit the increase of COVID-19 associated lethal cases in Moscow in the period from 10 February 2021 to 06 April 2021 we increase the pd value for Moscow within this period. This may originate from deficiency of bed capacity or increased fatality ratio due
7x105-
to 6x105-.2
5x105-
4x105 -
I 3x105-c .
2x105" <0 .
Q 1x105-
a
2020
to virus mutations. The values of the key probabilistic parameters employed in simulations are summarized in Table 1.
In our simulations the testing model is based on the available data of daily tests for the entire country (Fig. 4a), since the availability of the data for specific regions is much lower. Consequently, the number of daily tests for each region is either taken from official statistics [43] or determined from the daily number of cases for the entire Russian Federation (Fig. 4a) [44] in proportion to the particular region population. Fig. 4b shows the relevance of such an evaluation for available data for the Nizhny Novgorod region and compares estimated numbers with actual data obtained for the period 10.10.2020 - 10.04.2021. The fraction of tested agents with symptoms d was assumed to be equal to zero for first 25 days of simulation (26 February 2020 -23 March 2020) indicating the time required to develop a testing system, and then shifts to a non-zero value and increases with the number of day i according to the following law:
d(0 =
0,i < 25
dn
io + (df -do)-3-3ii
i +id
> 25
(46)
where d0 is initial value for i > 25, df is asymptotic value being reached for i » id, id is an empirically chosen parameter.
The parameter Ncont determining the average number of infected persons tested in the frames of the testing contacts of the agents Q i-1 who got the positive test on the previous day i - 1, was assumed to be either constant through the entire simulation period, or to have a step-wise decrease in the end of the second wave, since the first wave period was characterized by the maximal efforts to reveal all the cases.
It is worth mentioning that the tests accuracy ac was rather low in the beginning of the pandemic, and, in this respect, the test sensitivity is assumed to increase with time in the beginning of the pandemic. Similar to the dynamics of d, it was assumed that ac is described by the following law:
2021
a/ •*•*>.• ».
0
20.02.
2020
20.05.2020
20.08.2020
Date
20.11.2020
20.02.2021
1.8x104 1.6x104
V)
■§5 1.4x104
"5 1.2x104
jj 1.0x104
H 8.0x103 c
.S^6.0x103 Q 4.0x103 2.0x103
10.10.2020
10.12.2020 10.02.2021 Date
10.04.2021
Fig. 4 Daily number of tests in Russian Federation (RUS) in the period 20 February 2020 - 20 April 2021 (a) and daily number of tests in Nizhny Novgorod region (NN) versus evaluation in the period 10 October 2020 - 10 April 2021 (b). Ratio Pnn/Prus is the ratio of populations of Nizhny Novgorod region and Russian Federation.
e
4.03.5-_ 3.02.5-
<0° 2.0 -
ti
1.51.00.5-
d
- N™
\ /
V B/-
I
_!_
0.0-_
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
2.0x104
0.07 0.06 0.05 0.04 0.03 0.02
0.0
-AQ, MC d
-A/, MC
■ AQ, Official
k Jfci"^!
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
200-,
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
2 25 j 2.0
r /s ' d ■ ac Ncont Pd b
lH-"—
■
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
aD, MC g 25-
■ aD, Official 20-
A i 4: _ 15- Q < 10-
A / 5-
0-
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
3.53.021 25 " cc' 2.0
1.51.00.5-
" ac
Pd
- Nrn
T_Tr
0.0-(_......
10.02.2020 10.06.2020 10.10.2020 10.02.2021
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
18 16 14 12- 10 a « 8
64 2010.02.2020 10.06.2020 10.10.2020 10.02.2021
Date
Fig. 5 Empirically chosen simulation parameters (a, b, c) and comparison of simulated scenarios (MC) and official statistical data (Official) for daily numbers of newly revealed COVID-19 (d, e, f) and COVID-19-associated deaths (g, h, i) in Moscow (a, d, g), Nizhny Novgorod region (b, e, h) and Novosibirsk region (c, f, i) in the period of 26 February 2020 - 06 April 2021.
0,i < 25
ac(0 =
c,0 + (ac,/ - ac,o)p~3,£ > 25
(47)
where i is the number of simulation day starting from 26 February 2020, initial ac,o and asymptotic acj values are chosen empirically, as well as iT. The specificity of the tests is assumed to equal 100%, and no false positive test results are considered in simulations.
To provide a best-fit scenario, we manipulated with parameters of the model, namely, primarily the self-isolation index Is and testing strategy parameters d, Ncont and ac. Each scenario was constructed by simulation of 10 random scenarios with the same parameters and further averaging five scenarios that gave the minimal deviation from the official data. Due to the stochastic origin of the developed model, similar starting parameters may result in totally different scenarios of an epidemic progress.
3 Results
This study includes the analysis for three regions: Moscow, Nizhny Novgorod region, and Novosibirsk
region. Fig. 5 shows the best fit scenarios obtained for these regions. The dependencies of the key time-dependent parameters for the considered regions are shown in Figs. 5a-c. At the first step the parameters for the Moscow region (Fig. 5a) were determined. The variation of the self-isolation index is the most significant among other parameters and its dynamics is regulated by the dates of municipal and governmental decrees and social processes: introduction of the quarantine regime (1), release of the restrictions (2) followed by further cancellation of restrictive measures (3), restoration of restrictive measures after the second wave start (4), a winter-time isolation increase related to New Year holidays (5), and its release (6). The parameters of the testing system are quite hard to be accurately extracted from the official data, so they were empirically fitted based on general conclusions. The test accuracy ac was assumed to be quite low (ac,0 = 0.4) in the beginning of the testing period and rose up to 90% (acf = 0.9) during two months of the testing period (4). In Moscow the rise of the test accuracy (Fig. 5a) is assumed to be faster, since the newly developed tests arrive first of all to the capital, which is followed by their distribution to other regions. Similar trend is chosen for the percentage d of the
a
c
0.11
2
3
400
1000
200
0
150-
Q
< 100-
50-
0
a
symptomatic agents with symptoms manifestation that are being tested (3). It is assumed, that in the very beginning of the pandemic only few of them are tested, however, with its progression the testing strategy changes both due to governmental requirements and people's awareness. The latter decrease in the d parameter for Moscow (Fig. 5a) simulates smaller attention to early symptom manifestation.
Figs. 5d-f show the simulated dynamics of the new daily COVID-19 counts, while Figs. 5g-i present the simulated dynamics of the daily COVID-19 associated lethal cases for Moscow, Nizhny Novgorod region and Novosibirsk region, respectively. Official statistical data for daily counts and daily deaths that was fitted by manipulating the simulation parameters (Figs. 5a-c) are given in each plot for the reference in Figs. 5d-i; scatter data show the official statistical data, while solid lines show the simulated dynamics. The new daily counts dynamics AQi show good agreement with the official statistics, especially for Moscow and Novosibirsk region, while an insignificant discrepancy is observed for the beginning of the first wave for Nizhny Novgorod region. Agreement of the simulation results with the official statistics data (Figs. 5d-f) demonstrate the relevance of the developed model. Although the deviations may occur in particular days, the general trend is reproduced quite accurately for all three considered regions (Figs. 5d-f) with the same parameter r governing disease properties and associated with the reproducibility rate. An important advantage of the employed approach is the ability to predict the real number of new cases, including those that were not tested and thus did not contribute to the official statistics. In Figs. 5d-f solid green lines demonstrate calculated dynamics of daily detected cases, while red solid lines show the dynamics of the evaluated total new COVID-19 cases. The simulations show a 1.5-2 fold difference between the evaluated numbers of the officially confirmed new daily counts and the evaluated number of total daily cases.
It is important to note that in general the daily deaths time series ADi follow the trend of the daily new cases AQi with the delay equal to the average period between first manifestation and death, as can be seen from official statistics data. Since in accordance with official statistics the death rate does not exceed 5% [45, 46], the number of daily new COVID-19-associated deaths was simulated by varying the parameter pd only.
The observed trend is reproduced in all the considered regions except the second wave for Moscow, where the daily death statistics (Fig. 5g) shows a plateau, which is not observed in the simulated daily new cases graph (Fig. 5d). In this connection, a discrepancy occurs between the model predictions and official statistics (Fig. 5g), although for other considered regions (Figs. 5h, i) the agreement is much better. It is also worth mentioning that the time shift between the dependencies AQi and ADi is smaller for Moscow as compared to Novosibirsk and Nizhny Novgorod, which is most probably due to delay in the representation of COVID-19 associated death cases in official statistics in the regions with respect to Moscow,
since the difference in disease progression in different regions is not expected. This discrepancy could be compensated by introducing a corresponding delay to the model, however, it was avoided to maximally keep the basic parameters of the model for different regions.
Despite the similarity of the key model parameters for the first wave of the pandemic in the chosen regions, for the second wave different trends were revealed to govern the wave progression in different sites: for Moscow an increased death probability was introduced in the end of the second wave (Fig. 5a), while for Nizhny Novgorod region the number of tested contacts Ncont was decreased (Fig. 5c). In Moscow city the summer relaxation period followed by the beginning of the second wave was related to a higher decrease in isolation index as compared to other regions. In Supplementary materials we present a collection of the simulation results which do not show the good fit to the official data, however allow to see the effect of particular parameters on the model predictions.
1200 1000800-
O 600-
400-
200-
0
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
360030002400-
o
< 1800 A
1200-
600-
0
10.02.2020 10.06.2020
10.10.2020 Date
10.02.2021
Fig. 6 Comparison of simulated scenarios (MC) and real statistical data (Official) for daily numbers of newly revealed COVID-19 cases for different percentage of tested agents (different dynamics of fraction of tested agents with symptoms d) with symptoms manifestation (a) and different number of tested agents, which were in contact with an agent with positive COVID-19 test (b) in the period of 26 February 2020 - 06 April 2021.
Date
Date
Fig. 7 Comparison of simulated scenarios and real statistical data (Official) for daily numbers of newly revealed COVID-19 (a) and COVID-19-associated deaths (b) for different self-isolation index Is in Nizhny Novgorod region in the period of 26 February 2020 -06 April 2021.
In order to analyze the impact of testing strategy parameters on the epidemic process we performed simulations for Nizhny Novgorod region varying parameters d and Ncont. We limit ourselves by only one region in this part of the study, since we expect that the general trends will be similar for different regions owing to similarity of the model parameters, since the considered variations are higher than those between different regions. The results for the daily revealed cases AQi for different values of the parameters d and Ncont are shown in Figs. 6a and 6b, respectively.
In simulations presented in Fig. 6a, the evolution of the fraction of agents with symptoms manifestation d that are tested is described by Eq. 3, similar to previous simulations with a fixed value of d0 = 0.35, while df value is varied between 0.7 and 1.0. The increase in df value leads to higher detection rate and, hence, a faster isolation of the infected agents and a corresponding decrease in the number of the revealed and total COVID-19 cases. Note that even the total testing of all the agents with symptoms manifestation (df = 1) does not provide full elimination of the pandemic progression owing to both high percentage
of symptomless disease progression and limited accuracy of the tests.
Fig. 6b shows simulated dynamics of AQi for different number of tested infected agents Ncont that were previously contacting with the revealed cases. Similar to variations in d, the increase in Ncont leads to decrease in daily detected COVID-19 cases and vice versa. A pronounced rise in AQi for Ncont = 2 (red line) with a further sharp decrease indicate the herd immunity development, when most agents were infected and then recovered or died. The decrease of AQi to zero for Ncont = 3 and Ncont = 4 indicates the situation when all the infected agents were revealed and isolated thus preventing further development of the pandemic.
Fig. 7 demonstrates the effect of the basic level of the self-isolation index Is on the dynamics of daily revealed and lethal COVID-19 associated cases AQi and ADi calculated for Nizhny Novgorod region. The results are shown for basic self-isolation index dynamics /sase (as shown in Fig. 5c) and its variations (± 0.1 and ± 0.3 from base values for all days during considered period).
1.0x104-8.0x103-6.0x103-
o
< ■
4.0x103-2.0x103-
0.010.02.2020 10.06.2020 10.10.2020 10.02.2021
Date
300250200-
Q~ 150-
<
10050-
010.02
Fig. 8 Comparison of simulated scenarios and real statistical data (Official) for daily numbers of newly revealed COVID-19 (a) and lethal COVID-19 associated (b) cases in Nizhny Novgorod region for different probability of symptoms manifestation pm in the period of 26 February 2020 - 06 April 2021.
The simulated dynamics show clear evidence that the increase of the self-isolation index leads to the decrease
■ Official b
Pm = 0.1
pm = 0.3 (Best fit)
Pm = 05
-Pm = 07
A
• : ■
I-"-■<—---1-1-1-■-1-
.2020 10.06.2020 10.10.2020 10.02.2021
Date
in social activity and, hence, slows down the epidemic process and vice versa. Note that similar to the case of Ncont = 2, the decrease of Ishase for the value of 0.3 induces a sharp increase in the number of daily revealed cases followed by a sharp drop indicating the formation of herd immunity, while similarly to the case of Ncont = 4 (Fig. 6b), a suppression of epidemic progression is observed for the increase of the isolation level (ibase + 0.3).
Another important parameter of the model is the probability of symptoms manifestation pm. Fig. 8 shows the results of simulations for Nizhny Novgorod region
with basic parameters (Fig. 5c) and different manifestation probability pm.
As it can be seen form Fig. 8a, the increase in the percentage of symptoms manifestation leads to a larger number of potentially tested agents and, therefore, to the timely isolation of infected agents. With the significant decrease of symptoms manifestation probability, the pandemic scenario is represented by a quick spread of infection and the development of herd immunity (Fig. 8a, red curve), accompanied by huge number of lethal cases (Fig. 8b, red curve).
0.05
0.04
0.03
0.02
0.01
0.00
0.07
10.02.2020 10.06.2020
AQ, MC A/,- MC AQ, Official
10.10.2020 Date
10.02.2021
2x104-1x104-0
10.02.2020 10.06.2020
200
10.10.2020 Date
10.02.2021
150-
< 100-
50-
0
10.02.2020 10.06.2020
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
AQ,- MC A/, MC AQ, Official
1.4x104 -T
1.2x104-
1.0x104 -
< 8.0x103-
O
<1 6.0x103 -
4.0x103-
2.0x103 -
0.0-
10.02.2020 10.06.2020 10.10.2020 10.02.2021 Date
150-
125100-
Q 7CJ
< 75-
50 25
10.10.2020 Date
10.02.2021
0
10.02.2020 10.06.2020
10.10.2020 Date
10.02.2021
Fig. 9 Empirically chosen alternative simulation parameters (a, b) and comparison of corresponding simulated scenarios (MC) and official statistical data (Official) for daily numbers of newly revealed COVID-19 (c, d) and COVID-19 -associated deaths (e, f) in Moscow in the period of 26 February 2020 - 06 April 2021.
<B tf
10.06.2020
AQ, MC A/, MC AQ, Official
0.0 10.02.2020
1800 160014001200-<f 1000800-1 6004002000
10.02.2020 10.06.2020
10.10.2020 Date
10.02.2021
10.10.2020 Date
10.02.2021
10.02.2020 10.06.2020
10.10.2020 Date
10.02.2021
(B Ti
0.15
- 0.14
0.13
10.02.2020 10.06.2020
10.10.2020 Date
10.02.2021
500
400-
AQ, MC A/, MC AQ, Official
300-
<
cf
<
200-
100-
0
10.02.2020 10.06.2020
18
161412108642-
10.10.2020 Date
10.02.2021
010.02.2020
10.06.2020
10.10.2020 Date
10.02.2021
Fig. 10 Empirically chosen alternative simulation parameters (a, b) and comparison of corresponding simulated scenarios (MC) and official statistical data (Official) for daily numbers of newly revealed COVID-19 (c, d) and COVID-19-associated deaths (e, f) in Novosibirsk region in the period of 26 February 2020 - 06 April 2021.
Since the considered simulation of pandemic progression is a multiparameter problem, several solutions could be found to fit the given official statistics data. Abovementioned results were obtained for realistic parameters, e.g. relatively smaller self-isolation index Is during the second wave of COVID-19 in comparison to the first one. Nonetheless, different solutions were found for Moscow (Fig. 9) and Novosibirsk region (Fig. 10) fitting official statistics well, and suggesting comparable
isolation index during the first and the second waves of the pandemic (Figs. 9a, b and 10a, b) or even higher isolation index during the second wave (Fig. 9d, e and 10d, e). However, these scenarios implied a significant variation of pd and d parameters during the second wave (Figs. 9c and 9f).
As can be seen from Figs. 9 and 10, the major difference of these scenarios from those shown in Fig. 5 lies in the total number of daily COVID-19 cases, which is significantly
higher than routinely revealed COVID-19 cases (demonstrating from 2 to 9-fold difference). This situation indicates that although the proposed model provides good agreement with the reported official statistics, additional information regulating model parameters is to be employed for more accurate simulation and prediction of further pandemic progression.
4 Conclusion
In this paper we presented an agent-based model of COVID-19 epidemic spread with the employment of Monte Carlo simulation principles. The model is able to account for the age-dependent disease development, restrictive measures as well as testing system. The results are encouraging and were validated on the statistical data for daily new cases and deaths which were reported during the first and the second waves of COVID-19 pandemics from three representative regions of Russia, such as the Moscow city, Nizhny Novgorod region, and Novosibirsk region. The results of simulations show the possibilities of developed model for the simulation of disease spread with regard to analysis of restrictive measurements and different testing strategies.
It is worth noting that different combinations of simulations parameters could provide similar results in the dynamics of newly revealed cases and deaths, although the dynamics of total cases in such scenarios may differ significantly. This indicates, that some parameters of the model should be defined independently based on available information on social activity/testing strategy, etc., which could increase the accuracy of the prognosis provided by the model.
The reported results demonstrate the abilities of the developed model in simulating the progression of the epidemics, however, owing to a multiparametric nature of the considered problem, the choice of the parameters remains an important issue in performing simulations. The approach used in this study is based on employing available parameters, such as statistics on population, percentage of symptomless cases, number of daily tests etc., and finding the empirical values for other parameters that govern social activity issues, that cannot be extracted directly. It is worth noting that the parameters found for particular regions could be used as starting sets for another regions providing fast solution.
As it was mentioned in the model description, the general pool concept could be connected with the classical compartmental model, however, consideration of individual agents gives more flexibility to the model allowing to avoid considering numerous small compartments that would be required in case of application of a compartmental model in similar settings and detalization.
Disclosures
The authors have no relevant financial interest in this article and no conflict of interest to disclose.
Acknowledgements
The study is supported by RFBR (project no. 20-51-80004), CNPq (project no. 441016/2020-0), and NSFC (project no. 82161148012).
References
1. T. Carletti, D. Fanelli, and F. Piazza, "COVID-19: The unreasonable effectiveness of simple models," Chaos, Solitons & Fractals: X 5, 100034 (2020).
2. E. Pelinovsky, A. Kurkin, O. Kurkina, M. Kokoulina, and A. Epifanova, "Logistic equation and COVID-19," Chaos, Solitons & Fractals 140, 110241 (2020).
3. K. Wu, D. Darcet, Q. Wang, and D. Sornette, "Generalized logistic growth modeling of the COVID-19 outbreak: comparing the dynamics in the 29 provinces in China and in the rest of the world," Nonlinear Dynamics 101, 1561— 1581 (2020).
4. A. Cunha, F. da C. Batista, P. R. de L. Gianfelice, R. S. Oyarzabal, J. M. V. Grzybowski, and E. E. N. Macau, "epidWaves: A code for fitting multi-wave epidemic models," Software Impacts 14, 100391 (2022).
5. E. Pelinovsky, M. Kokoulina, A. Epifanova, A. Kurkin, O. Kurkina, M. Tang, E. Macau, and M. Kirillin, "Gompertz model in COVID-19 spreading simulation," Chaos, Solitons & Fractals 154, 111699 (2022).
6. H. S. Burkom, S. P. Murphy, and G. Shmueli, "Automated time series forecasting for biosurveillance," Statistics in Medicine 26(22), 4202-4218 (2007).
7. J. C. Brillman, T. Burr, D. Forslund, E. Joyce, R. Picard, and E. Umland, "Modeling emergency department visit patterns for infectious disease complaints: results and application to disease surveillance," BMC Medical Informatics and Decision Making 5(1), 4 (2005).
8. D. Lai, "Monitoring the SARS epidemic in China: a time series analysis," Journal of Data Science 3(3), 279-293 (2005).
9. J. Díaz-Hierro, J. J. Martín Martín, Á. Vilches Arenas, M. P. López Del Amo González, J. M. Patón Arévalo, and C. Varo González, "Evaluation of time-series models for forecasting demand for emergency health care services," Emergencias 24(3), 181-188 (2012).
10. S. Unkel, C. P. Farrington, P. H. Garthwaite, C. Robertson, and N. Andrews, "Statistical methods for the prospective detection of infectious disease outbreaks: a review," Journal of the Royal Statistical Society: Series A (Statistics in Society) 175(1), 49-82 (2012).
11. R. Kiang, F. Adimi, V. Soika, J. Nigro, P. Singhasivanon, J. Sirichaisinthop, S. Leemingsawat, C. Apiwathnasorn, and S. Looareesuwan, "Meteorological, environmental remote sensing and neural network analysis of the epidemiology of malaria transmission in Thailand," Geospatial Health 1(1), 71-84 (2006).
12. W. O. Kermack, A. G. McKendrick, "A contribution to the mathematical theory of epidemics," Proceedings of the Royal Society of London. Series A 115(772), 700-721 (1927).
13. Z. Zhao, X. Li, F. Liu, G. Zhu, C. Ma, and L. Wang, "Prediction of the COVID-19 spread in African countries and implications for prevention and controls: A case study in South Africa, Egypt, Algeria, Nigeria, Senegal and Kenya," Science of The Total Environment 729, 138959 (2020).
14. N. Ghaffarzadegan, H. Rahmandad, "Simulation-based estimation of the early spread of COVID-19 in Iran: actual versus confirmed cases," System Dynamics Review 36(1), 101-129 (2020).
15. S. Annas, M. I. Pratama, M. Rifandi, W. Sanusi, and S. Side, "Stability analysis and numerical simulation of SEIR model for pandemic COVID-19 spread in Indonesia," Chaos, Solitons & Fractals 139, 110072 (2020).
16. L. López, X. Rodo, "A modified SEIR model to predict the COVID-19 outbreak in Spain and Italy: simulating control scenarios and multi-scale epidemics," Results in Physics 21, 103746 (2021).
17. X. Yu, G. Qi, and J. Hu, "Analysis of second outbreak of COVID-19 after relaxation of control measures in India," Nonlinear Dynamics 106, 1149-1167 (2020).
18. S. Kurmi, U. Chouhan, "A multicompartment mathematical model to study the dynamic behaviour of COVID-19 using vaccination as control parameter," Nonlinear Dynamics 109(3), 2185-2201 (2022).
19. Y. Fang, Y. Nie, and M. Penny, "Transmission dynamics of the COVID-19 outbreak and effectiveness of government interventions: A data-driven analysis," Journal of Medical Virology 92(6), 645-659 (2020).
20. A. Temerev, L. Rozanova, O. Keiser, and J. Estill, "Geospatial model of COVID-19 spreading and vaccination with event Gillespie algorithm," Nonlinear Dynamics 109(3), 239-248 (2022).
21. C. Siettos, C. Anastassopoulou, L. Russo, C. Grigoras, and E. Mylonakis, "Modeling the 2014 Ebola virus epidemic-agent-based simulations, temporal analysis and future predictions for Liberia and Sierra-Leone," PLoS Currents 7, (2015).
22. H. Arduin, M. Domenech de Celles, D. Guillemot, L. Watier, and L. Opatowski, "An agent-based model simulation of influenza interactions at the host level: insight into the influenza-related burden of pneumococcal infections," BMC Infectious Diseases 17(1), 382 (2017).
23. J. T. Tuomisto, J. Yijóla, M. Kolehmainen, J. Bonsdorff, J. Pekkanen, and T. Tikkanen, "An agent-based epidemic model REINA for COVID-19 to identify destructive policies," medRxiv: 20047498 (2020).
24. N. Hoertel, M. Blachier, C. Blanco, M. Olfson, M. Massetti, F. Limosin, and H. Leleu, "Facing the COVID-19 epidemic in NYC: a stochastic agent-based model of various intervention strategies," medRxiv: 20076885 (2020).
25. J. R. Koo, A. R. Cook, M. Park, Y. Sun, H. Sun, J. T. Lim, C. Tam, and B. L. Dickens, "Interventions to mitigate early spread of SARS-CoV-2 in Singapore: A modelling study," The Lancet Infectious Diseases 20(6), 678-688 (2020).
26. G. España, S. Cavany, "NotreDame-FRED COVID-19 forecasts," GitHub (accessed 18 November 2022). [https://github. com/confunguido/covid 19_ND_forecasting].
27. N. M. Ferguson, D. Laydon, G. Nedjati-Gilani, N. Imai, K. Ainslie, M. Baguelin, S. Bhatia, A. Boonyasiri, Z. Cucunubá, G. Cuomo-Dannenburg, A. Dighe, I. Dorigatti, H. Fu, K. Gaythorpe, W. Green, A. Hamlet, W. Hinsley, L. C. Okell, S. van Elsland, H. Thompson, R. Verity, E. Volz, H. Wang, Y. Wang, P. G. T. Walker, C. Walters, P. Winskill, C. Whittaker, C. A. Donnelly, S. Riley, and A. C. Ghani, "Impact of non-pharmaceutical interventions (NPIs) to reduce COVID-19 mortality and healthcare demand," Imperial College London (2020).
28. S. L. Chang, N. Harding, C. Zachreson, O. M. Cliff, and M. Prokopenko, "Modelling transmission and control of the COVID-19 pandemic in Australia," Nature Communications 11, 5710 (2020).
29. D. Kai, G.-P. Goldstein, A. Morgunov, V. Nangalia, and A. Rotkirch, "Universal masking is urgent in the COVID-19 pandemic: SEIR and agent based models, empirical validation, policy recommendations, arXiv:2004.13553 (2020).
30. J. A. Moreno López, B. Arregui García, P. Bentkowski, L. Bioglio, F. Pinotti, P.-Y. Boelle, A. Barrat, V. Colizza, and C. Poletto, "Anatomy of digital contact tracing: Role of age, transmission setting, adoption, and case detection," Science Advances 7(15), eabd8750 (2021).
31. Md. S. Shamil, F. Farheen, N. Ibtehaz, I. M. Khan, and M. S Rahman, "An Agent-Based Modeling of COVID-19: Validation, Analysis, and Recommendations," Cognitive Computation (2021).
32. A. Aleta, D. Martin-Corral, A. Pastore y Piontti, M. Ajelli, M. Litvinova, M. Chinazzi, N. E. Dean, M. E. Halloran, I. M. Longini, S. Merler, A. Pentland, A. Vespignani, E. Moro, and Y. Moreno, "Modeling the impact of social distancing, testing, contact tracing and household quarantine on second-wave scenarios of the COVID-19 epidemic," medRxiv:20092841 (2020).
33. P. C. L. Silva, P. V. C. Batista, H. S. Lima, M. A. Alves, F. G. Guimaräes, and R. C. P. Silva, "COVID-ABS: An agent-based model of COVID-19 epidemic to simulate health and economic effects of social distancing interventions," Chaos, Solitons & Fractals 139, 110088 (2020).
34. P. Keskinocak, B. E. Oruc, A. Baxter, J. Asplund, and N. Serban, "The impact of social distancing on COVID19 spread: State of Georgia case study," PLoS One 15(10), e0239798 (2020).
35. P. T. Gressman, J. R. Peck, "Simulating COVID-19 in a university environment," Mathematical Biosciences 328, 108436 (2020).
36. F. Ying, N. O'Clery, "Modelling COVID-19 transmission in supermarkets using an agent-based model," PLoS One 16(4), e0249821 (2021).
37. A. Truszkowska, B. Behring, J. Hasanyan, L. Zino, S. Butail, E. Caroppo, Z. Jiang, A. Rizzo, and M. Porfiri, "HighResolution Agent-Based Modeling of COVID-19 Spreading in a Small Town," Advanced Theory and Simulations 4(3), 2000277 (2021).
38. C. Cheng, D. Zhang, D. Dang, J. Geng, P. Zhu, M. Yuan, R. Liang, H. Yang, Y. Jin, J. Xie, S. Chen, and G. Duan, "The incubation period of COVID-19: a global meta-analysis of 53 studies and a Chinese observation study of 11 545 patients," Infectious Diseases of Poverty 10(1), 119 (2021).
39. W. Dhouib, J. Maatoug, I. Ayouni, N. Zammit, R. Ghammem, S. B. Fredj, and H. Ghannem, "The incubation period during the pandemic of COVID-19: a systematic review and meta-analysis," Systematic Reviews 10, 101 (2021).
40. S. Paul, E. Lorin, "Distribution of incubation periods of COVID-19 in the Canadian context," Scientific Reports 11(1), 12569 (2021).
41. H. Xin, J. Y. Wong, C. Murphy, A. Yeung, S. Taslim Ali, P. Wu, and B. J. Cowling, "The Incubation Period Distribution of Coronavirus Disease 2019: A Systematic Review and Meta-analysis," Clinical Infectious Diseases 73(12), 2344-2352 (2021).
42. M. Kirillin, E. Sergeeva, A. Khilov, D. Kurakina, and N. Saperkin, "Monte Carlo simulation of the covid-19 spread in early and peak stages in different regions of the Russian Federation using an agent-based modelling," in Saratov Fall Meeting, Chinese-Russian workshop on Biophotonics and Bioimaging-2020 1, 71-74 (2020).
43. Coronavirus-monitor - interactive map of COVID-19 spread and statistics (accessed 06 April 2021). [https://coronavirus-monitor.info/country/russia/nizhegorodskaya-oblast/, in Russian].
44. Newsletter on the situation and measures taken to prevent the spread of diseases caused by the new coronavirus, (accessed 06 April 2021). [https://rospotrebnadzor.ru/about/info/news/ news_details.php?ELEMENT_ID=19030&sphrase_id=4390258, in Russian].
45. E. Mathieu, H. Ritchie, L. Rodes-Guirao, C. Appel, D. Gavrilov, C. Giattino, J. Hasell, B. Macdonald, S. Dattani, D. Beltekian, E. Ortiz-Ospina, and M. Roser, "Coronavirus Pandemic (COVID-19)," Our World in Data (accessed 19 April 2022). [https://ourworldindata.org/mortality-risk-covid].
46. E. Dong, H. Du, and L. Gardner, "An interactive web-based dashboard to track COVID-19 in real time," Lancet Infectious Diseases 20(5), 533-534 (2020).
