key: cord-0727856-t6zboemd
authors: Dimou, Argyris; Maragakis, Michael; Argyrakis, Panos
title: A network SIRX model for the spreading of COVID-19
date: 2021-12-08
journal: Physica A
DOI: 10.1016/j.physa.2021.126746
sha: 1167ca123b65581cbedb4b58e171d8df8d88d914
doc_id: 727856
cord_uid: t6zboemd

Infectious diseases, such as the current COVID-19, have a huge economic and societal impact. The ability to model its transmission characteristics is critical to minimize its impact. In fact, predicting how fast an infection is spreading could be a major factor in deciding on the severity, extent and strictness of the applied mitigation measures, such as the recent lockdowns. Even though modelling epidemics is a well studied subject, usually models do not include quarantine or other social measures, such as those imposed in the recent pandemic. The current work builds upon a recent paper by Maier and Brockmann (2020), where a compartmental SIRX model was implemented. That model included social or individual behavioral changes during quarantine, by introducing state [Formula: see text] , in which symptomatic quarantined individuals are not transmitting the infection anymore, and described well the transmission in the initial stages of the infection. The results of the model were applied to real data from several provinces in China, quite successfully. In our approach we use a Monte-Carlo simulation model on networks. Individuals are network nodes and the links are their contacts. We use a spreading mechanism from the initially infected nodes to their nearest neighbors, as has been done previously. Initially, we find the values of the rate constants (parameters) the same way as in Maier and Brockmann (2020) for the confirmed cases of a country, on a daily basis, as given by the Johns Hopkins University. We then use different types of networks (random Erdős-Rényi, Small World, and Barabási-Albert Scale-Free) with various characteristics in an effort to find the best fit with the real data for the same geographical regions as reported in Maier and Brockmann (2020). Our simulations show that the best fit comes with the Erdős-Rényi random networks. We then apply this method to several other countries, both for large-size countries, and small size ones. In all cases investigated we find the same result, i.e. best agreement for the evolution of the pandemic with time is for the Erdős-Rényi networks. Furthermore, our results indicate that the best fit occurs for a random network with an average degree of the order of <k> [Formula: see text] 10–25, for all countries tested. Scale Free and Small World networks fail to fit the real data convincingly.

Spreading phenomena have been an active area of research due to the multitude of applications that they have historically been associated with. One prominent such application is the spreading of viruses in a population, whether the virus is a software produced entity being spread in a computer network [2] , or a physical organism that is causing a disease to a living species [3] . One such virus, COVID-19, produced the pandemic which currently affects the entire planet. It started up in a single point (presumably the region of Wuhan in China where the first case was reported), but was quickly spread throughout all continents with devastating results, both in the health sector, and the global economy, as reported by the World Health Organization [4] . Therefore, understanding the spreading mechanism is of utmost importance in an effort to initially contain the pandemic, and to subsequently eliminate it. A large number of potential strategies have been reported, especially during the past 20 years, aiming to mitigate the spreading process. Networks played a prominent role in this kind of studies, since interactions between people can be easily replicated in such kind of network structures [5] .

The SIR model and its variations has been the most commonly used scenario for disease spreading [6] [7] [8] [9] [10] . The variants include amongst others the SIS [11] , the SIRD SEIR [7] , and the SIRX model [1] . These models have been solved in the past both by Mean Field approaches, by writing down a set of coupled differential equations which take into account all paths of the disease with the proper rate constants, as well as by Monte-Carlo computer simulations. The properties of interest include the distribution of the infected population [12] , the methodology of immunization in order to contain the spreading [13] [14] , the inclusion and effects of several social measures, such as quarantine [15] [16] , social distancing [17] , etc. These theoretical models are sometimes coupled with the medical efforts to find a cure [18] , or find the most effective vaccine [19] . Thus, they will help establish the optimum approach to mitigate the negative impact of the disease to society. Even though epidemiological models exist for a very long time [20] , what needs to be described and examined in detail are the societal changes which the pandemic brings about, such as social distancing and lockdown measures. Such measures, despite being simple and existing for more than a century, have to be coupled assiduously with the models of the physics of spreading and medical/pharmaceutical treatments available to the population.

A model that coupled the above factors was recently given in [1] , describing successfully the spread of the infection in China for the first few days, which happened to be the point of origin of this specific virus. The methodology used was that of solving a set of four (4) differential equations. These were chosen to include containment policies, such as quarantine measures of symptomatic people and accounted for isolation of populations in geographical areas. The model shows a universal scaling law for the confirmed cases of different areas in China when compared with the real data, implying that the underlying spreading mechanism may include the social measures undertaken. Such measures have not been sufficiently captured in the past by the usual epidemiological models.

In the present study we use Monte-Carlo simulations performed on network structures to tackle the same problems. Networks offer a wide variety of connectivity options, depending on their type. We expect that the spreading of a virus on a network, whether this is Erdős-Rényi, Small World, or Scale-Free, will be different depending on its structure, as has been shown in the past [21] . Thus, the question we pose here is the following: is there any type of network to agree both qualitatively and quantitatively well with the real data of the virus spreading in China. The work described in [1] was able to capture the overall (macro) picture, but did not provide any information as to J o u r n a l P r e -p r o o f Journal Pre-proof what types of interactions between individual people (micro) are responsible for the spreading of the infection. Here, we use the three types of networks mentioned above and investigate which one best fits the real data.

In order to proceed with the different probabilities for the occurrence of the different events we borrow the rate constants used in [1] , which were based on comparison with the real data. Thus, we solve the same differential equations as described therein, and use the values of the rate constants found there. We then apply the model to a variety of different countries, other than China, using the same approach. These include large countries, such as Turkey, Romania, as well as small countries, such as Andorra, Barbados, etc. For all these cases we try to find the type of network that agrees best with the real data corresponding to each country, and subsequently to produce the values of parameters that give the best fit, such as the average degree, < k >, etc.

Our results show very good agreement with the results of the original works for China, when comparing the two different methods. We are also in agreement with the real data for all countries tested and mentioned above. All our calculations point to the fact that a substrate based on an Erdős-Rényi network can be reliably used to provide results similar to those of real data. Other types of networks examined show significant deviations from the real data. This implies that a connectivity of Erdős-Rényi type can describe with an sufficient approximation the outcome of disease spreading in a country. By using this methodology we are unable to identify the existence of "hub" type of superspreaders, which have been hypothesized in the past [12] .

Our current work is based on a modified version of the SIR model in networks, that of SIRX, as presented and implemented in [1] . SIRX is a compartmental model, which divides the nodes of the network into different categories according to their properties. A node could be susceptible to the infection (S), infected (I), recovered (R), or quarantined (X). State X denotes a less common variant of the more standard SIR epidemiological models, and it symbolizes all quarantined individuals. X is assumed to approach the cumulative number of confirmed COVID-19 cases. All data are tracked by the Center for Systems Science and Engineering (CSSE) at Johns Hopkins University, Baltimore [22] , and refer to the period from January 20 to June 4, 2020.

The analytical model is structured as follows: Initially, for each given region or country we extract from the data, the number of confirmed cases of the first day that any cases are publicly announced and we set it as the starting value of X. It must be noted that since it is assumed that not all infected individuals are symptomatic or have exhibited symptoms yet, the total number of infected individuals can in reality be larger than X in some occasions. In order to estimate the total number of infected individuals in the very first timestep, we use the relation I = I 0 f actor * X (similarly to the notation used in [1] ). I 0 f actor is the required factor that when multiplied by X it produces I, the estimated number of infected individuals for the first time step. Moreover, by assuming that the R nodes at t = 1 are still 0, we conclude that all the remaining individuals, those not I, and not X, will be S. To describe the virus transmission process, we allow the I individuals to infect their S neighbours with a certain rate. This procedure is defined by a constant, α, which is the number of contacts that an infected node would on average infect in a single timestep, if all its neighbours were S. An infected individual will become R with probability β. The values of α and β allow us to understand how fast the infection is being transmitted, and define the quantity R 0,f ree = α β , which will be constant in our results and equal to 6.2, same as in [1] . The effective R 0 in SIRX however is smaller, defined as R 0,ef f = α β+κ+κ 0 [1] .

These are the main processes that take place in most epidemiological models. However, it is also important to consider how lockdown measures would affect the dynamics of the system. Behavioural changes, like wearing masks and hand washing can also slow down the spread of the infection in the case of COVID-19. It is assumed that due to these changes S and I individuals will be removed from the transmission process with a rate κ 0 . Due to this fact, an S individual may become R with rate κ 0 and will never get infected, while (not yet symptomatic) I individuals can become quarantined (X) with the same rate κ 0 , as measures like social distancing are assumed to have the same impact to the population of both S and I. Subsequently, a lot of I individuals become at some point symptomatic, which leads them to either be hospitalized or self-quarantined. To explain this procedure, it is assumed in the model that I individuals will become X with a rate of κ. Thus, the total rate that I individuals become X is κ + κ 0 .

The analytical model can eventually be described by the following figure (Fig. 1 ) and a set of differential equations (as found in [1] ):

S I X R α β κ + κ 0 I 0 f actor κ 0 Figure 1 : Schematic of the SIRX model. The parameters α, β, κ 0 , and κ are the rate constants for the transitions described in the text. I 0 f actor = I X at the initial timestep.

During the time associated with a single step, which equals to a single day of data, all possible transmissions between any 2 states of the model, as it was described above, will be treated as possible event outcomes. The occurrence of any of those will be tested by creating random numbers following a Monte Carlo procedure to implement the simulation model. i) An infection coming from an I node and targeted towards a nearest neighbour which is in state S, turning it into I, with a certain probability τ , which is given by:

in which equation τ represents the transmissibility, and < k > the average degree of the network. ii) an S becoming an R, when removed because of social distancing, lockdown or behavioural changes. S nodes become R with probability κ 0 , iii) an I becoming X for the same reasons. I nodes become X with the same probability κ 0 , iv) an I becoming X after exhibiting symptoms of the disease and becoming quarantined or hospitalized. I nodes become X with probability κ, and v) an I becoming R when recovering from the virus, with probability β.

We note that 4 out of 5 of these event types (all except for the first one) remove nodes from the transmission process. As the spreading process starts at time t=1, we first process the S-I infection (step 1 above), and subsequently we proceed with the other four (4) events in random order.

Time increases by a single step after evaluating and processing all possible events for all active nodes once. This procedure is repeated in the same fashion for a certain length of time. To account for noise from the stochastic phenomena imposed by these processes, we average over many independent realizations the values of the total number of S, I, R, and X nodes in each timestep.

Finally, we create a network on top of which the above processes occur, and allow for an SIRX spreading process to occur on it. The network nodes represent the individuals and the network links represent their contacts for a given region (or country). The type of the network defines the way the contacts between nodes are formed. We test several different kinds of networks, including Erdős-Rényi random networks, Watts-Strogatz Small World networks, and Barabasi Albert, preferential attachment Scale-Free networks [23] .

Given the fact that our results fit better on Erdős-Rényi networks, we focused on studying the characteristics that would best fit. In random networks, which have been used in many such studies [24] , every possible link between two nodes has the same probability to be formed. In such networks, the relation between the average degree < k > of the network and the probability p of the creation of a link is:

where N is the total population of the nodes. Thus, by setting the < k > value of the network, p can be computed and the network contacts are easily formed. The contacts will remain the same throughout the entire simulation, as we assume that the nodes have the same contacts each day.

We show here the results of the simulation model along with those of the analytical model and the real data. Following the work of [1] we start by solving the differential equations (Eq. 1), and vary the values of the parameters (κ, κ 0 , I 0 f actor) until we find the best fit of X to the real data. As in the analytical model, we keep α and β constant (α = 0.775, β = 0.125). To accommodate for J o u r n a l P r e -p r o o f Journal Pre-proof the fact that the differential equations produce in some cases unrealistically small κ + κ 0 values, which themselves lead to a need for unrealistically high < k > values for the networks, we set a minimum to the parameter κ. We chose that minimum to be κ = 0.2, since it is a plausible assumption that a rate equal to, or even slightly larger than β + κ 0 , can be used to justify that the number of asymptomatic individuals is comparable to that of the symptomatic ones, as referred on all published research so far [25] [26] [27] . Increased testing of individuals that were in close contact with confirmed cases, which was applied in some countries more fervently than in others, can also help explain why further increasing the κ value is a reasonable choice. The set of derived parameters from this solution is stored, so as to be used later on in the simulations. In Fig. 2 we plot the value of X as a function of time together with the actual number of confirmed cases for several provinces in China. Our results do not show any palpable difference from those in [1] after having set a minimum value for κ. Indeed, we still find a satisfactory fit for the analytical model with the real data. However, the parameters (I 0 f actor, κ, and κ 0 ) for the best fit in the real data, which are presented in Table 1 , are different from the ones shown in the original work.

Given that our simulations would be computationally very demanding for large countries (or typical Chinese regions), we use random networks of smaller size, typically of size N = 10000 nodes. By acquiring the values of the model for such N , we can directly compare the results of the model with those of well averaged simulations. The result of such a comparison is shown in Figure  3 where the regions are kept the same as before. To obtain these results we examined several different < k > values for each region, until we found the best fit with the results of the model. Thus, from this best fit we surmise which < k > values pertain for each region in China. We see that all < k > values are in the range < k >= 12 − 18. We then use the same modified analytical model with data from several different countries and we show the results in Figure 4 . The parameters (κ, κ 0 , and I 0 f actor) that best fit the model to the real data are presented in Table 2 . We get an equally good agreement of the modified analytical model with the real data of confirmed cases, as we did for China. Figure 5 shows the numerical results of the differential equations with N = 10000 together with the result of the simulation. We also observe here good agreement for data in several different countries, with < k > in the range of 5 − 15. In Figure 6 we present the number of S, I, R as a function of time for Romania and Serbia. We show the calculation of the analytical model and compare it to the calculation derived from the simulation on Erdős-Rényi networks, for the < k > presented in Figure 5 . We observe good agreement. Finally, we investigate some countries with a relatively small population (N < 0.5 million people), and compare with their confirmed cases data. Regarding these countries we carry out some more extensive simulations using their actual population figures. In Figure 7 we produce the same calculations as above with the modified model for these small countries, using their entire population. The results of the model and the real data are in good agreement. The parameters (κ, κ 0 , I 0 f actor) that best fit the model to real data are presented in Table 3 . A comparison between the results of both the differential equations and the simulations is also presented in the same figure. We observe that, as before, both methods agree well, with < k > values of the random networks used in the simulations being in the range of 13 − 26.

Apart from the random networks that we used for the calculations we also investigated the behavior of other network structures, in order to find their best fit values to the real data. Next, we show the results for small-world and scale-free networks in the number of confirmed cases, as a function of time.

The results for a small-world network are shown in Figure 8a where we observe lower values of infected individuals even for values of < k > significantly greater than the ones produced by the best fit of an Erdős-Rényi network. It should be noted that small world networks can have different characteristics according to their rewiring probability, p, which is a property that defines the networks characteristics and its small-world-ness [28] . To this purpose, we have used various values for the rewiring probability. Even those with larger rewiring values (up to 0.5) produce qualitatively similar results. Thus, we reach the conclusion that a simulation with a small-world network would result in getting a much smaller disease outbreak than that expected by other means (such as the differential equations), and the ones ultimately observed in the given country. We then used a scale-free network and carried out the same simulations, where we observe that even for very small values of < k >, we eventually obtain much larger values of confirmed cases (Figure 8b) . Thus, even with < k > values as small as 2, we obtain values larger than the real number of confirmed cases. As we observe from the figure, the presence of hubs results to an unrealistically high spread of the disease, which does not agree with the real data Table 3 : Parameter values of κ, κ 0 , I 0 f actor of the modified model for the 4 countries that correspond to those presented in Figure 7 Country 

In the present work we developed a simulation model for the spread of COVID-19 that is based on the analytical approach of Maier and Brockmann [1] . We varied slightly the original model, to make it more realistic, and added a lower threshold for the values of the parameter κ, the rate by which an I becomes X when being symptomatic. An additional reason why we believe that this threshold is needed is that soon after the initial outbreak most countries instigated a policy where confirmed contacts of infected individuals were advised to get tested. This policy increased the number of infected individuals that were confirmed as being infected, even if they exhibited no obvious symptoms. We use a Monte-Carlo simulation model and find that it agrees nicely with the results of the differential equations of the original work. While the previous work was applied only to regions of China, we applied our simulation model to several other countries besides China. Simulations for other large countries were carried out on networks of size N = 10000 nodes, and the results agreed nicely with the real data. Finally, we investigated some small countries for which it was possible to carry out simulations, even for their entire population (smaller than half a million). We found no signs in theses systems that our simulations would not be applicable for even larger systems (larger countries), although at a high computational cost.

The simulations were carried out using three common network topologies, i.e. Erdős-Rényi, small-world and scale-free networks. The structure appears to greatly affect the results of the disease spreading process. When trying to simulate a model of SIRX type, we found that only the Erdős-Rényi results agreed well with the real data, while the other two topologies showed considerable deviations from the real data. Our conclusion agrees nicely with [29] where the results of the implemented model fit real data for a random type of network. Thus, despite our expectations, networks where "hub" type superspreaders exist (scale free) or are in theory more J o u r n a l P r e -p r o o f Journal Pre-proof descriptive of a societal structure (small-world), did not produce convincingly accurate results.

The Erdős-Rényi networks used were probed for several values of < k >, and we kept only those that fit best with the real data. We concluded that the COVID-19 outbreak progression, in the regions and countries examined, can be well described by a spreading process on a random network whose peers are connected on average with 10-25 other people. Values of < k > outside this range do not agree well with the real data. Thus, by performing simulations on a given country and finding the best < k > for its confirmed cases, we could predict the spreading and progression of this virus An interesting question is how far out in the future can one predict the evolution of the disease. Simply put, this is possible only in the short run, when the parameters used do not change much in their values, in the same spirit as was done and explained in [1] . Once other external factors enter the system, such as changes in quarantine measures, social distancing measures etc. then one would need considerable more information and parameters to predict the disease evolution.

It is possible to extent this work and include vaccination strategies in our work, thus, account in greater detail for the decreasing population of susceptible people in some cases in a better way. In addition to this, the model presented in the original paper seems to be unable to account on its own for the loosening of quarantine measures and the re-population of susceptible individuals with individuals that underwent quarantine without having been infected and no longer are in quarantine. A re-adjustment of the differential equations presented in Eq. 1, so as to account for these two events, would allow for a better representation of the phenomenon and more accurate simulation results.

Effective containment explains sub-exponential growth in confirmed cases of recent COVID-19 outbreak in Mainland China. medRxiv

Computer viruses and malware

Infectious Diseases of humans: Dynamics and control

A critical review of the impacts of COVID-19 on the global economy and ecosystems and opportunities for circular economy strategies. Resources, Conservation and Recycling

Piet Van Mieghem, and Alessandro Vespignani. Epidemic processes in complex networks

The mathematical theory of infectious diseases and its applications

The mathematics of infectious diseases

Forecast and control of epidemics in a globalized world

Explicit formulae for the peak time of an epidemic from the sir model. which approximant to use?

Explicit formulae for the peak time of an epidemic from the sir model

A reaction-diffusion SIS epidemic model in a time-periodic environment

Distribution of infected mass in disease spreading in scale-free networks

Quantifying the effect of quarantine control in Covid-19 infectious spread using machine learning. medRxiv

Quarantine-generated phase transition in epidemic spreading

Social distancing measures to control the COVID-19 pandemic: Potential impacts and challenges in Brazil

Virological and clinical cure in COVID-19 patients treated with hydroxychloroquine: A systematic review and meta-analysis

A strategic approach to covid-19 vaccine r&d

A Contribution to the Mathematical Theory of Epidemics

Vaccination strategies against COVID-19 and the diffusion of anti-vaccination views

An interactive web-based dashboard to track covid-19 in real time

Emergence of scaling in random networks

Random graph models of social networks

Universal testing for covid-19 in essential orthopaedic surgery reveals a high percentage of asymptomatic infections

Comparison of transmissibility of coronavirus between symptomatic and asymptomatic patients: Reanalysis of the ningbo COVID-19 data

Estimating the asymptomatic proportion of coronavirus disease 2019 (COVID-19) cases on board the Diamond Princess cruise ship

Collective dynamics of 'small-world' networks

Spreading of infections on random graphs: A percolation-type model for covid-19

Results presented in this work have been produced using the Aristotle University of Thessaloniki (AUTh) High Performance Computing Infrastructure and Resources. Funding: This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.