In this paper we study the efect of targeted immunization on the peak number of infections in an epidemic outbreak. For this we extend a previously developed python-based dashboard environment for the time eficient simulation-based study of SIR epidemics spread on complex network topologies, using realistic Continuous Time Markov Chain (CTMC) simulations by means of Gillespie's stochastic simulation algorithm. The new components make it possible to study targeted immunization by means of state-of-the-art methods and to visualize typical paths of infection during the temporal evolution of an epidemic. We show results obtained with diferent centrality measures (eigenvalue centrality of the adjacency and non-backtracking matrix, degree centrality, and average path length centrality), used in targeted immunization. In the results we focus on studying the peak number of infections (PNI). The PNI is very relevant when it comes to the practical management of an epidemic, as it determines, for instance, the number of intensive care units that are needed to ofer an appropriate treatment of the disease in critical cases. However, the PNI has received much less attention in studies than the epidemic threshold, which is more relevant in the early stage of an epidemic. Our example study on classical network topologies reveal that the choice of the centrality measure for targeted immunization as well as the number of targeted nodes will have a strong impact on the drop of the PNI. Our simulation-based results on scale-free Barabasi-Albert networks show that the PNI reduction that can be achieved by using modern centrality metrics such as the non-backtracking eigenvalue drop, can lead to up to 40% lower peaks than those achieved with naïve methods such as degree based immunization (immunization of the biggest node(s)) in case of immunization of 2.5% nodes. These results underpins the crucial importance of the correct choice of the centrality metric in targeted immunization.
Efective resistance against virus disease spreading could be organized by taking into account all
or at least most of the influenced conditions. The level of such influence should be estimated by
modeling and simulation as realistically and eficiently as possible. Of particular interest is the
influence of the topology of the contact network and the impact of targeted immunization on the
dynamics of an epidemic. In the COVID-19 pandemic, it became very clear that immunization
should not just target the total spread of a virus in the population, or the epidemic threshold,
but also seek to decrease peak numbers of infections (PNI)[
However, in previous studies, the reduction of the epidemic threshold was mainly the focus
of studies. The epidemic threshold should be as high as possible and it determines whether an
epidemic is likely to invade the network, which is the case when the contagiousness exceeds
the epidemic threshold, or, otherwise, tends to disappear before invading the network. Various
measures have been suggested in the literature to increase the epidemic threshold by means
of targeted immunization. As was argued in [
Despite the big interest of the research community in methods that target the increase of the
epidemic threshold, there are relatively few studies that focus directly on the peak number of
infections [
The PNI for such moderate size networks with a power law degree distribution can, however,
be studied with Gillespie’s stochastic simulation algorithm (SSA) which was originally designed
for chemical reaction systems with a moderate, finite number of molecules [
To simulate the event’s PNI, in previous work we developed an eficient method of stochastic
simulation of epidemic outbreaks adopting Gillespie’s stochastic simulation algorithm and
exploiting specifics of the epidemic outbreak to avoid computing the full transition matrix.
Whereas the dimension of the full transition matrix (or generator matrix) of the CTMC for
the state of infections in complex networks grows exponentially with the number of nodes,
it sufices to focus on non-zero elements in this matrix, the number of which is only linearly
depending on the number of nodes in the network. This dashboard can simulate networks
with several thousand nodes as outlined in [
Next, we will briefly outline the details of the simulation method before proceeding with the topic of targeted immunization.
There are three main epidemiological models: SI, SIS and SIR. The models allow every node of a network to be in one of the following states respectively: - S - susceptible; - I - infected. Such node has symptoms of a disease and is able to infect other nodes; - R - removed (immunized). An immunized node is not a disease carrier.
In addition, our system allows taking into account simply infected(people who can be treated at home or who have no noticeable symptoms) and critically infected(people who need medical assistance). There are a number of other epidemiological models that take into account the incubation period, detected and undetected infected nodes, the possibility to return from immunized (removed) state to susceptible, etc. The main problem of such models to accurately simulate their impact on complex networks is lack of data about a disease. It is assumed here that the study of SIS and SIR dynamics can already reveal major trends in the epidemic dynamics efectiveness of targeted immunization. It remains, however, an important topic for future research whether or not the inclusion of further details of the disease progression have a major impact on the resulting dynamics.
Furthermore, in this study we show results on diseases that, after a successful recovery, result in the immunity of the infected individuals. Although in some diseases such as the COVID-19 outbreak (that was the motivating example for our study), this is debated among experts, we assume that for networks of moderate size the immunity to the disease after recovering is long enough to assume that nodes do not get infected again until the end of the epidemic outbreak. Therefore the SIR model is more suitable for this study than SIS.
For modeling the epidemic process we use Continuous Time Markov Chains (CTMCs). The
advantage of using CTMS over Markov Chains with discrete time is that CTMC allow realistic
modeling of time. Moreover, they are not based on asymptotical simplifications and stability of
mean values as do the classical epidemiological models based on diferential equations [
In a straightforward implementation of the CTMC model of the contact process all rates
at which transitions occur between network states would be described in a generator matrix
. Here the state of the entire population forms a row-column of that matrix. There are 2
possible population states ( nodes that can be infected (0) or non-infected/susceptible (1)),
leading to a 2 × 2 matrix. where denotes the number of nodes in the contact network.
A matrix component ∈ R , ̸= denotes the transtition rates., and = − ∑︀̸= is
the rate of leaving state (on the diagonal). To avoid exponential state-space, we proposed in
[
0, if node is not susceptible.
Here is the contagiousness of the virus, () denotes the number of infected neighbors of the
vertex in state . In the theory of the stochastic simulation of CTMC [
Next to the time duration between two events, we need to also simulate the type of event
that occurs in the next transition, which is decoupled from the duration[
Targeted immunization strategies select nodes that are to be quarantined, immunized or removed from the network in order to achieve a maximum desirable efect on the dynamics of the epidemic, in our case the maximum reduction of the peak number of infections.
The most common method of immunization is the immunization of nodes with the highest degree of a network (). The degree of a node is considered to be the number of its connections (edges) with other nodes. In the case of a directed graph, the degree of a node is divided into input and output , which can also be used for a certain type of contacts of network individuals.
Using this method of targeting, it is possible to select nodes with the highest degrees, or apply iterative selection with the calculation of the highest degree after removing each of the nodes, since the degrees of the nodes change after each immunization of a node.
Relatively early, and based on diferential equation models, it was argued that the number of contacts of a node in the network (degree) is not a reliable indicator for selecting nodes for targeted immunization[12]. Various alternative centrality measures for node-selection were subsequently suggested:
The largest eigenvalue of a network is a measure that originated from the spectral theory of systems where it is also known as the spectral radius of a graph. From the point of view of the epidemiology, the eigenvalue characterizes the possibility of the successful occurrence of an epidemiological outbreak in a network and can be roughly characterized at the speed of expansion of a linear system along the direction of its principal eigenvector. The relation of the dominant eigenvalue to the epidemic threshold was described in [12] and it is determined by the following formula: < = 1 , where is the virus birth rate and is the virus death rate, and is the highest eigenvalue of the adjacency matrix. It is important to pay attention to the highest eigenvalue, because the measure takes into account not only the degree of a node, but also the entire connectivity of the network, including the degree of its neighbors. Therefore, immunization by the method increases the distance between hubs and breaks the network into clusters from which it is more dificult for the disease to escape. Depending on the type of epidemic outbreak - SIS or SIR - the eigenvalue is either computed for the adjacency matrix or, respectively, for the HNB matrix.
In [12] it was shown that for SIS dynamics the maximum eigenvalue of the adjacency matrix of
the contact network is inversely proportional to the epidemic threshold [12], [
The HNB matrix was considered first by Hashimoto for message-passing systems in complex
networks and later found its application in the study of epidemic outbreaks of the SIR type
[
The HNB matrix is defined by following formula: (→),(→ℓ) = {︃1, if = and ̸= ℓ Thus the HNBM is the adjacency matrix of a network where the edges of the original network serve as nodes. A pair of edges and is adjacent (connected), if and only if the edges are adjacent that the node for which is the incoming edge has the outgoing edge , and does not point back to the node from which is the outgoing edge. In some cases the HNB matrix has no positive eigenvalues. Based on our observations this does not occur if every node of the original network has at least one incoming edge and one outgoing edge (the start node of the incoming edge is not equal end node of the outcoming edge in the case of the non-backtracking matrix). In the current software we demand the network to have positive non-backtracking matrix eigenvalues, but further analysis is needed to figure out when exactly this is the case.
The shortest path length is the number of edges that need to be traversed to get from node to node . Average Shortest Path Length (ASPL) is the average of the shortest path lengths from each node to each.
The average shortest path length turned out to be an efective measure for reducing the PNI in COVID-19 dynamics, when considering multiple-node immunization strategies. This has been shown in our previous work on the COVID-19 outbreak in Ukraine [9]. This also seems plausible, since in networks with a small average shortest path length between pairs of nodes the propagation of the virus can proceed very fast. Therefore, we include increasing the average shortest path length here as an immunization strategy.
In order to better understand the process of the spread of a virus in the network we created a dashboard, an earlier version of which was described in detail in [10]. This dashboard provides a means to run a simulation with a user-defined network and parameters of the particular epidemic, to visualize the virus propagation process over time, and plot the (averaged) results of the simulation(s). Three types of settings are distinguished: • Simulation parameters: maximum simulation time (simulation days), time(simmulational) step for collecting data from a simulation, and the number of repeated simulations with given parameters. Repeated simulations are required to test for significant efects and understand the average, best-case and worst-case behavior of the stochastic process at hand. • Network parameters: network model (complete graph, Erdős–Rényi, Watts-Strogatz, and
Barabási–Albert models), number of nodes in the network, and number of edges per node • Epidemiological parameters: epidemiological model, infection rate, period of being in infected and critically infected states. As a default data on demographic distribution and mortality statistics is the statistic on COVID-19 in Ukraine from 19 December 2020 and can be changed with replacing files by path ./scripts/datasets/.
New in this version of the dashboard are options for nodes immunization. For this, a user is given the opportunity to enter the number of nodes for immunization and one of the following methods of immunization: - Decreasing of the highest degree; - Increasing of the average shortest path length; - Maximum eigenvalue-drop of Adjacency matrix; - Maximum eigenvalue-drop of Non-backtracking matrix;
According to chosen method, the user-defined number of nodes will be immunized using a greedy algorithm.
The network is presented in the form of an interactive three-dimensional graph (Fig.2), where each of the nodes is a person and each edge is a contact. The initial color of the nodes is white and the edges are black. When a node is infected, it turns green, and if the node is critically infected, it turns red. After recovery, the node turns blue, or purple if it was critically infected. When infection occurs, the edge through which it was performed turns red. In the case of multiple simulations for given parameters, the thickness of the edge increases each time an infection is performed through it. The colors of all elements at each step are taken from the last simulation run. Also, it is possible to go to each of the steps using a slider, which is thinned under the graph.
In addition to the visualization of the network, for the last simulation, a plot is built that displays the number of nodes in each of the states (susceptible, infected, critically infected, infectious and dead) during the simulation with a given time interval. The plot, like the graph is interactive, therefore user can move, change scale and disable curves that they do not need.
For the experiments, 4 methods of immunization were chosen: decreasing the highest degree (HD), increasing the highest average shortest path length (ASPL), decreasing the highest eigenvalue of the adjacency matrix (EAM), and decreasing the highest eigenvalue of the nonbacktracking matrix (ENBM). By every method of targeted immunization, the nodes are removed in a greedy manner that is nodes are removed subsequently and after each removal, the centrality measures are recomputed.
Note that greedy removal strategies do not always lead to the best result in multiple node removal when compared to an optimization of the removed subset. However, the optimization is usually very time expensive and limits the size of the networks that can be considered. The assessment of the gap between greedy and optimization-based subset selection remains however an interesting topic for future research.
As models of networks, three classical models of complex networks were chosen: BarabásiAlbert, Watts-Strogatz (with reconnection rate of 0.3), and Erdős-Rényi model. Each network was undirected and contains 200 nodes and about 400 edges. The Erdős-Rényi model can be considered a base model in which each node has about the same number of contacts which is a common but unrealistic assumption in many previous studies on epidemic outbreaks. The Barabási-Albert network has a scale-free degree distribution and the Watts-Strogatz network features a small world structure with a small degree of separation [13].
For the infectious period, 10 simulation days were chosen. The results on the number of infected nodes on every timestep obtained and averaged from 50 simulations for each combination of the network model and method of nodes immunization on the Fig. 3 were compared with simulations in the corresponding network models without any immunization. The networks for immunization by each method were not exactly the same.
The observed in figure 3 results are quite controversial. Immunization of one node by none of the methods did not give significant results in reducing the PNI and delaying it. Only in the case of immunization by decreasing the eigenvalue of the non-backtracking matrix of ER network and also in the case of the eigenvalue of the adjacency matrix of the BA network we can observe reducing of PNI by about 10 percent. In addition, increasing the average shortest path length slightly decreased and postponed the PNI for all networks. According to facts that simulation is stochastic and for every simulation a new network was generated, there can be an error.
To understand the efect of immunization better, immunization of 5 nodes by greedy removal was performed on the same network models. The results collected on the same network models from at least 50 simulations are shown in figure 4. According to the results of the simulations, decreasing the highest eigenvalue of an adjacency matrix provides a significant decrease and postpones the PNI in the case of degree-centralized, free-scale Barabas-Albert network and only in the case. Immunization by decreasing the highest eigenvalue of a non-backtracking matrix showed a considerable efect in the case of every network. On the highly clusterized Watts-Strogatz network, the influence of immunization was the lowest. Also, we can make the conclusion that the greedy algorithm is not suitable for increasing the average shortest path length.
The parameters of networks with every method of immunization and without any are given in table 1. The parameters are averaged from at least 50 networks in every case of immunization 1 and 5 nodes by each method and without any immunization. nodes
Time to a new event(infection) generates randomly, according to the exponential distribution, the expected value of which is the number of contacts with infected nodes. The time is recalculated after every event(infection). Let’s assume that there is a moment when in a network are not many contacts with infected nodes in the network (as usual, it happens at the beginning of the simulation), in this case, the time to an event is relatively high, and might be the case, that the highest time of recovering a node is lower than . So, there will be not any infected nodes at the time of happening a new event and the simulation stops. Let’s call it failed outbreak.
Condition of failed outbreak: Total number of infected (TNI) nodes is relatively small, TNI → 0.
Fractions of failed outbreaks for the experiments with 1 and 5 immunized nodes by each method are given in table 2 with the condition TNI < 0.2* number of nodes in the network). nodes
This paper has reported progress in the simulative analysis of epidemic outbreaks following the SIR dynamics. The contributions reported in the paper were two-fold:
Firstly we provided a discussion of an extension of a dashboard for the stochastic simulation of epidemics complex networks. The extension includes the visualization of the trajectory of a single run, and the summarizing visualization for multiple runs. The latter highlights contacts in the network that are very likely to transmit the virus and the analysis of these could yield insights into some typical infection paths and how to block them.
Secondly we extended the dashboard with the possibility to assess targeted immunization strategies using state of the art centrality measures. This new feature was used to study on realistic network topologies the outbreak of a epidemic before and after targeted immunization. In contrast to existing studies that focus on the epidemic threshold, our analysis focused on the peak number of infections (PNI), which is a crucial indicator when it comes to the management of a health crisis caused by an epidemic. Interestingly, the choice of the targeted immunization strategy has a strong impact on the drop of the PNI. Naïve strategies that immunize merely by targeting "big" nodes with many contacts fall short in performance as compared to those that use Eigenvalue based immunization. The diference is up to 40% less infections when using the right immunization strategy as compared to no-immunization or degree based immunization in the averaged results in case of immunization of 2.5% nodes.
This paper is a rather theoretical study, but it has an important message for the practical research: Even if common sense might tell us that we should focus on nodes with many contacts ifrst in immunization, the study shows that a significantly better efect can be achieved on the peak number of infections by using more sophisticated centrality measures than degree, namely global centrality measures such as eigenvector based or path-based centralities that take into account the topology of the entire network. It will be an important topic for future research to validate this result also in real world studies. However, it is unethical to experiment with real situations, but instead we can make the simulation more realistic if data from real contact networks and on the virus/disease dynamics becomes available.
In this paper a shortcut is used to simulate the treatment of nodes, which does not use rates but a deterministic threshold duration after which nodes are considered to have been treated (and are removed from the network). It will be a topic of future research to make the treatment also a stochastic process governed by a rate in the CTMC model. Moreover temporal changes of the network, e.g., the establishment or disappearance of contacts in the network could provid a meaningful extension which also can be governed by using rates in the CTMC model. In all these extensions the CTMC modeling approach can be eficiently used, because of the state-space reduction that makes use of the fact that in a single iteration only a small number of state transitions occur at non-zero rates.
The source code of the system can be found by the link: https://github.com/YulianKuryliak/ EpidemicOutbreakPredictor
Michael Emmerich ackknowledges support for research on spread dynamics in complex networks in the scheme EU Horizon 2020 RISE SMA. [9] Y. Kuryliak, M. Emmerich, D. Dosyn, On the efect of complex network topology in managing epidemic outbreaks, in: 3rd International Workshop on Modern Machine Learning Technologies and Data Science Workshop, MoMLeT and DS 2021, 2021, pp. 1–15. [10] Y. Kuryliak, System of modeling the spreading of infectious diseases, Proceedings of the International Conference on Applied Informatics ICDD, Sibiu, Romania, November 4-6 (2021) 106–114. [11] A. Lipowski, D. Lipowska, Roulette-wheel selection via stochastic acceptance, Physica A: Statistical Mechanics and its Applications 391 (2012) 2193–2196. URL: https://ideas.repec. org/a/eee/phsmap/v391y2012i6p2193-2196.html. doi:10.1016/j.physa.2011.12.0. [12] D. Chakrabarti, Y. Wang, C. Wang, J. Leskovec, C. Faloutsos, Epidemic thresholds in real networks, ACM Transactions on Information and System Security (TISSEC) 10 (2008) 1–26. [13] D. J. Watts, Six degrees: The science of a connected age, WW Norton & Company, 2004.