key: cord-0637520-qdz5rsdy authors: Zhou, Yunxiu; Levin, Simon A.; Leonard, Naomi E. title: Active Control and Sustained Oscillations in actSIS Epidemic Dynamics date: 2020-07-02 journal: nan DOI: nan sha: 514621e124f2c37c5cf5f6f72cbcb81c8f623af8 doc_id: 637520 cord_uid: qdz5rsdy An actively controlled Susceptible-Infected-Susceptible (actSIS) contagion model is presented for studying epidemic dynamics with continuous-time feedback control of infection rates. Our work is inspired by the observation that epidemics can be controlled through decentralized disease-control strategies such as quarantining, sheltering in place, social distancing, etc., where individuals actively modify their contact rates with others in response to observations of infection levels in the population. Accounting for a time lag in observations and categorizing individuals into distinct sub-populations based on their risk profiles, we show that the actSIS model manifests qualitatively different features as compared with the SIS model. In a homogeneous population of risk-averters, the endemic equilibrium is always reduced, although the transient infection level can exhibit overshoot or undershoot. In a homogeneous population of risk-tolerating individuals, the system exhibits bistability, which can also lead to reduced infection. For a heterogeneous population comprised of risk-tolerators and risk-averters, we prove conditions on model parameters for the existence of a Hopf bifurcation and sustained oscillations in the infected population. Deterministic compartmental models in epidemiology have provided valuable insights for the understanding of evolutionary dynamics of infectious disease spread in a host population (Kermack and McKendrick (1927) ). These models have also been widely applied to study various other spreading dynamics, including but not limited to information dissemination in social networks (Jin et al. (2013) ), sentiment contagion in human society (Zhao et al. (2014) ) and the propagation of systemic risks in financial market (Demiris et al. (2014) ). Although these deterministic models unavoidably ignore some important details such as individual heterogeneity and network structure, they adequately capture the qualitative features of the spreading dynamics, including transient system behavior and stability of solutions. The models have been shown to be close approximations of certain Markov chain models of the underlying stochastic dynamics, see, e.g., Sahneh et al. (2013) . The Susceptible-Infected-Susceptible (SIS) model has been widely studied and applied in epidemiological modeling. While its assumption that individuals acquire no immunity after recovery may not be suitable for certain diseases, it provides a worst case scenario, which is valuable for a large class of contagious diseases in general. The SIS model in its simplest form assumes a constant infection rate. However, it is well acknowledged that the rate can vary ⋆ The research was supported in part by Army Research Office grant W911NF-18-1-0325, Office of Naval Research grant N00014-19-1-2556, National Science Foundation grants CMMI-1635056, CNS-2027908, CCF-1917819, and C3.ai Digital Transformation Institute. over time due to the different control strategies taken by individuals, which in turn affects the contagion dynamics. Indeed, understanding such models and their extensions is fundamental to developing disease control strategies. A review of analysis and control of epidemics is provided in Nowzari et al. (2016) . In this paper, we propose a model with feedback controlled infection rates to account for active control strategies. Individuals modify their contact rates with others based on what they observe about the level of infection in the population. We model a time lag in the observations and distinguish individuals as riskaverters, risk-tolerators, and risk-ignorers. Risk-averters represent those who change their contact rate with others in the opposite direction as the change in the observed infected population level, e.g., those who could and did stay at home and practiced increased social distancing during the COVID-19 pandemic as they saw the infected population grow. Risk-tolerators represent those whose contact rate with others changes in the same direction as the change in the observed infected population level, e.g., health care workers, delivery workers, and other essential workers, who were obliged to work during the COVID-19 pandemic. Risk-ignorers represent those who do not actively modify their contact rates. Our work contributes to the literature in the following ways. First, while active and passive spreading was first distinguished in the social economics literature (see Hartmann et al. (2008) ), our model rigorously demonstrates the differences between them, and we prove new results on the dynamics of contagion with active control. Further, our model serves as one form of the state-dependent approaches discussed in Rands (2010) that offer evolutionarily grounded ways for studying social contagion in collective processes. Second, our feedback model uses a low-pass filter of the measured infected population level. This models the observation delay and introduces an important robustness to uncertainty. This is in contrast to contagion models which feed back infected population level directly, as in Baker (2020); Franco (2020) . Third, we prove a new and relatively simple scenario under which sustained oscillations appear within epidemiological frameworks without external forcing. Understanding mechanisms that can lead to oscillations is critically important in the context of infectious disease spread (Lin et al. (1999) ; Dushoff et al. (2004); Camacho et al. (2011 ), Xu et al. (2020 ). It is likewise of great interest in many other socioeconomic processes, e,g., the rise and fall of business cycles (Mishchenko (2014) ) and fluctuation of behavioral preferences in social networks (Pais et al. (2012) ). Fourth, in contrast to control strategies that require global knowledge about the exact underlying spreading dynamics, e.g., Nowzari et al. (2016) , our work provides evidence for the promise of tunable decentralized active control strategies to manage the dynamics of epidemics. The paper is organized as follows. In Section 2, we review the Susceptible-Infected-Susceptible (SIS) model. We introduce the actively controlled Susceptible-Infected-Susceptible (act-SIS) model in homogeneous populations of risk-tolerators and risk-averters respectively in Section 3, where we analyze the equilibrium solutions and stability conditions and show their qualitatively different features as compared with the SIS model. We examine the actSIS model in a heterogeneous population comprised of risk-tolerators and risk-averters in Section 4, and prove conditions for a Hopf bifurcation with a stable limit cycle. We conclude in Section 5. Consider a disease spreading in a large, randomly-mixed population, where individuals are divided into either susceptible (S) or infected (I) classes. Susceptible individuals get infected at rate β while infected individuals recover at rate δ. Let p(t) be the fraction of infected individuals and s(t) the fraction of susceptible individuals at time t. The SIS model iṡ The force of infection term βp in (1) describes the rate at which susceptible individuals get infected (Muench (1934) ). It can be decomposed into three terms: βp =β × α × p, whereβ is the transmission rate of the disease and an intrinsic property of the disease, and α is the effective number of contacts per unit time. Since s(t) = 1 − p(t), (1) can be rewritten aṡ The steady-state behavior of solutions to the well-mixed SIS model (2) is characterized by the basic reproduction number R 0 = β/δ, a key concept in epidemiology that defines the epidemic threshold of a particular infection (Diekmann et al. (1990) ). If R 0 > 1, the disease persists and a nonzero fraction of the population is infected at steady state: as t → ∞, p(t) → 1 − δ/β, the endemic equilibrium (EE). If R 0 ≤ 1, the disease dies out at steady state: as t → ∞, p(t) → 0, the infection free equilibrium (IFE). At R 0 = 1, there is a transcritical bifurcation. A natural extension of the homogeneous population setting is the introduction of heterogeneities of various kinds. Characterizing population heterogeneity in terms of infection and/or recovery rates has been discussed both in population subgroups (Anderson and May (1992) ) and in networks (Hethcote and Yorke (1984) , Pagliara and Leonard (2020) ). In the following, we review the network SIS model that was originally introduced as the multi-group SIS model in Lajmanovich and Yorke (1976) . Consider a heterogeneous population of n sub-populations, each large, well-mixed and homogeneous. Let susceptible individuals in sub-population i get infected through contact with infected individuals in sub-population j at rate β ij ≥ 0 and infected individuals in sub-population i recover at rate δ i ≥ 0. The rate β ij can be decomposed as β ij =β × α ij whereβ is the transmission rate and α ij represents the effective contact rate between sub-population i and j. The network SIS model iṡ where p i (t) denotes the fraction of infected individuals in the ith sub-population, or equivalently, the probability that a typical individual in sub-population i is infected at time t. Let B = {β jk } and Γ = diag (δ 1 , . . . , δ N ) be the infection matrix and the recovery matrix, respectively. For the network SIS model (3), the basic reproduction number Our proposal of the actively controlled Susceptible-Infected-Susceptible (actSIS) model is based on the following observations. First, individuals conduct behavioral changes as they acquire information of the epidemics which consequently affect their contact rates with others (Funk et al. (2009) ). Such information, however, often involves estimations that are delayed and unavoidably omits details in the finer time scale. We are motivated in part by the model and study of change in susceptibility after first infection as presented Pagliara and Leonard (2020) . Accordingly, we present the actSIS model for studying epidemic dynamics with continuous-time feedback control of infection rates. We let the infection rate be the product of the intrinsic transmission rateβ and an effective contact rate α(·) that is actively modified by individuals based on their observations of the system state. To account for the uncertainty and delay in measurements of the infection level in the population, we let the feedback responses depend on the filtered state p s of the infected fraction p, where p s tracks p and possibly some external stimulus r(t) with a time constant τ s . The actSIS model iṡ Acknowledging that people conduct social-behavioral changes in a soft-threshold manner (Smaldino et al. (2018) ), we consider sigmoidal-shaped functions for the feedback response α(·). Similar feedback mechanisms in neuronal dynamics have been shown to exhibit ultra-sensitivity and robustness to inputs and variability (Sepulchre et al. (2019) ). Let φ : [0, 1] → [0, 1] be a monotonically increasing saturating function with location parameter µ ∈ (0, 1) and slope parameter ν ∈ (0, 1). 1 For risk-tolerators we define α to vary directly with p s : It is not surprising to find that the EE of the actSIS model for both risk-tolerators and risk-averters is upper bounded by that of the SIS model, since the incorporation of feedback responses α(·) ∈ [0, 1] always decreases the effective infection rates. However, the underlying structure resulting in these lower endemic solutions differs among types of individuals, and the actSIS model shows qualitatively different dynamical features as compared to the SIS model. For homogeneous risk-tolerators, (4) undergoes a saddle node bifurcation and exhibits bistability as illustrated in Fig. 1 and Fig. 2(a) . The saddle node bifurcation point is greater than the transcritical bifurcation point in the SIS model, which has implications for control design since it implies that it is more difficult for the disease to spread in a population of risk-averters than in a population of risk-ignorers. For homogeneous risk-averters, (4) undergoes a transcritical bifurcation, as in the SIS model, but with a reduced EE, as illustrated in Fig. 2(b) . Further, the EE becomes a stable focus under certain conditions, resulting in large overshoot and/or undershoot in the transient dynamics, as illustrated in Fig. 3 . We devote the rest of this section to the detailed description and proof of these rich dynamics. For simplicity of exposition, we set δ = 1, τ s = 10 throughout the rest of this paper. 2003)). with r(t) = 0 and α(·) = φ T (·). Then the following hold: (i) There are two types of equilibrium solutions: (a) IFE: p = p s = 0. It is always stable; (b) EE: p = p s = p * satisfying φ T (p * )(1 − p * ) = δ β . There are zero, one, or two such solutions. A solution is stable if Φ T ′ (p * ) · (1 − p * ) 2 < δ β . (ii) The system undergoes a saddle node bifurcation at the bifurcation point (β c , p * ) with stable upper branch (stable EE) and unstable lower branch (unstable EE). The epidemic thresholdβ c /δ is always larger than that of the SIS model, i.e.,β c /δ > 1. (iii) Forβ >β c the system exhibits bistability of the IFE and the EE. Forβ ≤β c the IFE is the only stable equilibrium. (iv) Asβ c <β → ∞, p * → 1 − δ/β, the EE of SIS. Proof. (i) The equilibrium solutions are straightforward to compute. For stability, we compute the Jacobian as For the IFE, (7) reduces to implying that the IFE is always stable. For the EE, Stability requires the equilibrium solution to satisfy δ 1−p * − β(1 − p * )φ T ′ (p * ) > 0. Since p * ∈ (0, 1) , it is equivalent to requiring that φ T ′ (p * ) · (1 − p * ) 2 < δ β . (ii) We start by computing the critical valueβ c . Since φ T (·) is a monotonically increasing function taking values between 0 and 1, g(p) := φ T (p)(1 − p) takes value between 0 and 1 and it first increases from 0 (since g(0) = 0) and then decreases to 0 (g(1) = 0). Depending on parameter values (µ T , ν T ,β, δ), the EE has either zero, one or two solutions. Let h(p) := φ T ′ (p) · (1 − p) 2 . It first increases and then decreases for p ∈ [0, 1]. Denotingp := argmax p g(p), one can check that h(p) intersects with g(p) at three points: 0,p and 1. This implies that when the EE has two solutions (when δ/β < g(p)), the smaller solution is unstable and the larger one is stable.β c can be solved analytically by solvingβ c = δ/g(p). To prove the existence of a saddle-node bifurcation, we use the classification of equilibria for a two-dimensional system presented in Section 4.2.5 in Izhikevich (2007) . Specifically, we show that one of the eigenvalues of the Jacobian at (β c , p * ) becomes zero. Using equationsβ c = δ/g(p) and g ′ (p) = 0, one can easily verify the determinant of the Jacobian at (β c , p * ) equals zero thus the existence of a saddle node bifurcation. Since g(p) ∈ (0, 1), we haveβ c /δ is always larger than the epidemic threshold in the SIS model. (iii) This follows directly from (i) and (ii). (iv) Comparing the equations satisfied by the EE solutions for the SIS and the homogeneous actSIS (risk-tolerators) model, we observe that φ T (p * ) and p * increase withβ. As φ T (p * ) approaches 1, p * approaches 1 − δ/β. Theorem 3.2 (Steady-state behavior of the actSIS model with homogeneous risk-averters). Consider the actSIS dynamics for a homogeneous population of risk-averters given by (4) with r(t) = 0 and α(·) = 1 − φ A (·). Then the following hold: (i) There are two equilibrium solutions: (a) IFE: p s = p = 0. It is stable ifβ < δ; (iii) The EE is always upper bounded by 1 − δ/β, the SIS EE. Proof. (i) Since 1 − Φ A (·) is monotone decreasing, there exists at most one EE solution when δ/β < 1. The Jacobian of (4) is which simplifies to for the IFE, and for the EE. Therefore, the IFE is stable whenβ < δ. For the EE, stability requires δp * 1−p * +β(1 − p * )p * Φ A′ (p * ) > 0. Since Φ A ′ (p * ) is always non-negative, the stability conditions always holds when such a solution exists. (ii) We prove when EE is a stable focus, by proving when (12) has a pair of complex-conjugate eigenvalues with negative real part. We use the classification of equilibria for a twodimensional system according to the trace and determinant of the Jacobian (Izhikevich (2007)). Denoting A := J| EE and using the definitions of a, b and c, we have Tr(A) = −(a + c) and det(A) = ac(1 + b). Therefore the condition b > (a−c) 2 4ac , guarantees that Tr(A) 2 − 4 det(A) < 0, and together with Tr(A) < 0, that the EE is a stable focus. (iii) Asβ increases, the right hand side of (1 − Φ A (p * ))(1 − p * ) = δ/β decreases. As a result, p * increases withβ. On the other hand, 1 − φ A (p * ) decreases away from 1 thus the EE are always bounded above by 1 − δ/β. We introduce the heterogeneous network actSIS model. The transmission rates, recovery rates and feedback responses can all be distinct. However, we restrict to the following set up as a first step in exploring the role of heterogeneity. Let the population comprise two homogeneous sub-populations, risktolerators and risk-averters, which differ only in their feedback responses to the infection. We assume the disease transmission occurs across sub-populations but not within. For the generalization of the network SIS model (3) to the network actSIS model, this translates as n = 2,β i =β, δ i = δ, β ii = 0, α(·) = α 1 (·) for all risk-tolerators and α(·) = α 2 (·) for all risk-averters. Let p 1 (t) (p 2 (t)) denote the fraction of risk-tolerators (risk-averters) that are infected at time t. Let β 12 =βα 1 (p s2 ) =βφ T (p s2 ) be the effective infection rate from risk-averters to risk-tolerators, and β 21 =βα 2 (p s1 ) =β 1 − φ A (p s1 ) from risk-tolerators to risk-averters. The heterogeneous network actSIS model in the case of these two subpopulations iṡ The equilibrium solutions satisfy p * 1 = p * s1 , p * 2 = p * s2 and As shown in Fig. 4(a) , there are model parameters for which system (13) undergoes a Hopf bifurcation with a stable limit cycle. Asβ increases from zero, system (13) undergoes a saddlenode bifurcation from a single stable IFE to bistability of the EE and the IFE. Asβ increases further across the Hopf bifurcation pointβ * , the system exhibits a stable limit cycle about the EE. For completeness, we present the following theorem (Theorem 3.4.2 in Guckenheimer and Holmes (1983) ), which we use to prove the existence of stable limit cycles for (13). Theorem 4.1 (Guckenheimer and Holmes) . Suppose that the heterogeneous actSIS model (13) expressed asṗ = f (p,β), p = (p 1 , p 2 , p s1 , p s2 ),β ∈ R, has an equilibrium at p * ,β * and the following properties are satisfied: • (H1) The Jacobian D p f | (p * ,β * ) has a simple pair of pure imaginary eigenvalues λ β * and λ β * and no other eigenvalues with zero real parts, Then the dynamics undergo a Hopf bifurcation at p * ,β * resulting in periodic solutions. The stability of the periodic solutions is given by the sign of the first Lyapunov coefficient of the dynamics ℓ 1 | (p * ,β * ) . If ℓ 1 < 0, then these solutions are stable limit cycles and the Hopf bifurcation is supercritical, while if ℓ 1 > 0, the periodic solutions are repelling. We show in Proposition 4.2 conditions on the model parameters that guarantee that the non-hyperbolicity condition (H1) is satisfied for the heterogeneous actSIS model (13). We verify numerically that condition (H2) is satisfied and ℓ 1 < 0. Matching coefficients of (15) and (16), we derive the expressions for a, b, and c, in terms of the model parameters (β, µ T , γ T , µ A , γ A ) for (13) to satisfy the first part of (H1). We have λ β * = √ ci, a pure imaginary eigenvalue. a = 0 guarantees there are no other eigenvalues with zero real part. A proof of conditions guaranteeing (H2) of Theorem 4.1 is the subject of ongoing work. Fig 4(b) shows numerically that (H2) is satisfied for the parameters selected. We also checked that ℓ 1 < 0. An illustration of the sustained oscillations corresponding to the stable limit cycle of (13) is depicted in Fig 5. The actSIS model incorporates two novel mechanisms as compared with the SIS model: a feedback mechanism for the effective infection rates, and a time scale separation between the state of the system and the state used in the feedback law. The qualitative differences we have shown for the actSIS model are due to both of the mechanisms. We have observed sustained oscillations even when some of the individuals are risk-ignorers and under relaxed assumptions on interconnections. We will examine the broader set of possibilities in future work and consider applications in other biological and socio-ecological processes. (13) exhibits stable oscillations. States p 1 and p 2 initially decrease quickly due to different mechanisms: p 1 decreases because p s2 (0) < µ T , which drives β 12 down and p 2 decreases because p s1 (0) > µ A , which drives β 21 down. The low p 1 leads to a decrease in p s1 , which in turn drives up β 21 , thus increasing p 2 . The increase in p 2 leads to an increase in p s2 , which in turn drives up p 1 . The process repeats resulting in sustained oscillations. b-c): The time evolution of the effective infection rates. Parameters used for the simulation:β = 8.8, µ T = 0.45, ν T = 0.8, µ A = 0.3, ν A = 20. Initial conditions: p 1 (0) = 0.55, p 2 (0) = 0.45, p s1 (0) = 0.45, p s2 (0) = 0.66. On the epidemic of financial crises MAT-CONT: a MATLAB package for numerical bifurcation analysis of ODEs On the definition and the computation of the basic reproduction ratio R0 in models for infectious diseases in heterogeneous populations Dynamical resonance can account for seasonality of influenza epidemics Epidemiological models and lyapunov functions A feedback SIR (fSIR) model highlights advantages and limitations of infection-based social distancing The spread of awareness and its impact on epidemic outbreaks Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields Modeling social interactions: Identification, empirical methods and policy implications Lecture Notes in Biomathematics Dynamical Systems in Neuroscience Epidemiological modeling of news and rumors on twitter A contribution to the mathematical theory of epidemics A deterministic model for gonorrhea in a nonhomogeneous population Dynamics of influenza A drift: the linear three-strain model On the dynamics of deterministic epidemic propagation over networks Oscillations in rational economies Derivation of rates from summation data by the catalytic curve Analysis and control of epidemics: A survey of spreading processes on complex networks Adaptive susceptibility and heterogeneity in contagion models on networks Hopf bifurcations and limit cycles in evolutionary network dynamics Group movement 'initiation' and statedependent decision-making Generalized epidemic mean-field model for spreading processes over multilayer complex networks Control across scales by positive and negative feedback Sigmoidal acquisition curves are good indicators of conformist transmission Mechanistic modelling of multiple waves in an influenza epidemic or pandemic Sentiment contagion in complex networks Proposition 4.2. Denote m = (1 − p * 1 )(1 − p * 2 ), q = p * 1 β * 21 /p * 2 , s = p * 2 β * 12 /p * 1 , v = p * 1 β ′ * 21 /p * 2 , w = p * 2 β ′ * 12 /p * 1 , β * 12 =βα 1 (p * 2 ) , β * 21 =βα 2 (p * 1 )). Then, for system (13), the non-hyperbolicity condition (H1) in Theorem 4.1 is satisfied ifProof. For a four-dimensional system to satisfy (H1) the eigenvalues of the Jacobian D p f | (p * ,β * ) must satisfyfor some a = 0,b ∈ R,and c > 0. We compute the Jacobianwhich has eigenvalues that satisfy