key: cord-0104829-7qoifljf
authors: Par'e, Philip E.; Liu, Ji; Beck, Carolyn L.; Nedi'c, Angelia; Bacsar, Tamer
title: Multi-Competitive Viruses over Static and Time--Varying Networks
date: 2017-02-24
journal: nan
DOI: nan
sha: 1491b1994a59738353532cc5faba1767f289e11b
doc_id: 104829
cord_uid: 7qoifljf

Epidemic processes are used commonly for modeling and analysis of biological networks, computer networks, and human contact networks. The idea of competing viruses has been explored recently, motivated by the spread of different ideas along different social networks. Previous studies of competitive viruses have focused only on two viruses and on static graph structures. In this paper, we consider multiple competing viruses over static and dynamic graph structures, and investigate the eradication and propagation of diseases in these systems. Stability analysis for the class of models we consider is performed and an antidote control technique is proposed.

Spread dynamics have been studied for hundreds of years. Bernoulli developed one of the first known models inspired by the smallpox virus [1] . In this paper we focus exclusively on susceptible-infected-susceptible (SIS) models, which have been developed for both continuous [2] - [5] and discrete time domains [5] - [7] . SIS models consist of a number of agents that are either infected or healthy (susceptible), which may cycle (aperiodically) between these two states. The infection rate combined with the connectivity of the ith agent with infected neighbors j (denoted by β ij ) positively affects the probability of being infected, while the healing rate δ i negatively affects the infection probability. This is depicted in Figure 1a .

The idea of two competing SIS viruses, namely the bi-virus model, has been recently pursued in [8] - [13] . The main motivation for such systems is that of competing ideas spreading on different social networks. However these models can have broader applications to political stances, adaptation of competing products, competing practices in farming, etc. and can be generalized to more than two viruses. Consider, for example, the case of three competing viruses; then each state has four possible states: susceptible, infected with virus 1, 2, or 3. The idea of information diffusion on two layered networks has also been explored for a susceptible-infected-recovered (SIR) model in [14] .

In [15] , a different multi-virus model is considered.

Further, all previous work on competing viruses has focused on viruses over static graph structures. There are recent results for the single virus model over time-varying networks [16] - [20] . Some of the ideas from [19] , [20] will be employed in this paper and applied to a more general model. Various control techniques have been applied to SIS virus systems [13] , [21] - [23] . These techniques assume the healing rate is a control variable. In [13] , it is shown that there exists no distributed linear feedback control that can stabilize the system, and in fact, will destabilize the system. Alternative approaches focus on reducing the maximum eigenvalue of the linearized system using the healing rate and/or the infection rate. In [21] , [22] , distributed control techniques for setting healing rate and quarantine protocols are proposed and implemented on a severe acute respiratory syndrome (SARS) simulation model. In [23] , a bound is provided for the cost of fairness of mitigating the spread of disease, that is, the difference between the optimal solution and the fair or homogeneous solution, for several classes of graphs. In [25] , geometric programming ideas are used to control single SIS virus systems and the authors present a polynomial time algorithm illustrated on an air transportation network. In [24] , similar ideas to [25] are applied to the bi-virus model.

In this paper we present a generalization of the bi-virus model to an arbitrary number, m, of competing viruses.

We provide conditions for stability of the disease-free equilibrium (DFE) for static as well as time-varying graph structures. We also provide sufficient conditions for stability of the non-disease free equilibrium (NDFE). We provide two control techniques based on minimizing the maximum eigenvalue of the linearized system, appealing to some of the theorems presented herein. These control techniques, which are different from other approaches in the literature, allow every agent to have a base healing rate and an additive control term.

The paper is organized as follows: we first introduce, in Section II, the SIS model and the competing virus model for m viruses. In Sections III and IV we analyze the model, providing conditions for stability of the DFE and the NDFE, and in Section V we provide an antidote control formulation. In Section VI, we present a set of illuminating simulations of various competing virus models over time-varying networks, and we conclude with some discussion in Section VII.

For any positive integer n, we use [n] to denote the set {1, 2, . . . , n}. We view vectors as column vectors. We use x T to denote the transpose of a vector x. The ith entry of a vector x will be denoted by x i . The ijth entry of a matrix A will be denoted by a ij and, also, by [A] ij when convenient. We use 0 and 1 to denote the vectors whose entries are all equal to 0 and 1, respectively, and I to denote the identity matrix, while the dimensions of the vectors and matrices are to be understood from the context. For any vector x ∈ IR n , we use diag(x) to denote the n × n diagonal matrix whose ith diagonal entry equals x i . For any two sets A and B, we use A \ B to denote the set of elements in A but not in B. 

where p i is the probability that agent i is infected, the β ij 's are (possibly asymmetric) infection rates incorporating the nearest-neighbor graph structure, and δ i is the healing rate. Neighbor relationships among the n agents are described by a directed graph G on n vertices with an arc from vertex j to vertex i whenever agent i can be infected by agent j. The agents can also be thought of as groups of people and p i 's as the percentages of the groups that are infected, and therefore the neighbor graph G can have self-arcs at all n vertices. Hence, β ij equals zero if there is not an edge in G from node j to node i. The model in (1) is more general because the underlying graph G can be directed and the weights given by β ij can be any non-negative number. The representation in (1) can be put into matrix form:ṗ

where p is the vector of the p i 's, B is the matrix of the β ij 's, P = diag(p), and D = diag(δ 1 , . . . , δ n ). In the analysis that follows, as stated above, B is not assumed to be symmetric unless explicitly stated so.

This model has been extended to have two viruses, providing a generalization of the model introduced in [10] ,

where p 1 i (t) and p 2 i (t) are the probabilities that agent i has virus 1 and 2 respectively, and each virus has its own infection rates and healing rates. Each virus spreads over a (possibly different) spanning subgraph of G, where their union is the neighbor graph G. It will be assumed that both of the two subgraphs are strongly connected and, thus, so is G. 1 We need not restrict ourselves to two viruses, however. A direct generalization leads to the following multi-virus

for all k ∈ [m]. This representation can be written in matrix form as:

where the matrices are the same as in (2), but now they are dependent on which virus they correspond to. Since the subgraph for each virus k is strongly connected, it follows that B k is irreducible, meaning that it cannot be permuted into block triangular matrix form. The assumption that B k is bounded means ∀i, j, β k ij < ∞. The set

is invariant with respect to the system defined by (5) . If p k i denotes the probability of agent i being infected by virus k and 1 − m k=1 p k i denotes the probability of agent i being healthy, it is natural to assume that their initial values are in [0, 1], since otherwise the values will lack any physical meaning for the epidemic model considered herein. Similarly, if the states were representative of the density of infected members of a sub-population, they would also be bounded between zero and one. 

. Consider an index i ∈ [n]. If p k i (τ ) = 0, then from (4) and the assumption that the matrices B k are non-negative,ṗ k i (τ ) ≥ 0. The same holds for p 1 i (τ ) + · · · + p m i (τ ). If p k i (τ ) = 1, then from (4) and the assumption that the matrices B k are non-negative,

For the rest of the paper we assume p k

. It has been shown that there are disease-free equilibrium and non-disease free equilibria for the single virus system [19] , [20] , [26] , [27] , as well as for the two-virus system [13] ; the same applies to multi-virus systems as well. However, in this case the scenario becomes slightly more complicated because all viruses can reach the DFE, 1 A directed graph is strongly connected if for any two distinct vertices i and j, there is a directed path from i to j. or a NDFE, or there may be some viruses at a DFE and some at a NDFE. We will explore several conditions for convergence to these different equilibria.

First, we explore stability of the DFE for both the static and dynamic graph cases.

We first give conditions under which the DFE is asymptotically stable.

, then the healthy state is the unique equilibrium of (5), which is asymptotically stable with domain of attraction D, as defined in (6).

Proof. To prove the theorem, it is sufficient to show that for all k ∈ [m], p k (t) will asymptotically converge to 0 as t → ∞ for any initial condition.

Since for all k ∈ [m], p k i (t) is always non-negative by Lemma 1, from (4),

which implies that the trajectories of p k i (t) are bounded above by a single-virus model. Since the B k 's are nonnegative and irreducible, by Proposition 3 in [13] , p k i (t) will asymptotically converge to 0 as t → ∞ for all k ∈ [m], and thus the healthy state is the unique equilibrium of (5).

We next state a result on global exponential stability for the case when the underlying subgraphs are undirected and the infection rates are symmetric.

Theorem 2. Suppose B k is symmetric, and the maximum eigenvalue of B k − D k is less than zero, that is λ 1 (B k − D k ) < 0. Then the DFE is exponentially stable for virus k, with domain of attraction D, in (6) .

The first inequality holds because (P l B k ) ij ≥ 0, ∀l, i, j by construction since each p l i (t) is a probability. The second inequality holds by the Rayleigh-Ritz Theorem because B k − D k is symmetric. Therefore, the system converges exponentially fast to the origin by Theorem 8.5 in [28] .

Note that this is a generalization of the result in [19] , [20] .

We can state that the condition in Theorem 1 is necessary and sufficient for eradication of all viruses. Proof. Sufficiency has been shown in Theorem 1. Therefore, to prove the theorem, all that needs to be shown is that if for any j ∈ [m] s(B j − D j ) > 0, the system (5) admits a NDFE.

Without loss of generality, suppose that s(B 1 − D 1 ) > 0. Set p k = 0 for all k = 2, . . . , m. Then, the dynamics of p 1 simplifies to a single-virus system, which admits a NDFE by Proposition 4 in [13] . Therefore, in the case when s(B 1 − D 1 ) > 0, the system (5) always admits an equilibrium of the form (p 1 , 0, . . . , 0) withp 1 0.

We can generalize the model from (4) to have dynamic graph structure aṡ

where β k ij (t) is a function of time and the equation holds for k = 1, . . . , m. We now provide a sufficient condition for global exponential stability of the DFE. 

Then the DFE is exponentially stable for virus k, with domain of attraction D, in (6).

The first inequality holds because (P l B k (t)) ij ≥ 0, ∀l, i, j, t by construction since each p l i (t) is a probability. The second inequality holds by the Rayleigh-Ritz Theorem because B k (t) − D k is symmetric. The last inequality holds by definition of the supremum. Therefore, the system converges exponentially fast to the origin by Theorem 8.5 in [28] .

This result is a generalization of Theorem 1 in [19] , [20] .

We can also show exponential stability for the case when the infection rates are not symmetric and the underlying subgraphs are undirected, with some added assumptions.

Definition 1. For a given virus k, assume that for all t ≥ 0, there exist c k (t), λ k (t) > 0 such that

We then define

Note that

Theorem 5. Consider the dynamics for virus k in

enough µ > 0, then the DFE is exponentially stable for virus k, with domain of attraction D, in (6).

Proof. Note that since (P l (t)B k (t)) ij ≥ 0 ∀l, i, j, by construction,

Therefore, by Grönwall's Inequality ( [28] ), the solution of the original system will be bounded above by the solution of the linear system. Thus by Lemma 2 in [20] , the DFE is exponentially stable for virus k.

Note that this theorem is a generalization of a single virus result provided in [20] , where Lemma 2 in [20] is for a less general model; however the same arguments hold by replacing BA(t) with B k (t) and BȦ(t) withḂ k (t).

Theorem 6. Consider the dynamics for virus k:

for all t 0 ≥ 0, and for some ν > 0 there exists an h > 0 such that

for all t ≥ 0 and some γ, 0 < γ ≤ 1. Assume further that

for some negative scalarᾱ and for all t 0 ≥ 0,

for all t 0 ≥ 0, and for all i, j and t ≥ 0 the perturbation

Then the origin is exponentially stable for virus k.

Proof. Since (P (t)(B k (t) + ∆ k (t))) ij ≥ 0 ∀i, j by (18) and Lemma 1,

Therefore, by Grönwall's Inequality ( [28] ), the solution of the original system will be bounded above by the solution of the linear system. Thus by Lemma 2 in [20] , the origin is exponentially stable for virus k. Note that this result is an extension of Theorem 3 in [13] , and we present a sketch of the proof.

Sketched proof of Theorem 7: From the proof of Theorem 1, p k (t) will asymptotically converge to 0 as t → ∞ for all initial values (p 1 (0), . . . p m (0)) ∈ {(p 1 , . . . , p m )|p i = 0 and p k ∈ [0, 1] n ∀k = i}, for k = i. From (5),

Thus, we can regard the dynamics of p i (t) as an autonomous systeṁ

with a vanishing perturbation − k =i P k (t)B i p i (t), which converges to 0 as t → ∞. From Proposition 5 in [13] , the autonomous system (19) will asymptotically converge to a unique epidemic state (0, . . . , 0,p i , 0, . . . , 0) for any (p 1 (0), . . . , p m (0)) ∈ D \ {(p 1 , . . . , p m )|p i = 0 and p k ∈ [0, 1] n ∀k = i}, with D defined in (6) .

Another possible NDFE is that of coexisting equilibrium, which is where more than one virus survives. We have the following interesting result similar to Theorem 7 in [13] .

Theorem 8. Consider the model in (4) with each virus propagating over the same strongly connected graph G with the corresponding adjacency matrix A, and each virus homogeneous in healing and infection rates, that is,

for some constant α ik > 0.

Proof: The homogeneity assumption on the infection rates allows us to factor B k = β k A, for each virus k ∈ [m].

To be an equilibrium of (4) the following must hold for all k ∈ [m]

in which (I −P 1 − · · · −P m )A is an irreducible Metzler matrix 2 , sinceP 1 + · · · +P m is diagonal and [P 1 + · · · +P m ] ii < 1 for all i ∈ [n], by Lemma 1. From Lemma 2 in [13] , it must be true thatp k 0 ∀k ∈ [m] and

While the stability of the time-varying case has been explored in [20] , the time-varying NDFE or epidemic limit

cycle is an open problem, even for the single virus case. Some work has been done to show the existence of a periodic NDFE for a single virus switching system in [18] .

Let us assume that for each agent, in addition to the healing rate, there is a control input u i (t) that acts as an additive boost to the healing rate. This implies that the controller can increase the agents' ability to recover from the virus, which can be thought of as the administration of an antidote or some other type of treatment. This effect is portrayed in the model aṡ

We define U (t) = diag(u) with u = [u 1 (t), . . . , u n (t)] T . To simplify the discussion in this section, we assume that B k (t) is symmetric, piecewise continuous in t, and bounded ∀t ≥ 0. Similar to the approaches in [21] - [23] , [25] , [29] , we focus on minimizing the maximum eigenvalue of B k (t) − (D k + U k (t)). Even though these control techniques are generally effective, we believe the approaches herein are more general and simpler, and therefore more scalable. Also, the assumption that our control input is additive to the base healing rate is novel and more sensible for the main motivating example, that is, every agent should have some inherent healing rate that should not be affected by the controller.

While the solutions to the following posed problems may not meet the conditions of Theorems 2 and 4, that is, they may not result in the maximum eigenvalues being less than zero, they push the system towards those conditions, consistent with the principle of the average being less than zero, presented in Theorem 6. And in practice, illustrated by simulation in the next section, these techniques reduce the spread of the epidemics. Under the aforementioned assumptions we can formulate the following optimization problem for each virus k, appealing to Theorems 2 and 4 depending on whether B k is constant or time dependent:

given that B k (t) is symmetric for all t ≥ 0.

From the Gershgorin Disc Theorem [30] it is clear that by sufficiently increasing the u k i 's, the conditions of Theorems 2 and 4 will be satisfied. Therefore we can relax the above optimization problem to obtain the following:

This is clearly a linear program and can easily be solved.

To make this a more compelling and realistic problem, we can impose a constraint on the number of agents that can be affected, which is a reasonable assumption because the cost of providing a low-dose treatment to all agents is higher than providing that same treatment dose to a few select members of the population (such as the sickest or most susceptible agents). Define the sparsity metric · 0 as the number of the non-zero entries in its argument.

Employing the sparsity metric, we have the following problem, with a capacity constraint and a sparsity constraint:

where d k is the maximum number of agents that can be treated for virus k. At first glance, the second and third constraints may seem redundant; however, the 1 constraint limits the total amount of antidote that can be used while the sparsity constraint limits the number of agents that can be treated. The inclusion of the 1 constraint prevents an infinite amount of antidote being administered to the limited number of agents allowed by the sparsity constraint.

It is well known that · 0 is highly non-convex [31] , making the above problem difficult to solve. Therefore, to solve it we employ another relaxation using the reweighted 1 norm [32] .

Definition 2. The weighted 1 norm is

where w i 's are positive and can be a constant or depend on time.

In view of this, we can rewrite the above problem as the following:

where κ is a constant weighting factor.

An effective heuristic for the selection of the weights w k i 's in (21) , proposed in [32] , is, for some small > 0,

For completeness, we include Algorithm 1, which explains the implementation of this heuristic to solve Problem 2. The notation Problem 2(w k−1 ) indicates that w k−1 is used for the weighted 1 norm in the objective function of Problem 2 in the kth iteration. Employing this heuristic yields a good solution to Problem 2 but clearly is expensive, since it requires the calculation of multiple solutions. The effectiveness of this approach is illustrated in the following section via simulation. (23) and (24) and the graph structure follows (25)- (27) . A video of this simulation can be found at youtu.be/j_MHm08dA_o.

In this section we present a set of illuminating simulations of various competing virus models over static and time-varying graph structure networks. Due to limit of dimensions in color and size, for the simulations we will only have three competing viruses. Virus 1 is depicted by the color red (r), virus 2 is depicted by the color blue (b), and virus 3 is depicted by the color green (g). For all i ∈ [n], the color at each time t for agent i is given by

When p 1 i (t) + p 2 i (t) + p 3 i (t) = 0, the color goes to black, indicating completely healthy, susceptible. These are used to facilitate the depiction of the parallel equilibrium (p 1 = α 2p2 = α 3p3 ), which will be shown by all nodes converging to the same color. For all i ∈ [n], the diameter of the node representing agent i is given by

with d 0 being the default/smallest diameter and r 0 being the scaling factor depending on the total sickness of agent i. Therefore the color indicates the sickness each agent has and the diameter indicates how sick each agent is.

For systems that have three different subgraphs, viruses 1, 2, and 3 spread on the graphs depicted by gray, green, and pink edges, respectively. If all viruses spread on the same graph, the edges are gray.

The simulation in Figure 2 has three viruses spreading over the same time-varying graph. Similar to [20] , the graph structure is determined by where z i (t) ∈ R 2 is the position of agent i, withr = 10. The agents have piece-wise constant drifts, that is,

where φ(t) ∈ R 2 and is determined, for each dimension l ∈ [2], by

where the agents hover around a square, centered at some point z c . The initial positions and φ's are chosen randomly.

Each virus is homogeneous in infection rate. The first two viruses meet the assumptions of Theorem 4, while the maximum eigenvalue of the third virus fluctuates between being positive and negative. See Figure 3 for a plot of the maximum eigenvalues of the three-virus dynamics. Consistent with the theorem, the first two viruses are eradicated quite quickly. The third virus is also eliminated, but it takes a little longer. This eradication is illustrated in Figure   2b .

The simulation in Figure The simulation shown in Figure 5 meets the assumptions of Theorem 8, that is, the three viruses are each homogeneous, with δ 1 β 1 = δ 2 β 2 = δ 3 β 3 , and propagate over the same graph structure. There are 15 agents and the initial conditions are given in Figure 5a . Consistent with the theorem, the system converges to a co-existing parallel equilibria.

We conclude with a simulation that implements the control techniques presented in Section V. Consider the single virus system in Figure 6 . This system is homogeneous in infection rate, with β = 0.492. We compare the system with no controller (on the left), a controller using Problem 1 (in the middle), and a controller that uses Algorithm 1 to solve Problem 2 iteratively with κ = .05 (on the right). The sum of the final probabilities of infection for all (23) and (24) . A video of this simulation can be found at youtu.be/zCRiLr8sWEM. (23) and (24) . A video of this simulation can be found at youtu.be/sy_RoUP7qUs. agents ( n i=1 p i (100)) for the three plots are 10.7, 4.92, and 3.6, respectively. Therefore, Algorithm 1 performed the best, however both had significant improvements over the uncontrolled simulation. The maximum eigenvalues of the three linearized systems are, from left to right, 1.893, 0.557, and 0.421; so none of the linearized systems are Hurwitz. Therefore, consistent with Theorem 3, the systems are all at NDFE. However, even though the control efforts do not completely eradicate the virus, they do mitigate its effect.

We have explored the competing multi-virus SIS model with several theorems exploring stability of the equilibria of the model for the static and time-varying graph cases. We have also proposed several control techniques that In future work we would like to explore more generic cases of co-existing epidemic states. Further, we would like to compare the techniques in Section V to other existing techniques. We would also like to implement the control techniques on large scale systems with at least tens of thousands of nodes.

Essai d'une nouvelle analyse de la mortalité causée par la petite vérole et des avantages de l'inoculation pour la prévenir

Contributions to the mathematical theory of epidemics. ii. the problem of endemicity

Epidemiological models and lyapunov functions

Virus spread in networks

Global dynamics of epidemic spread over complex networks

Epidemic spreading in real networks: An eigenvalue viewpoint

Epidemic threshold and immunization on generalized networks

Winner takes all: Competing viruses or ideas on fair-play networks

Competing memes propagation on networks: A network science perspective

Competitive epidemic spreading over arbitrary multilayer networks

Bi-virus sis epidemics over networks: Qualitative analysis

On the analysis of a continuous-time bi-virus model

On the analysis of a continuous-time bi-virus model

Conjoining speeds up information diffusion in overlaying social-physical networks

A stochastic model of multivirus dynamics

Virus propagation on time-varying networks: Theory and immunization algorithms

Spread of epidemics in time-dependent networks

Stability criteria for SIS epidemiological models under switching policies

Stability analysis and control of virus spread over time-varying networks

Epidemic processes over time-varying networks

Network design problems for controlling virus spread

Designing spatially heterogeneous strategies for control of virus spread

Cost of fairness in disease spread control

Optimal resource allocation for competitive spreading processes on bilayer networks

Optimal resource allocation for network protection against spreading processes

In-homogeneous virus spread in networks

Stability properties of infected networks with low curing rates

Nonlinear systems

Controllability metrics, limitations and algorithms for complex networks

Matrix analysis

Robust principal component analysis: Exact recovery of corrupted low-rank matrices via convex optimization

Enhancing sparsity by reweighted 1 minimization