# Coupled Evolutionary Behavioral and Disease Dynamics under Reinfection Risk

Abhisek Satapathi, Narendra Kumar Dhar, Ashish R. Hota and Vaibhav Srivastava
Abhisek Satapathi and Ashish R. Hota are with the Department of Electrical
Engineering, IIT Kharagpur, India. Narendra Kumar Dhar is with the School of
Computing and Electrical Engineering, IIT Mandi, India. Vaibhav Srivastava is
with the Department of Electrical and Computer Engineering, Michigan State
University, East Lansing, MI 48824. (Email<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>vaibhav@msu.edu). This work
was supported in part by ARO grant W911NF-18-1-0325. A preliminary version of
this work Satapathi et al. (2022) appeared at the American Control Conference,
2022.

###### Abstract

We study the interplay between epidemic dynamics and human decision making for
epidemics that involve reinfection risk; in particular, the susceptible-
infected-susceptible (SIS) and the susceptible-infected-recovered-infected
(SIRI) epidemic models. In the proposed game-theoretic setting, individuals
choose whether to adopt protection or not based on the trade-off between the
cost of adopting protection and the risk of infection; the latter depends on
the current prevalence of the epidemic and the fraction of individuals who
adopt protection in the entire population. We define the coupled epidemic-
behavioral dynamics by modeling the evolution of individual protection
adoption behavior according to the replicator dynamics. For the SIS epidemic,
we fully characterize the equilibria and their stability properties. We
further analyze the coupled dynamics under timescale separation when
individual behavior evolves faster than the epidemic, and characterize the
equilibria of the resulting discontinuous hybrid dynamical system for both SIS
and SIRI models. Numerical results illustrate how the coupled dynamics
exhibits oscillatory behavior and convergence to sliding mode solutions under
suitable parameter regimes.

##  1 Introduction

Infectious diseases or epidemics spread through human society via social
interactions among infected and healthy individuals. In the absence of
suitable medicines or vaccines, particularly at the onset of an epidemic,
individuals adopt protective measures (such as wearing masks or maintaining
social distancing) in order to avoid becoming infected, often in a strategic
and decentralized manner. Thus, adoption of protective behavior both affects
and is affected by the level of infection prevalence. Therefore, it is
critical to examine the interplay between epidemic dynamics and human decision
making for improved preparedness against future pandemics.

There is a rich body of work devoted to developing mathematical models that
capture the dynamic evolution of epidemics Mei et al. (2017); Nowzari et al.
(2016). There are two broad classes of epidemic models: (i) epidemics where an
individual upon recovery is no longer at risk of future infection (e.g., the
susceptible-infected-recovered (SIR) epidemic model) and (ii) epidemics where
recovery from infection does not lead to immunity from future infections.
Motivated by a large number of infectious diseases where recovered individuals
are still at risk of future infection, we focus on the second class of
epidemics in this work; in particular, the susceptible-infected-susceptible
(SIS) and the susceptible-infected-recovered-infected (SIRI) Pagliara et al.
(2018); Pagliara and Leonard (2020) epidemic models. In the SIRI model, the
infection rate of individuals who have been infected in the past is different
from those who have never been infected which captures diseases that impart
compromised or strengthened immunity after the initial infection. This is in
contrast with the SIS epidemic where a recovered individual encounters an
identical infection risk as an individual who has never been infected.

Selfish adoption of protective behavior against epidemics has been studied in
the framework of game theory Chang et al. (2020); Huang and Zhu (2021); both
static or single-shot games to model vaccination decisions Hota and Sundaram
(2019); Trajanovski et al. (2015) and dynamic games that model evolution of
protective decisions Theodorakopoulos et al. (2013); Eksin et al. (2017);
Huang and Zhu (2020); Hota et al. (2022) have been analyzed for the SIS
epidemic in recent past. Most of these works explore structural properties of
the equilibrium strategies and epidemic evolution under those strategies.
Specifically, in Theodorakopoulos et al. (2013), individuals in a well-mixed
population decide whether or not to adopt protection (which eliminates the
risk of infection) at a certain cost, and epidemic evolution under optimal
protection strategies is studied. This setting is generalized in Eksin et al.
(2017) and Hota et al. (2022) to include different types of networked
interactions.

However, the problem of a large group of strategic individuals finding or
learning equilibrium strategies is quite challenging Daskalakis et al. (2009);
Sandholm (2010). In settings with a large population of agents, evolutionary
dynamics have been proposed and their convergence behavior to the set of
equilibria are analyzed Hofbauer and Sigmund (2003); Sandholm (2010). This
class of learning dynamics involves repeated play of a static game where the
payoff functions remain unchanged. It is only recently that evolutionary
learning with dynamic payoff functions have been investigated Arcak and
Martins (2020); Park et al. (2019).

In contrast with the above settings, epidemic games exhibit a dynamically
changing proportion of healthy and infected individuals due to epidemic
dynamics. A recent work Elokda et al. (2021a) investigated this dynamic
evolution of infection dynamics and strategic decisions for a class of
susceptible-asymptomatic-infected-recovered (SAIR) epidemic setting in the
framework of dynamic population games Elokda et al. (2021b). Coupled evolution
of disease and behavior have recently been investigated in Martins et al.
(2022); Khazaei et al. (2021) for SIRS and SEIR epidemic models, respectively.
While the setting in Martins et al. (2022) allows for reinfection, the payoff
function of an individual does not directly depend on the infection
prevalence, rather the payoffs are governed by a dynamic reward term likely
designed by a central authority. In a contemporary work Frieswijk et al.
(2022), the authors consider a coupled SIS epidemic-behavioral setting in
which the payoff contains a social influence factor and a cost proportional to
the current level of infection, and imitation dynamics is used to model the
evolution of behavior.

Few recent works have also explored the notion of timescale separation to
study the interaction of evolutionary learning dynamics with other dynamical
systems Weitz et al. (2016); Gong and Cao (2021). In the context of epidemics,
Al-Radhawi et al. (2021) explores slow-fast dynamics in the context of SIR
epidemics. SIS dynamics with risk-tolerator or risk-averter individuals that
respond differently is studied in Zhou et al. (2020). Neither works consider
any behavioral or game-theoretic response to epidemics. Furthermore, to the
best of our knowledge, game-theoretic analysis of protection strategies for
the SIRI epidemic model has yet not been reported.

In this paper, we build upon the above line of work, and consider the SIS and
SIRI epidemic settings where a large population of individuals choose whether
to adopt protection or to remain unprotected as the epidemic evolves. For a
healthy (susceptible or recovered) individual, being protected reduces the
probability of becoming infected. Similarly, an infected individual who adopts
protection transmits or causes a new infection with a smaller probability
compared to an infected individual who is unprotected.111Thus, our
formulation, with partially effective protection, generalizes the setting in
Theodorakopoulos et al. (2013) where adopting protection completely eliminated
the risk of becoming infected for an individual. For a susceptible individual,
adopting protection comes with a cost (such as discomfort due to wearing masks
or adhering to social distancing protocols) which is weighed with the
instantaneous risk of becoming infected; the latter depends on the current
epidemic prevalence and the proportions of individuals in different epidemic
states adopting protection. In the SIRI model, recovered individuals behave
similar to susceptible individuals, albeit with different infection
probabilities.

In order to lay the foundations for a rigorous investigation of coupled
epidemic and behavioral dynamics, we focus on the replicator dynamics Hofbauer
and Sigmund (2003); Sandholm (2010) to model the evolution of protection
decisions and study their interaction with SIS and SIRI epidemic dynamics. The
major contributions of this work are fourfold. First, for the coupled SIS
epidemic and replicator dynamics, we completely characterize the equilibria
(existence and local stability) and show how the stability of different
equilibrium points get exchanged as certain parameters change. Second, we
explore the behavior of the coupled dynamics under timescale separation
leading to a slow-fast dynamical system Berglund and Gentz (2006). Third, for
the coupled SIRI epidemic and replicator dynamics, we determine all possible
equilibria and analyze the coupled dynamics under timescale separation. For
both (SIS and SIRI) coupled dynamics, we specifically focus on the case where
the replicator dynamics are faster and characterize the evolution of the
epidemic by analyzing the resulting slow dynamics which takes the form of a
variable structure system in the limiting regime. We show that the infected
proportion may converge to an endemic sliding mode of the hybrid dynamics for
a certain range of parameters. Fourth, we illustrate additional insights into
the transient behavior of these coupled dynamics using numerical simulations.
Particularly, we illustrate the impact of the timescale separation parameter
and oscillatory convergence to the equilibrium point.

The remainder of the paper is organized as the following. In Section 2, we
introduce the coupled SIS epidemic and replicator dynamics, characterize the
associated equilibria and their bifurcations, and study the transient behavior
of the system trajectories under timescale separation. In Section 3, we
introduce the coupled SIRI epidemic and replicator dynamics, determine the
associated equilibria, and study these dynamics under timescale separation. We
present additional insights into these coupled dynamics using numerical
simulations in Section 4. Conclusions are drawn in Section 5, and the proofs
of the technical results are presented in the Appendix.

##  2 Coupled SIS Epidemic and Evolutionary Behavioral Model

In this section, we formally introduce the coupled evolution of the SIS
epidemic and protection adoption behavior in a homogeneous large-population
setting. Let the proportion of susceptible and infected individuals be $s(t)$
and $y(t)$, respectively. Both $s(t),y(t)\in[0,1]$ with $s(t)+y(t)=1$ for all
$t\geq 0$. Individuals choose whether to adopt protection against the epidemic
or not; these actions are denoted by $\mathtt{P}$ and $\mathtt{U}$,
respectively. Consequently, the population state at time $t$ is defined as
$x(t):=[x_{\mathtt{SU}}(t)\quad x_{\mathtt{SP}}(t)\quad
x_{\mathtt{IU}}(t)\quad x_{\mathtt{IP}}(t)]^{\intercal}\in\Delta_{4}$, where
$\Delta_{4}$ is the probability simplex in 4, $x_{\mathtt{SU}}$ denotes the
proportion of individuals who are susceptible and choose to remain
unprotected, $x_{\mathtt{IP}}$ denotes the proportion of infected individuals
who adopt protection, etc. At time $t$, we have
$x_{\mathtt{SU}}(t)+x_{\mathtt{SP}}(t)=s(t)$,
$x_{\mathtt{IU}}(t)+x_{\mathtt{IP}}(t)=y(t)$, and
$\mathds{1}^{\intercal}x(t)=1$ where $\mathds{1}$ is a vector of appropriate
dimension with all entries $1$.

Individuals choose their action to maximize their payoffs which depends on the
population state $x(t)$ (which includes information regarding infection
prevalence since $x_{\mathtt{IU}}(t)+x_{\mathtt{IP}}(t)=y(t)$). For an
infected individual, there is no further risk of infection, and as a result,
we define its payoff to be constant parameters given by $-c_{\mathtt{IU}}$ if
it remains unprotected and $-c_{\mathtt{IP}}$ if it adopts protection,
respectively. In particular, $c_{\mathtt{IU}}>0$ captures the cost of breaking
isolation/quarantine protocols for an infected individual or could represent
worsening of the disease if an infected individual continues to be
irresponsible. Thus, we assume $c_{\mathtt{IU}}>c_{\mathtt{IP}}\geq 0$.

A susceptible individual trades off the cost of adopting protection, denoted
by $c_{P}>0$, and expected cost of becoming infected. The latter is the
product of the loss upon infection $L>0$ and the instantaneous probability of
becoming infected which depends on its action as well as the population state.
Specifically, let $\beta_{u}$ and $\beta_{p}$ denote the probabilities of an
infected individual causing a new infection if it is unprotected and
protected, respectively. We impose the natural assumption
$\beta_{u}>\beta_{p}\geq 0$ throughout the paper. Consequently, the
instantaneous probability of infection for an unprotected susceptible
individual at population state $x$ is given by
$\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}}$. Similarly, let a
protected susceptible individual be $\alpha\in(0,1)$ times (less) likely to
become infected compared to an unprotected susceptible individual. Thus, the
instantaneous probability of infection for a protected susceptible individual
at population state $x$ is given by
$\alpha(\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}})$. The payoff vector
of individuals at population state $x$ is now defined as

$F(x)\\!:=\begin{bmatrix}F_{\mathtt{SU}}(x)\\\ F_{\mathtt{SP}}(x)\\\
F_{\mathtt{IU}}(x)\\\
F_{\mathtt{IP}}(x)\end{bmatrix}\\!=\\!\begin{bmatrix}-\\!L(\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}})\\\
-c_{P}-\\!L\alpha(\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}})\\\
-c_{\mathtt{IU}}\\\ -c_{\mathtt{IP}}\end{bmatrix},$ (1)

where $F_{\mathtt{SU}}(x)$ denote the payoff for an individual who is
susceptible and unprotected at population state $x$, and so on. Susceptible
individuals who adopt protection pay a cost $c_{P}$ but experience a reduced
risk of becoming infected scaled by factor $\alpha$ as discussed above.

The population state $x$ is dynamic and time-varying, both due to the epidemic
dynamics as well as the evolutionary dynamics of strategically adopting
protective behavior. Since the proportion of susceptible and infected
individuals is time-varying, the framework of classical population games is
not applicable. We now introduce the coupled epidemic and evolutionary
dynamics. We first introduce some notation. Let $z_{S}(t)\in[0,1]$ denote the
fraction of susceptible individuals who remain unprotected, i.e.,
$x_{\mathtt{SU}}(t)=z_{S}(t)s(t)$ and $x_{\mathtt{SP}}(t)=(1-z_{S}(t))s(t)$.
Similarly, let $z_{I}(t)\in[0,1]$ denote the fraction of infected individuals
who remain unprotected.

Due to the presence of both unprotected and protected individuals with
different infection probabilities, the evolution of infected proportion is
given by

$\displaystyle\dot{y}(t)$ $\displaystyle=(x_{\mathtt{SU}}(t)+\alpha
x_{\mathtt{SP}}(t))(\beta_{u}x_{\mathtt{IU}}(t)+\beta_{p}x_{\mathtt{IP}}(t))-\gamma
y(t)$
$\displaystyle=\big{[}(1-y(t))(z_{S}(t)+\alpha(1-z_{S}(t)))\cdot(\beta_{u}z_{I}(t)+\beta_{p}(1-z_{I}(t)))-\gamma\big{]}y(t)$
$\displaystyle=:f_{y}(y(t),z_{S}(t),z_{I}(t)).$ (2)

The above dynamics is analogous to the conventional scalar SIS epidemic
dynamics with effective infection rate
$\beta_{\textup{eff}}(z_{S},z_{I})=(z_{S}+\alpha(1-z_{S}))(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))$
which now depends on the efficacy of protection and the fractions that adopt
protection. In particular, if $\alpha=1$ and $\beta_{p}=\beta_{u}=\beta$,
i.e., protection is not effective, the effective infection rate is $\beta$.

We now define the evolutionary dynamics to capture the evolution of
unprotected fractions in both susceptible and infected subpopulations. We
focus on the class of replicator dynamics Sandholm (2010) in this work and
assume that susceptible individuals only replicate the strategies of other
susceptible individuals (likewise for infected individuals). Consequently for
susceptible individuals, we obtain

$\displaystyle\dot{z}_{S}(t)$
$\displaystyle={z}_{S}(t)(1-{z}_{S}(t))\big{[}[F]_{\mathtt{SU}}-[F]_{\mathtt{SP}}\big{]}$
$\displaystyle={z}_{S}(t)(1-{z}_{S}(t))\big{[}c_{P}-L(1-\alpha)(\beta_{u}z_{I}(t)+\beta_{p}(1-z_{I}(t)))y(t)\big{]}$
$\displaystyle=:f_{S}(y(t),z_{S}(t),z_{I}(t)).$ (3)

Similarly, for infected individuals, we have

$\displaystyle\dot{z}_{I}(t)$
$\displaystyle={z}_{I}(t)\big{[}-c_{\mathtt{IU}}+(c_{\mathtt{IU}}{z}_{I}(t)+c_{\mathtt{IP}}(1-{z}_{I}(t)))\big{]}$
$\displaystyle={z}_{I}(t)(1-{z}_{I}(t))(c_{\mathtt{IP}}-c_{\mathtt{IU}})$
$\displaystyle=:f_{I}(y(t),z_{S}(t),z_{I}(t)).$ (4)

Thus, equations (2), (3) and (4) characterize the coupled evolution of the
epidemic and population states.

######  Lemma 1 (Invariant set of Coupled SIS Epidemic).

For the coupled SIS epidemic and evolutionary behavioral dynamics defined by
(2), (3) and (4), the set
$\\{(y,z_{S},z_{I})\;|\;(y,z_{S},z_{I})\in[0,1]^{3}\\}$ is invariant.

###### Proof.

For $y(t)\in\\{0,1\\}$, $\dot{y}$ is either zero or negative. Likewise, for
$z_{I},z_{S}\in\\{0,1\\}$, $\dot{z}_{I}=\dot{z}_{S}=0$. The result follows
from Nagumo’s theorem (Blanchini and Miani, 2008, Theorem 4.7). ∎

###  2.1 Equilibrium Characterization and Stability Analysis

We now examine the equilibrium points of the above coupled epidemic-replicator
dynamics and their stability properties. First, we consider the evolution of
$z_{I}$ in (4) which does not depend on $y$ and $z_{S}$. There are two
stationary points $z_{I}=0$ and $z_{I}=1$, and it is easy to see that, for
$c_{\mathtt{IU}}>c_{\mathtt{IP}}$, $z_{I}=1$ is unstable and $z_{I}=0$ is
exponentially stable with basin of attraction $[0,1)$. It is also quite
intuitive that if $c_{\mathtt{IP}}<c_{\mathtt{IU}}$, infected individuals
prefer to use protection and the strategy for infected individuals should
converge to adopting protection.

Thus, in the remainder of this paper, we only focus on equilibria with
$z_{I}=0$. We begin with introducing a few variables that will be used to
define the equilibrium points:

$y^{*}_{u}:=1-\frac{\gamma}{\beta_{p}},\quad
y^{*}_{\mathrm{int}}:=\frac{c_{P}}{L(1-\alpha)\beta_{p}},\quad
y^{*}_{p}:=1-\frac{\gamma}{\alpha\beta_{p}},\quad\text{and}\quad
z^{*}_{S,\mathrm{int}}:=\frac{1}{1-\alpha}\left[\frac{\gamma}{\beta_{p}(1-y^{*}_{\mathrm{int}})}-\alpha\right].$

We now define all possible equilibria of the coupled SIS epidemic and
evolutionary behavior dynamics (2–4) corresponding to $z_{I}=0$:

$\mathbf{E0}=(0,0,0),\quad\mathbf{E1}=(0,1,0),\quad\mathbf{E2}=(y^{*}_{u},1,0),\quad\mathbf{E3}=(y^{*}_{\mathrm{int}},z^{*}_{S,\mathrm{int}},0),\text{
and}\quad\mathbf{E4}=(y_{p}^{*},0,0).$

At $\mathbf{E0}$ everyone adopts protection and there is no infection. At
$\mathbf{E1}$, there is no infection, and susceptible individuals do not adopt
protection. $\mathbf{E2}$ is an endemic equilibrium, i.e., a fraction of the
population is infected, at which susceptible individuals continue to remain
unprotected. $\mathbf{E3}$ is an endemic equilibrium at which some susceptible
individuals adopt protection. Finally, $\mathbf{E4}$ is an endemic equilibrium
at which all susceptible individuals adopt protection. The existence and
stability of these equilibria are summarized in Table 1 and established in the
following proposition whose proof is presented in Appendix A.

Table 1: Existence and stability of equilibria of the coupled dynamics (2), (3) and (4). Epidemic | Endemic | Equilibria
---|---|---
Parameters | Infection Level | $\mathbf{E1}:(0,1,0)$ | $\mathbf{E2}:(y^{*}_{u},1,0)$ | $\mathbf{E3}:(y^{*}_{\mathrm{int}},z^{*}_{S,int},0)$ | $\mathbf{E4}:(y^{*}_{p},0,0)$
$\gamma>\beta_{p}$ | $-$ | $\checkmark$, stable | $-$ | $-$ | $-$
$\alpha\beta_{p}<\gamma<\beta_{p}$ | $y^{*}_{u}<y^{*}_{\mathrm{int}}$ |  $\checkmark$, unstable | $\checkmark$, stable | $-$ | $-$
$y^{*}_{\mathrm{int}}<y^{*}_{u}$ |  $\checkmark$, unstable |  $\checkmark$, unstable | $\checkmark$, stable | $-$
$\gamma<\alpha\beta_{p}$ | $y^{*}_{u}<y^{*}_{\mathrm{int}}$ |  $\checkmark$, unstable | $\checkmark$, stable | $-$ |  $\checkmark$, unstable
$y^{*}_{p}<y^{*}_{\mathrm{int}}<y^{*}_{u}$ |  $\checkmark$, unstable |  $\checkmark$, unstable | $\checkmark$, stable |  $\checkmark$, unstable
$y^{*}_{\mathrm{int}}<y^{*}_{p}$ |  $\checkmark$, unstable |  $\checkmark$, unstable | $-$ | $\checkmark$, stable

######  Proposition 1 (Equilibria and Stability).

For the equilibrium points of the coupled SIS epidemic and evolutionary
behavioral dynamics (2–4) corresponding to $z_{I}=0$, the following statements
hold:

  1. 1.

$\mathbf{E0}$ exists for all parameter regimes, and is unstable;

  2. 2.

$\mathbf{E1}$ exists for all parameter regimes, is stable if
$\beta_{p}<\gamma$, and is unstable, otherwise;

  3. 3.

$\mathbf{E2}$ exists only when $\beta_{p}>\gamma$, is stable when
$y^{*}_{u}<y^{*}_{\mathrm{int}}$, and is unstable otherwise;

  4. 4.

$\mathbf{E3}$ exists only when $y^{*}_{p}<y^{*}_{\mathrm{int}}<y^{*}_{u}$, and
is stable;

  5. 5.

$\mathbf{E4}$ exists only when $\gamma<\alpha\beta_{p}$, is stable when
$y^{*}_{p}>y^{*}_{\mathrm{int}}$, and is unstable otherwise.

######  Remark 1.

Note that as the recovery rate $\gamma$ increases, the infected proportion at
the stable equilibrium decreases. In practice, $\gamma$ is improved via direct
intervention of authorities by augmenting healthcare facilities. When further
resource augmentation is not possible, a more effective protection scheme with
a smaller value of $\alpha$ would result in a smaller value of $y^{*}_{p}$. In
addition, if $\alpha$ is reduced further such that $\gamma>\alpha\beta_{p}$,
the dynamics exhibits a new stable endemic infection level
$y^{*}_{\mathrm{int}}$ which also decreases as $\alpha$ is reduced. Thus, our
result shows that indirect intervention by facilitating availability of more
effective protection schemes and letting individuals choose whether to adopt
protection in a decentralized manner would significantly contribute towards
ensuring a smaller endemic infection level.

###  2.2 Bifurcation Analysis

The above proposition shows that the equilibrium points exchange stability
properties as certain parameters, e.g., $\gamma,\beta_{p},\alpha$, vary. We
now numerically explore the bifurcations associated with the transition of
stability among the equilibria. For the numerical illustration we choose the
parameter values in equations (2), (3) and (4) as summarized below:

$c_{P}$ | $\alpha$ | $\beta_{u}$ | $c_{\mathtt{IU}}$ | $L$ | $\beta_{p}$ | $c_{\mathtt{IP}}$
---|---|---|---|---|---|---
$1$ | $0.5$ | $0.3$ | $2$ | $80$ | $0.15$ | $1$

We adopt the recovery rate $\gamma$ as a bifurcation parameter and use the
numerical continuation package MATCONT Dhooge et al. (2003) to compute the
bifurcation diagram shown in Fig. 1.

Figure 1: Bifurcation diagram of the equilibria of the coupled SIS epidemic
and evolutionary behavioral dynamics. (solid blue: stable branch and dashed
red: unstable branch of equilibria)

For $\gamma\to 0^{+}$, $\mathbf{E4}$ is the stable equilibrium, while
$\mathbf{E0}$, $\mathbf{E1}$, and $\mathbf{E2}$ are unstable. As the value of
$\gamma$ is increased at point T3 in Fig. 1, $\mathbf{E4}$ exchanges stability
to $\mathbf{E3}$ in a transcritical bifurcation. Note that unstable branch of
$\mathbf{E3}$ is not visible since it is associated with negative values of
$z_{S}$. As $\gamma$ is increased, the fraction of susceptible population
adopting the protection decreases and at T2, $\mathbf{E3}$ exchanges stability
to $\mathbf{E2}$ in another transcritical bifurcation. Again, the unstable
branch of $\mathbf{E3}$ at T2 corresponds to $z_{s}>1$ and is not visible in
the diagram. Upon further increasing $\gamma$, the fraction of infected
population continues to decrease and at T1, $\mathbf{E2}$ exchanges stability
with the disease-free equilibria $\mathbf{E1}$. The unstable branch of
$\mathbf{E2}$ at T1 corresponds to negative values of $y$ and is not visible.

Another transcritical bifurcation takes place at T0
$(\gamma=\alpha\beta_{p})$, where $\mathbf{E0}$ and $\mathbf{E4}$ cross. For
$\gamma<\alpha\beta_{p}$ (resp., $\gamma>\alpha\beta_{p}$), $\mathbf{E0}$ has
two (resp., one) eigenvalues in the right-half plane. As $\mathbf{E4}$
approaches T0 from $y>0$, it has one eigenvalue in the right-half plane, while
for $y<0$ near T0, $\mathbf{E4}$ is stable. Thus, the transcritical
bifurcation at T0 corresponds to exchange of stable and unstable eigenvalues
of two unstable equilibria.

###  2.3 Coupled Epidemic-Behavioral Dynamics under Timescale Separation

We have thus far assumed that the the epidemic and the replicator dynamics
evolve at the same time-scale. However, it is not strictly necessary for the
coupled dynamics. In order to obtain further insights into their behavior, we
now study the coupled epidemic-behavioral dynamics (2), (3) and (4) under
timescale separation. In particular, we focus on the case in which the
replicator dynamics evolves faster than the epidemic dynamics; indeed in the
modern era, there is an increased awareness about infectious diseases due to
publicly available testing data at regular intervals, awareness campaigns by
public health authorities, spread of information via social media, etc., which
shapes human response at a much faster time-scale. To this end, we model the
coupled dynamics as a slow-fast system:

$\displaystyle\begin{split}\dot{y}(t)&=f_{y}(y(t),z_{S}(t),z_{I}(t))\\\
\epsilon\dot{z}_{S}(t)&=f_{S}(y(t),z_{S}(t),z_{I}(t))\\\
\epsilon\dot{z}_{I}(t)&=f_{I}(y(t),z_{S}(t),z_{I}(t)),\end{split}$ (5)

where $\epsilon\in(0,1]$ is a timescale separation variable Berglund and Gentz
(2006).

At a given epidemic prevalence $y$, we characterize the (stable) equilibria of
the fast system involving the replicator dynamics with states $(z_{S},z_{I})$.
For reasons discussed earlier, we focus on equilibria with $z_{I}=0$. It is
now easy to see that if $y\neq y^{*}_{\mathrm{int}}$, there are two
equilibrium points: $(0,0)$ and $(1,0)$. If $y=y^{*}_{\mathrm{int}}$, then
$(z_{S},0)$ is an equilibrium point of the fast system for any
$z_{S}\in[0,1]$. Following analogous arguments as in the proof of Proposition
1, it follows that $(0,0)$ is locally stable for the fast system when
$y>y^{*}_{\mathrm{int}}$ and $(1,0)$ is locally stable for the fast system
when $y<y^{*}_{\mathrm{int}}$. Consequently, we obtain the following reduced
dynamics for the slow system which approximates the coupled dynamics (5) in
the limit $\epsilon\to 0$ as

$\displaystyle\dot{y}(t)$
$\displaystyle=\begin{cases}\big{[}(1-y(t))\beta_{p}-\gamma\big{]}y(t),\quad&\text{if
}y(t)<y^{*}_{\mathrm{int}},\\\
\big{[}(1-y(t))\alpha\beta_{p}-\gamma\big{]}y(t),\quad&\text{if
}y(t)>y^{*}_{\mathrm{int}},\end{cases}$ (6) and if $\displaystyle
y(t)=y^{*}_{\mathrm{int}},\text{we have}$ $\displaystyle\dot{y}(t)$
$\displaystyle\in\Big{\\{}\big{[}(1-y^{*}_{\mathrm{int}})\beta_{p}(z_{S}+\alpha(1-z_{S}))-\gamma\big{]}y^{*}_{\mathrm{int}}\;\Big{|}\;z_{S}\in[0,1]\Big{\\}}.$

In particular, since the reduced dynamics is an instance of a discontinuous
dynamical system, we define the dynamics at the point of discontinuity
$y=y^{*}_{\mathrm{int}}$ as a differential inclusion which is also the convex
combination of the dynamics on both sides of $y=y^{*}_{\mathrm{int}}$. It is
easy to see that the right hand side of the dynamics is measurable and is
locally essentially bounded, and therefore (6) admits a Filippov solution
Cortes (2008). We now establish the convergence behavior of the infected
fraction under (6).

######  Proposition 2 (Trajectories under fast behavioral response).

For the epidemic dynamics (6) with $y(0)\neq 0$, the following statements hold

  1. 1.

if $y^{*}_{u}\leq 0$, then $y(t)$ monotonically decreases and converges to the
origin;

  2. 2.

if $0<y^{*}_{u}<y^{*}_{\mathrm{int}}$, then $y(t)$ monotonically converges to
$y_{u}^{*}$;

  3. 3.

if $y^{*}_{p}<y^{*}_{\mathrm{int}}<y^{*}_{u}$, then $y(t)$ converges to
$y^{*}_{\mathrm{int}}$ which acts as a sliding mode of the dynamics;

  4. 4.

if $y^{*}_{\mathrm{int}}\leq y^{*}_{p}$, then $y(t)$ monotonically converges
to $y_{p}^{*}$.

The proof of the above proposition is presented in Appendix A.

######  Remark 2.

The dynamics in (6) potentially represents a class of non-pharmaceutical
interventions where authorities impose social distancing measures that reduces
the infection rate by a factor $\alpha$ when the infection prevalence exceeds
a threshold. Thus, the result in Proposition 2 is potentially of independent
interest. Further, the above proposition generalizes analogous results
obtained in prior works Theodorakopoulos et al. (2013); Hota et al. (2022)
which assumed $\alpha=0$, i.e., adopting protection completely eliminates risk
of infection.

##  3 Coupled SIRI Epidemic and Evolutionary Behavioral Model

In the SIS epidemic model, a recovered individual encounters the same
infection rate as an individual who has never been infected. However, in many
infectious diseases, initial infection could lead to compromised immunity
Marshall (2020) or it might even lead to reduced risk of future infection
Gomes et al. (2004); Crellen et al. (2021). The susceptible-infected-
recovered-infected (SIRI) epidemic model Pagliara et al. (2018), captures the
above characteristics. In this section, we investigate the implications of
game-theoretic protection decisions on the evolution of the SIRI epidemic
dynamics.

In the SIRI epidemic model, an individual belongs to one of three possible
compartments or states: susceptible, infected and recovered; while the
proportion of individuals in each of the above states at time $t$ is denoted
by $s(t),y(t)$ and $r(t)$, respectively. The evolution of these proportions is
given by

$\displaystyle\dot{s}(t)$ $\displaystyle=-\beta s(t)y(t)$ (7a)
$\displaystyle\dot{y}(t)$ $\displaystyle=\beta
s(t)y(t)+\hat{\beta}r(t)y(t)-\gamma y(t)$ (7b) $\displaystyle\dot{r}(t)$
$\displaystyle=\gamma y(t)-\hat{\beta}r(t)y(t),$ (7c)

where $\beta>0$ is the rate at which susceptible individuals become infected
if they encounter an infected individual, $\hat{\beta}>0$ is the rate at which
recovered individuals become infected and $\gamma$ is the rate of recovery for
infected individuals. The above transitions are depicted in Figure 2.

Thus, the above model captures settings where an individual develops immunity
upon becoming infected ($\beta>\hat{\beta}$) and when the immunity of an
individual is compromised upon infection ($\beta<\hat{\beta}$). When
$\beta=\hat{\beta}$ and susceptible and recovered states are combined, we
recover the SIS epidemic model. The classical SIR model is also obtained as a
special case when $\hat{\beta}=0$.

$\mathtt{S}$$\mathtt{I}$$\mathtt{R}$$\beta$$\gamma$$\hat{\beta}$ Figure 2:
Evolution of states in the SIRI epidemic model. Self-loops are omitted for
better clarity.

###  3.1 Equilibria and Stability without Protective Behavior

The dynamics in (7) admits a continuum of equilibria which are infection free
(IFE) and may also have an isolated endemic equilibrium (EE) where infection
level is nonzero. Mathematically, at the IFE, we have $s=s^{*},y=0,r=1-s^{*}$
where $s^{*}\in[0,1]$. At the EE, we have
$s=0,y=1-\frac{\gamma}{\hat{\beta}},r=\frac{\gamma}{\hat{\beta}}$. Before we
analyze the implications of strategic adoption of protective behavior, we
summarize the main result from Pagliara et al. (2018) on the existence and
stability of the above equilibria.

######  Theorem 1 (Theorem 2, Lemma 1, Lemma 2 Pagliara et al. (2018)).

Let $R_{0}:=\frac{\beta}{\gamma}$, $R_{1}:=\frac{\hat{\beta}}{\gamma}$, and
$M:=\frac{1-R_{1}}{R_{0}-R_{1}}$. Then, the SIRI dynamics exhibits the
following behavior.

  1. 1.

Infection-free: If $R_{0}<1$ and $R_{1}<1$, EE is not an equilibrium, while
all points in the IFE are locally stable. Further, $y(t)$ decays monotonically
to $0$.

  2. 2.

Endemic: If $R_{0}>1$ and $R_{1}>1$, all points in the IFE are unstable while
the EE exists and is locally stable. Further, all solutions reach EE as
$t\to\infty$.

  3. 3.

Epidemic: If $R_{0}>1$ and $R_{1}\leq 1$, then IFE is locally stable if
$s^{*}<M$ and unstable if $s^{*}>M$. Further, as $t\to\infty$, all solutions
reach IFE with $s^{*}<M$.

  4. 4.

Bistable: If $R_{0}\leq 1$ and $R_{1}>1$, EE exists and is locally stable. IFE
with $s^{*}>M$ is locally stable while IFE with $s^{*}<M$ is unstable. If
$y(0)<1-M(R_{0}M)^{-R_{0}/R_{1}}$, then the dynamics reach a point at the IFE
as $t\to\infty$; otherwise it reaches the EE.

We now introduce a population game model to capture how individuals adopt
protective behavior, followed by analyzing the resulting evolution of the SIRI
epidemic.

###  3.2 Population Game Model of Protection Adoption

We build upon the formulation in the previous section. Consider a large
population of individuals, each in one of three possible infection states, who
choose whether to adopt protection or remain unprotected. The population state
$x:=[x_{\mathtt{SU}}\quad x_{\mathtt{SP}}\quad x_{\mathtt{IU}}\quad
x_{\mathtt{IP}}\quad x_{\mathtt{RU}}\quad
x_{\mathtt{RP}}]^{\intercal}\in[0,1]^{6}$, where $x_{\mathtt{RU}}$ denotes the
proportion of recovered individuals who remain unprotected, $x_{\mathtt{RP}}$
denotes the proportion of recovered individuals who adopt protection, and the
other states as described earlier. At time $t$, we have
$x_{\mathtt{SU}}(t)+x_{\mathtt{SP}}(t)=s(t)$,
$x_{\mathtt{IU}}(t)+x_{\mathtt{IP}}(t)=y(t)$,
$x_{\mathtt{RU}}(t)+x_{\mathtt{RP}}(t)=r(t)$ and
$\mathds{1}^{\intercal}x(t)=1$ as before. In the subsequent analysis, we omit
the dependence on $t$ for brevity of notation.

A susceptible individual becomes infected at a rate $\beta_{p}\geq 0$
($\beta_{u}\geq 0$) if it comes in contact with an infected individual who
adopts protection (remains unprotected). Similarly, a recovered individual
becomes reinfected at rate $\hat{\beta}_{u}\geq 0$ ($\hat{\beta}_{p}\geq 0$)
if it comes in contact with an infected and unprotected (protected)
individual. A susceptible or recovered individual is $\alpha\in(0,1)$ times
less likely to become infected (or reinfected) if it adopts protection
compared to an unprotected individual with same disease status. Building upon
the discussion in the previous section, we define the payoff vector as

$\displaystyle F(x)$ $\displaystyle\\!:=\begin{bmatrix}F_{\mathtt{SU}}(x)\\\
F_{\mathtt{SP}}(x)\\\ F_{\mathtt{IU}}(x)\\\ F_{\mathtt{IP}}(x)\\\
F_{\mathtt{RU}}(x)\\\
F_{\mathtt{RP}}(x)\end{bmatrix}\\!=\\!\begin{bmatrix}-L(\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}})\\\
-c_{P}\\!-\\!L\alpha(\beta_{u}x_{\mathtt{IU}}+\beta_{p}x_{\mathtt{IP}})\\\
-c_{\mathtt{IU}}\\\ -c_{\mathtt{IP}}\\\
-L(\hat{\beta}_{u}x_{\mathtt{IU}}+\hat{\beta}_{p}x_{\mathtt{IP}})\\\
-c_{P}\\!-\\!L\alpha(\hat{\beta}_{u}x_{\mathtt{IU}}+\hat{\beta}_{p}x_{\mathtt{IP}})\end{bmatrix},$
(8)

where $L,c_{P},c_{\mathtt{IU}}$ and $c_{\mathtt{IP}}$ are as defined earlier.

The payoffs for susceptible and infected individuals defined above coincide
with the payoffs in case of the SIS epidemic stated in (1). Due to risk of
reinfection, recovered individuals behave in a similar manner as susceptible
individuals and evaluate the trade-off between cost of adopting protection and
the instantaneous infection risk while choosing their protection status. As
before, we assume $c_{\mathtt{IU}}>c_{\mathtt{IP}}\geq 0$.

We denote the time-varying proportions of susceptible, infected and recovered
individuals who choose to remain unprotected by $z_{S}$, $z_{I}$, and $z_{R}$,
respectively. In particular, we have $x_{\mathtt{RU}}=z_{R}r$,
$x_{\mathtt{RP}}=(1-z_{R})r$ and so on. These proportions evolve according to
the replicator dynamics with payoffs defined in (8). We now state the coupled
evolution of disease and evolutionary behavioral dynamics as follows:

$\displaystyle\dot{s}$
$\displaystyle=-(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))(z_{S}+\alpha(1-z_{S}))sy$
(9a) $\displaystyle\dot{y}$
$\displaystyle=(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))(z_{S}+\alpha(1-z_{S}))sy$
$\displaystyle\qquad\quad+\\!(\hat{\beta}_{u}z_{I}+\hat{\beta}_{p}(1\\!-\\!z_{I}))(z_{R}+\alpha(1-\\!z_{R}))ry\\!-\\!\gamma
y$ (9b) $\displaystyle\dot{r}$
$\displaystyle=-\\!(\hat{\beta}_{u}z_{I}+\hat{\beta}_{p}(1\\!-\\!z_{I}))(z_{R}+\alpha(1-z_{R}))ry\\!+\gamma
y$ (9c) $\displaystyle\dot{z}_{S}$
$\displaystyle=z_{S}(1\\!-\\!z_{S})(c_{P}-\\!L(1-\\!\alpha)(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))y)$
(9d) $\displaystyle\dot{z}_{I}$
$\displaystyle=z_{I}(1\\!-\\!z_{I})(c_{IP}-c_{IU})$ (9e)
$\displaystyle\dot{z}_{R}$
$\displaystyle=z_{R}(1\\!-\\!z_{R})(c_{P}\\!-\\!L(1\\!-\\!\alpha)(\hat{\beta}_{u}z_{I}\\!+\\!\hat{\beta}_{p}(1\\!-\\!z_{I}))y).$
(9f)

######  Lemma 2 (Invariant Set for Coupled SIRI Dynamics).

For the coupled SIRI epidemic and evolutionary behavioral dynamics defined by
(9), the set
$\\{(s,y,r,z_{S},z_{I},z_{R})\;|\;(s,y,r,z_{S},z_{I},z_{R})\in[0,1]^{6}\\}$ is
invariant.

The proof follows from identical arguments as the proof of Lemma 1.

While $\dot{s}=-(\dot{y}+\dot{r})$ above, a complete characterization of the
equilibria and their stability properties is prohibitive due to the dynamics
being high dimensional. In order to gain insights into epidemic evolution
under game-theoretic decision-making, we analyze the coupled dynamics via
timescale separation arguments. Specifically, we analyze in detail the case
when the behavioral dynamics is much faster than the disease dynamics, derive
the stable equilibria of the fast system (consisting of variables
$z_{S},z_{I},z_{R}$) followed by analyzing the equilibria and convergence
behavior of the reduced dynamics for the slow system (consisting of variables
$s,y,r$).

###  3.3 Equilibria of the Replicator Dynamics Evolving on the Faster
Timescale

At a given infection prevalence $y$, we now characterize the equilibria of the
behavioral dynamics (9d), (9e), (9f).

######  Proposition 3.

Let $y^{*}_{\mathrm{int}}:=\frac{c_{P}}{L(1-\alpha)\beta_{p}}$ and
$\hat{y}^{*}_{\mathrm{int}}:=\frac{c_{P}}{L(1-\alpha)\hat{\beta}_{p}}$. Then,
the following are true.

  * •

Any equilibrium with $z_{I}=1$ is unstable,

  * •

For the dynamics (9d) and (9e),

    1. 1.

$(z_{S}=1,z_{I}=0)$ is the stable equilibrium when $y<y^{*}_{\mathrm{int}}$,

    2. 2.

$(z_{S}=0,z_{I}=0)$ is the stable equilibrium when $y>y^{*}_{\mathrm{int}}$,

    3. 3.

any $z_{S}\in[0,1]$, $z_{I}=0$ is an equilibrium when
$y=y^{*}_{\mathrm{int}}$.

  * •

For the dynamics (9e) and (9f),

    1. 1.

$(z_{I}=0,z_{R}=1)$ is the stable equilibrium when
$y\\!<\hat{y}^{*}_{\mathrm{int}}$,

    2. 2.

$(z_{I}=0,z_{R}=0)$ is the stable equilibrium when
$y\\!>\hat{y}^{*}_{\mathrm{int}}$,

    3. 3.

any $z_{R}\in[0,1]$, $z_{I}=0$ is an equilibrium when
$y=\hat{y}^{*}_{\mathrm{int}}$.

The proof is straightforward and is omitted. In fact, the stable equilibria of
the above behavioral dynamics corresponds to the Nash equilibrium strategies
of the underlying population game defined by the payoff vector (8) at a given
infection prevalence $y$. In particular, when $y=y^{*}_{\mathrm{int}}$,
susceptible individuals are indifferent between adopting protection or
remaining unprotected. When $y>y^{*}_{\mathrm{int}}$, it is optimal to adopt
protection which leads to $z_{S}=0$ being the stable equilibrium, and when
$y<y^{*}_{\mathrm{int}}$, it is optimal to remain unprotected which leads to
$z_{S}=1$ being the stable equilibrium. Stability of $z_{R}$ is analogous.

We now analyze the epidemic dynamics under timescale separation when the
population instantly converges to the Nash equilibrium strategies depending on
the current value of $y$. We consider two cases: $\beta_{p}>\hat{\beta}_{p}$
and $\beta_{p}<\hat{\beta}_{p}$ separately.

###  3.4 Epidemic Evolution under Strengthened Immunity upon Infection

We first consider the case when $\beta_{p}>\hat{\beta}_{p}$, i.e., the
infection rate for a susceptible individual adopting protection is higher than
the infection rate for a protected recovered individual. That is, initial
infection strengthens the immunity of the individual against subsequent
infections. In this case, we have
$y^{*}_{\mathrm{int}}<\hat{y}^{*}_{\mathrm{int}}$. The reduced epidemic
dynamics for infected and recovered subpopulations, where $z_{S}$ and $z_{R}$
are replaced by their stable equilibrium value in accordance with Proposition
3, is a hybrid system given by

$\displaystyle
y\in[0,y^{*}_{\mathrm{int}}):\begin{cases}&\\!\\!\\!\dot{y}=[\beta_{p}(1-r-y)+\hat{\beta}_{p}r-\gamma]y\\\
&\\!\\!\\!\dot{r}=[-\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (10a)
$\displaystyle
y=y^{*}_{\mathrm{int}}:\begin{cases}&\\!\\!\\!\dot{y}\in\Big{\\{}[(z_{S}+\alpha(1-z_{S})\beta_{p}(1-r-y)+\hat{\beta}_{p}r-\gamma]y\
\Big{|}\ z_{S}\in[0,1]\Big{\\}}\\\
&\\!\\!\\!\dot{r}=[-\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (10b)
$\displaystyle
y\in(y^{*}_{\mathrm{int}},\hat{y}^{*}_{\mathrm{int}}):\begin{cases}&\\!\\!\\!\\!\\!\dot{y}=\\![\alpha\beta_{p}(1\\!-\\!r-\\!y)\\!+\\!\hat{\beta}_{p}r\\!-\\!\gamma]y\\\
&\\!\\!\\!\\!\\!\dot{r}=\\![-\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (10c)
$\displaystyle
y=\hat{y}^{*}_{\mathrm{int}}:\begin{cases}&\\!\\!\\!\\!\dot{y}\in\Big{\\{}[\alpha\beta_{p}(1-r-y)+\hat{\beta}_{p}r(z_{R}+\\!\alpha(1-\\!z_{R}))-\\!\gamma]y\Big{|}\\!\
z_{R}\in[0,1]\Big{\\}}\\\
&\\!\\!\\!\\!\dot{r}\in\Big{\\{}[-\hat{\beta}_{p}(z_{R}+\alpha(1-z_{R})r+\gamma]y\
\Big{|}\ z_{R}\in[0,1]\Big{\\}},\end{cases}$ (10d) $\displaystyle
y\in(\hat{y}^{*}_{\mathrm{int}},1]:\begin{cases}&\\!\\!\\!\\!\\!\dot{y}=\\![\alpha\beta_{p}(1-\\!r-\\!y)\\!+\alpha\hat{\beta}_{p}r-\\!\gamma]y\\\
&\\!\\!\\!\\!\\!\dot{r}=\\![-\alpha\hat{\beta}_{p}r+\gamma]y.\end{cases}$
(10e)

In particular, when $y<y^{*}_{\mathrm{int}}<\hat{y}^{*}_{\mathrm{int}}$, it
follows from the previous subsection that the stable equilibrium of the fast
system (behavioral replicator dynamics (9d), (9e), (9f)) is
$z_{S}=1,z_{I}=0,z_{R}=1$. Setting these values in (9b) and (9c) together with
the observation $s=1-y-r$ yields the disease dynamics stated in
(LABEL:eq:siri_strong_lowy). The dynamics for moderate and high prevalence of
infection stated in (LABEL:eq:siri_strong_mediumy) and
(LABEL:eq:siri_strong_highy) are obtained in an analogous manner. The dynamics
at points of discontinuities are defined via differential inclusions. As
before, it is easy to see that (10) admits a Filippov solution.

We now characterize the equilibria of the dynamics in (10). The following
three types of equilibria are possible:

$\displaystyle\mathbf{E1}$ $\displaystyle=(0,r^{*}),\quad\text{for}\quad
r^{*}\in[0,1],$ $\displaystyle\mathbf{E2}$
$\displaystyle=\left(1-\frac{\gamma}{\hat{\beta}_{p}},\frac{\gamma}{\hat{\beta}_{p}}\right),\qquad\mathbf{E3}=\left(1-\frac{\gamma}{\alpha\hat{\beta}_{P}},\frac{\gamma}{\alpha\hat{\beta}_{P}}\right).$

$\mathbf{E1}$ corresponds to a continuum of infection free equilibria (IFE)
where the proportion of infected population is $0$ while the proportion of
recovered population depends on the initial condition and the value of the
other parameters. This set of equilibria always exist for the dynamics in
(10). $\mathbf{E2}$ and $\mathbf{E3}$ are endemic equilibrium points with a
nonzero proportion of infected population and the susceptible proportion being
$0$. The existence and stability of all these equilibria as well as a sliding
mode solution are established below.

######  Proposition 4 (Equilibria and Stability of (10)).

For the SIRI epidemic under game-theoretic protection stated in (10), the
following statements hold for any $y(0)\in(0,1]$.

  1. 1.

If $\gamma>{\beta}_{p}$, $\mathbf{E2}$ and $\mathbf{E3}$ are not equilibrium
points, all points at the IFE are locally stable, and $y(t)$ decays
monotonically to $0$.

  2. 2.

If $\hat{\beta}_{p}<\gamma<{\beta}_{p}$, $\mathbf{E2}$ and $\mathbf{E3}$ are
not equilibrium points. Points in the IFE with
$r^{*}>\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$ are locally stable,
and unstable otherwise, and $y(t)$ asymptotically converges to $0$.

  3. 3.

If $\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}$,
$\mathbf{E3}$ is not an equilibrium point, points in the IFE are unstable, and
$\mathbf{E2}$ is an asymptotically stable equilibrium point.

  4. 4.

If
$\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$,
then $\mathbf{E2}$ and $\mathbf{E3}$ are not equilibrium points, points in the
IFE are unstable and
$(\hat{y}^{*}_{\mathrm{int}},1-\hat{y}^{*}_{\mathrm{int}})$ acts as a sliding
mode of (10) with $y(t)\to\hat{y}^{*}_{\mathrm{int}}$ as $t\to\infty$.

  5. 5.

If $\gamma<\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$, $\mathbf{E2}$
is not an equilibrium point, $\mathbf{E3}$ is an asymptotically stable
equilibrium point, and points in the IFE are unstable.

The proof is presented in Appendix B and exploits the notion of input-to-state
stability. Note that in a certain parameter regime (Case 4 in the above
proposition), the infection free equilibria are not stable, endemic equilibria
do not exist and infected fraction converges to a sliding mode of the hybrid
dynamics. This sliding mode corresponds to an equilibrium of the original
coupled dynamics (9) given by
$(s=0,y=\hat{y}^{*}_{\mathrm{int}},r=1-\hat{y}^{*}_{\mathrm{int}},z_{S}=0,z_{I}=0,z_{R}=\frac{1}{1-\alpha}\left[\frac{\gamma}{\beta_{p}(1-\hat{y}^{*}_{\mathrm{int}})}-\alpha\right])$.
It is easy to see that for the range of $\gamma$ in Case 4, $z_{R}\in(0,1)$,
i.e., the outcome is an intermediate level of protection adoption by recovered
individuals and an endemic infection level.

###  3.5 Epidemic Evolution under Compromised Immunity

We now consider the case when $\beta_{p}<\hat{\beta}_{p}$, i.e., the infection
rate for a susceptible individual adopting protection is smaller than the
infection rate for a protected recovered individual. That is, initial
infection leads to compromised immunity against future infections. In this
case, we have $y^{*}_{\mathrm{int}}>\hat{y}^{*}_{\mathrm{int}}$. The reduced
epidemic dynamics for infected and recovered subpopulations is given by

$\displaystyle
y\in[0,\hat{y}^{*}_{\mathrm{int}}):\begin{cases}&\\!\\!\\!\dot{y}=[\beta_{p}(1-r-y)+\hat{\beta}_{p}r-\gamma]y\\\
&\\!\\!\\!\dot{r}=[-\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (11a)
$\displaystyle
y=\hat{y}^{*}_{\mathrm{int}}:\begin{cases}&\\!\\!\\!\dot{y}\in\Big{\\{}[\beta_{p}(1-r-y)+(z_{R}+\alpha(1-z_{R}))\hat{\beta}_{p}r-\gamma]y\
\Big{|}\ z_{S}\in[0,1]\Big{\\}}\\\
&\\!\\!\\!\dot{r}\in\Big{\\{}[-\hat{\beta}_{p}r(z_{R}+\alpha(1-z_{R})+\gamma]y\
\Big{|}\ z_{R}\in[0,1]\Big{\\}},\end{cases}$ (11b) $\displaystyle
y\in(\hat{y}^{*}_{\mathrm{int}},y^{*}_{\mathrm{int}})\\!:\\!\begin{cases}&\\!\\!\\!\\!\\!\dot{y}=[\beta_{p}(1-\\!r-\\!y)+\\!\alpha\hat{\beta}_{p}r-\\!\gamma]y\\\
&\\!\\!\\!\\!\\!\dot{r}=[-\alpha\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (11c)
$\displaystyle
y=y^{*}_{\mathrm{int}}\\!:\\!\begin{cases}&\\!\\!\\!\\!\dot{y}\in\Big{\\{}[\alpha(z_{S}+\\!\alpha(1-\\!z_{S}))\beta_{p}(1-\\!r-\\!y)+\alpha\hat{\beta}_{p}r-\gamma]y\
\Big{|}\ z_{S}\in[0,1]\Big{\\}}\\\
&\\!\\!\\!\\!\dot{r}=[-\alpha\hat{\beta}_{p}r+\gamma]y,\end{cases}$ (11d)
$\displaystyle
y\in(y^{*}_{\mathrm{int}},1]:\\!\\!\begin{cases}&\\!\\!\\!\\!\dot{y}=\\![\alpha\beta_{p}(1-r-y)+\\!\alpha\hat{\beta}_{p}r-\\!\gamma]y\\\
&\\!\\!\\!\\!\dot{r}=\\![-\\!\alpha\hat{\beta}_{p}r+\\!\gamma]y.\end{cases}$
(11e)

The above hybrid dynamics is obtained by setting $z_{S}$ and $z_{R}$ values to
their stable equilibrium values in accordance with Proposition 3. The dynamics
at points of discontinuities are defined via differential inclusions as
before.

The equilibria of (11) coincide with those of (10), i.e., the following three
equilibria exist:

$\displaystyle\mathbf{E1}$ $\displaystyle=(0,r^{*}),\quad\text{for}\quad
r^{*}\in[0,1],$ $\displaystyle\mathbf{E2}$
$\displaystyle=\left(1-\frac{\gamma}{\hat{\beta}_{p}},\frac{\gamma}{\hat{\beta}_{p}}\right),\qquad\mathbf{E3}=\left(1-\frac{\gamma}{\alpha\hat{\beta}_{P}},\frac{\gamma}{\alpha\hat{\beta}_{P}}\right).$

The existence and local stability of all these equilibria as well as a sliding
mode solution are established below.

######  Proposition 5 (Equilibria and Stability of (11)).

For the SIRI epidemic under game-theoretic protection stated in (11), the
following statements hold:

  1. 1.

if $\gamma>\hat{\beta}_{p}$, $\mathbf{E2}$ and $\mathbf{E3}$ are not
equilibrium points, all points at the IFE are locally stable, and $y(t)$
decays monotonically to $0$;

  2. 2.

if $\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}$,
$\mathbf{E2}$ is an equilibrium point which is locally stable, and
$\mathbf{E3}$ is not an equilibrium point;

  3. 3.

if
$\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$,
then $\mathbf{E2}$ and $\mathbf{E3}$ are not equilibrium points, and
$(\hat{y}^{*}_{\mathrm{int}},1-\hat{y}^{*}_{\mathrm{int}})$ acts as a sliding
mode of (11);

  4. 4.

if $\gamma<\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$, $\mathbf{E2}$
is not an equilibrium point, $\mathbf{E3}$ is an equilibrium point which is
locally stable.

In addition, for statements $2-4$ above,

  * •

if $\gamma<\beta_{p}$, then all points in the IFE $\mathbf{E1}$ are unstable;

  * •

if $\gamma>\beta_{p}$, then IFE with with
$r^{*}<\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$ are locally stable,
and IFE with $r^{*}>\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$ are
unstable.

######  Remark 3.

Note that when $\gamma>\beta_{p}$ and $\gamma$ satisfies any of the conditions
in statements $2-4$ of Proposition 5, the hybrid dynamics exhibits bistability
where both the IFE and an endemic equilibrium or attractive sliding mode
coexist and are stable. In contrast, when $\beta_{p}>\hat{\beta}_{p}$, we do
not observe coexistence of multiple stable equilibria. Due to the coexistence
of multiple stable equilibria, characterizing the respective regions of
attraction remains a direction for future research.

##  4 Numerical Results

We now provide further insights into the coupled epidemic-replicator dynamics
via numerical simulations.

###  4.1 SIS Epidemic Setting

Figure 3: Trajectories of the coupled SIS epidemic-behavioral dynamics when
the behavioral dynamics evolves on a faster timescale. The right panel shows
the oscillatory trajectories of infected proportion ($y$) and susceptible
proportion that does not adopt protection ($z_{S}$). The oscillations increase
as $\epsilon$ is decreased. The left panel shows trajectories in $z_{S}-\Delta
F$ plane. For $\Delta F<0$ (resp. $\Delta F>0$) $z_{S}=0$ (resp. $z_{S}=1$) is
the stable equilibrium point and is shown in solid blue line. Trajectories
reach an equilibrium point at which $\Delta F=0$.

We first investigate the SIS epidemic setting. We simulate the dynamics in (5)
for $\epsilon\in\\{0.01,0.1,1\\}$. The same model parameters as in Section 2.2
are selected and $\gamma$ is selected as $0.1$ so that $\mathbf{E3}$ is the
stable equilibrium point. Fig. 3 (right panel) shows the time-evolution of $y$
and $z_{S}$. More oscillatory behavior is observed as $\epsilon$ becomes
smaller, i.e., the behavioral dynamics becomes faster. To understand this, we
focus on $z_{S}$ dynamics (3) with

$\Delta F=c_{P}-L(1-\alpha)(\beta_{u}z_{I}(t)+\beta_{p}(1-z_{I}(t)))y(t)$

as a dynamic parameter. To this end, we illustrate $z_{S}$ trajectories in
$z_{S}-\Delta F$ plane in Fig. 3 (left panel).

Figure 4: Evolution of infected proportion in the SIRI model for different
values of recovery rate $\gamma$ and timescale separation parameter $\epsilon$
in accordance with Proposition 4 when initial infection leads to strengthened
immunity.

Recall that if $\Delta F$ is a positive (resp. negative) constant, then
$z_{S}=1$ (resp. $z_{S}=0$) is a stable equilibrium point. Accordingly,
$z_{S}=0$ and $z_{S}=1$ are marked blue and red in Fig. 3 (left panel), when
they are stable and unstable, respectively.

Since behavioral dynamics is fast, $y$ is quasi-stationary and $z_{I}$ very
quickly converges to zero. In Fig. 3 (left panel), the initial fraction of
infected population is sufficiently high such that $\Delta F<0$, then the fast
behavioral dynamics quickly converges to $z_{S}=0$ (the bottom solid blue
line), i.e., every susceptible individual adopts protection. This results in a
decrease in the fraction of infected population and increases $\Delta F$. As
$\Delta F$ becomes positive, $z_{S}=0$ becomes unstable and $z_{S}$ quickly
jumps to $z_{S}=1$ (the top solid blue line), and a similar process repeats
which again drives $\Delta F$ to negative values. This process leads to the
highly oscillatory behavior seen in Fig. 3. Eventually, trajectories converge
such that $\Delta F=0$ and in this regime $z_{S}$ settles to equilibrium value
$z^{*}_{S,\mathrm{int}}\in(0,1)$.

###  4.2 SIRI Epidemic Setting

We first consider the case where initial infection leads to strengthened
immunity. In order to highlight a wide range of transient behavior of the
coupled dynamics, we choose parameter values as
$\beta_{p}=0.3,\hat{\beta}_{p}=0.2,\beta_{u}=0.4,\hat{\beta}_{u}=0.25,L=75,\alpha=0.6,c_{P}=2,c_{\mathtt{IU}}=2,c_{\mathtt{IP}}=1$.
For this set of parameter values, we have $y^{*}_{\mathrm{int}}=0.2222$ and
$\hat{y}^{*}_{\mathrm{int}}=0.3333$. Fig. 4 shows the evolution of infected
proportion $y(t)$ for three different values of the recovery rate $\gamma$
with initial states as $y(0)=0.2,r(0)=0.01,z_{S}(0)=z_{I}(0)=z_{R}(0)=0.5$. In
particular, we compute the trajectories of the dynamics in (9) with $\epsilon$
multiplied to the R.H.S. of (9a), (9b), (9c) via Runge-Kutta fourth order
method with spacing $0.05$. When $\epsilon=1$, epidemic and replicator
dynamics evolve at the same timescale while a smaller value of $\epsilon$
signifies a faster evolution of behavior compared to disease evolution.

When
$\gamma=0.15\in(\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}],\hat{\beta}_{p})$,
it follows from Proposition 4 that points at the infection free equilibria are
not stable, $\mathbf{E3}$ does not exist while $\mathbf{E2}$ is a locally
stable equilibrium. The plot in the left panel of Figure 4 shows indeed $y(t)$
converges to the endemic infection level at $\mathbf{E2}$ in our simulations.
Furthermore, the convergence is not monotonic. When $\gamma=0.1$, $y(t)$
converges to the endemic infection level $\hat{y}^{*}_{\mathrm{int}}$ in
accordance with Case $4$ of Proposition 4. The convergence is in an
oscillatory manner similar to the setting illustrated above for the SIS
epidemic. Finally, when
$\gamma=0.078<\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$, $y(t)$
converges to the endemic infection level at $\mathbf{E3}$.

Figure 5: Evolution of infected proportion and bistable behavior in the SIRI
model for different values of $\gamma$ and $\epsilon$ when initial infection
leads to compromised immunity in accordance with Proposition 5.

Finally, we consider the case where initial infection leads to compromised
immunity. We choose parameter values as
$\hat{\beta}_{p}=0.25,\hat{\beta}_{u}=0.4,\beta_{u}=0.35,L=125,\alpha=0.6,c_{P}=2,c_{\mathtt{IU}}=2,c_{\mathtt{IP}}=1$.
For this set of parameter values, we have $y^{*}_{\mathrm{int}}=0.25$ and
$\hat{y}^{*}_{\mathrm{int}}=0.16$. Recall that when
$\beta_{p}<\hat{\beta}_{p}$ and $\gamma>\beta_{p}$, infection free equilibria
as well as an endemic equilibrium point can be simultaneously locally stable.
Figure 5 shows the evolution of infected proportion $y(t)$ for three different
values $\gamma,\beta_{p}$ and $y(0)$. Initial values of other states are set
as stated above. When $\gamma=0.14$ and $\beta_{p}=0.12$, both the sliding
mode solution $\hat{y}^{*}_{\mathrm{int}}$ and the IFE are locally stable.
From the plot in the left panel of Figure 5, we observe that when $y(0)=0.05$,
$y(t)$ converges to the sliding mode solution for all three $\epsilon$ values.
However, when $y(0)=0.01$, $y(t)$ converges to the IFE for
$\epsilon=0.1,0.025$ indicating that certain equilibria at the IFE are also
locally stable when replicator dynamics is sufficiently fast in accordance
with Proposition 5. Further, when $\epsilon=1$, $y(t)$ converges to
$\hat{y}^{*}_{\mathrm{int}}$ indicating that the relative speed of epidemic
and behavioral dynamics plays a role in the convergence behavior of the
infected proportion. Finally, when $\beta_{p}=0.15>\gamma$, we no longer have
bistable behavior and $y(t)$ converges to $\hat{y}^{*}_{\mathrm{int}}$ for all
values of $\epsilon$ even when $y(0)=0.001$ as shown in the right panel of
Figure 5.

##  5 Conclusion

We proposed and analyzed a novel model that captures the interaction of
epidemic propagation dynamics with human protection adoption behavioral
dynamics. For the coupled SIS epidemic-replicator dynamics, we completely
characterized the equilibrium points, their stability properties, and the
associated bifurcations. For both SIS and SIRI epidemic models, we further
analyzed the coupled dynamics under timescale separation when the behavioral
response is faster than the epidemic evolution, and characterize the
equilibria of the reduced epidemic dynamics which takes the form of a hybrid
dynamical system. The endemic level of infection at the stable equilibrium
strongly depends on the effectiveness of the protection mechanism and the
disease recovery rate. Numerical results show that the coupled dynamics may
lead to highly oscillatory behavior and may take long time to converge to the
stable equilibrium.

This work is one of the first attempts towards understanding the joint
evolution of infectious diseases and human behavioral response. There are
several promising directions for future research.

  * •

We have considered payoff functions defined in terms of the instantaneous risk
of becoming infected. However, individuals often take decisions with a view of
the future risk of infection and are affected by factors such as fatigue or
procrastination. Incorporating multi-stage payoff functions would require
defining a stronger notion of equilibrium (e.g., mean-field equilibrium) and a
different class of evolutionary learning dynamics (in terms of policies)
leading to a challenging research problem.

  * •

The parameters in the payoff function (such as loss due to infection or cost
of being unprotected) are often dynamic; both due to external intervention or
due to evolving opinions about disease severity. Analyzing the joint evolution
of opinion, behavior and disease remains an interesting problem.

  * •

Finally, extending the present analysis to a networked environment remains a
challenging open problem.

We hope our work stimulates further investigations along the above avenues.

##  Appendix A Proofs omitted from Section 2

We first evaluate the entries of the Jacobian matrix of the coupled dynamics
(2), (3) and (4) required to determine the stability of the equilibrium
points. Specifically, we compute

$\displaystyle\frac{\partial f_{y}}{\partial y}$
$\displaystyle=\big{[}(1-\\!y)(z_{S}+\\!\alpha(1-\\!z_{S}))(\beta_{u}z_{I}\\!+\beta_{p}(1-\\!z_{I}))\\!-\\!\gamma\big{]}$
$\displaystyle\qquad-y(z_{S}+\alpha(1-z_{S}))(\beta_{u}z_{I}+\beta_{p}(1-z_{I})),$
$\displaystyle\frac{\partial f_{y}}{\partial z_{S}}$
$\displaystyle=y\big{[}(1-y)(1-\alpha)(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))\big{]},$
$\displaystyle\frac{\partial f_{y}}{\partial z_{I}}$
$\displaystyle=y\big{[}(1-y)(z_{S}+\alpha(1-z_{S}))(\beta_{u}-\beta_{p})\big{]}.$

Similarly,

$\displaystyle\frac{\partial f_{S}}{\partial y}$
$\displaystyle=-z_{S}(1-{z}_{S})L(1-\alpha)(\beta_{u}z_{I}+\beta_{p}(1-z_{I})),$
$\displaystyle\frac{\partial f_{S}}{\partial z_{S}}$
$\displaystyle=(1-2z_{S})\big{[}c_{P}-L(1-\alpha)(\beta_{u}z_{I}+\beta_{p}(1-z_{I}))y\big{]}$
$\displaystyle\frac{\partial f_{S}}{\partial z_{I}}$
$\displaystyle=-z_{S}(1-{z}_{S})L(1-\alpha)(\beta_{u}-\beta_{p})y$

Finally,

$\displaystyle\frac{\partial f_{I}}{\partial y}$
$\displaystyle=0,\quad\frac{\partial f_{I}}{\partial
z_{S}}=0,\quad\frac{\partial f_{I}}{\partial
z_{I}}=(1-2z_{I})(c_{\mathtt{IP}}-c_{\mathtt{IU}}).$

###  A.1 Proof of Proposition 1

###### Proof.

It can be verified that $\mathbf{E0}$ and $\mathbf{E1}$ are always equilibria
of the coupled dynamics. The Jacobian matrix (see above) at $\mathbf{E0}$ and
$\mathbf{E1}$ are

$J_{\mathbf{E0}}=\left[\begin{smallmatrix}\alpha\beta_{p}-\gamma&0&0\\\
0&c_{P}&0\\\ 0&0&c_{\mathtt{IP}}-c_{\mathtt{IU}}\end{smallmatrix}\right],\quad
J_{\mathbf{E1}}=\left[\begin{smallmatrix}\beta_{p}-\gamma&0&0\\\ 0&-c_{P}&0\\\
0&0&c_{\mathtt{IP}}-c_{\mathtt{IU}}\end{smallmatrix}\right].$

$J_{\mathbf{E0}}$ has a positive eigenvalue if $c_{P}>0$. Thus, for any
nonzero cost of adopting protection, $\mathbf{E0}$ is not a stable equilibrium
point of the coupled dynamics.222This is intuitive: in the absence of
infection, it is optimal to not choose costly protective measures. Likewise,
$\mathbf{E1}$ is stable if $\beta_{p}<\gamma$, and unstable, otherwise.333In
other words, $\mathbf{E1}$ is stable if the disease-free equilibrium (DFE) is
stable for the epidemic dynamic (2) when susceptible individuals do not adopt
protection ($z_{S}=1$) and infected individuals do ($z_{I}=0$).

We now analyze equilibrium points where infection is endemic. It can be
verified that $\mathbf{E2}$ exists only when $\beta_{p}>\gamma$ as otherwise
$y^{*}_{u}<0$. The Jacobian matrix at $\mathbf{E2}$ is

$J_{\mathbf{E2}}=\left[\begin{smallmatrix}d_{1}&d_{2}&d_{3}\\\
0&-[c_{P}-L(1-\alpha)\beta_{p}y^{*}_{u}]&0\\\
0&0&c_{\mathtt{IP}}-c_{\mathtt{IU}}\end{smallmatrix}\right],$

where $d_{1}=(1-2y^{*}_{u})\beta_{p}-\gamma$,
$d_{2}=y^{*}_{u}(1-y^{*}_{u})(1-\alpha)\beta_{p}>0$ and
$d_{3}=y^{*}_{u}(1-y^{*}_{u})(\beta_{u}-\beta_{p})>0$. Thus, $J_{\mathbf{E2}}$
is an upper triangular matrix. Furthermore, the first diagonal entry is

$\displaystyle(1-2y^{*}_{u})\beta_{p}-\gamma$
$\displaystyle=(1-2(1-\frac{\gamma}{\beta_{p}}))\beta_{p}-\gamma$
$\displaystyle=(-1+2\frac{\gamma}{\beta_{p}})\beta_{p}-\gamma=\gamma-\beta_{p}<0,$

in the regime where $\mathbf{E2}$ exists. Therefore, $\mathbf{E2}$ is stable
when

$\displaystyle c_{P}$ $\displaystyle>L(1-\alpha)\beta_{p}y^{*}_{u}\iff
y^{*}_{u}<y^{*}_{\mathrm{int}}.$

It can be verified that
$\mathbf{E3}=(y^{*}_{\mathrm{int}},z^{*}_{S,\mathrm{int}},0)$ is an
equilibrium point. We now examine the conditions under which
$y^{*}_{\mathrm{int}}\in(0,1)$ and $z^{*}_{S,\mathrm{int}}\in(0,1)$. By
definition, $y^{*}_{\mathrm{int}}>0$. We now observe that

$\displaystyle z^{*}_{S,\mathrm{int}}>0$
$\displaystyle\iff\frac{\gamma}{\alpha\beta_{p}}>1-y^{*}_{\mathrm{int}}\iff
y^{*}_{\mathrm{int}}>1-\frac{\gamma}{\alpha\beta_{p}},$ $\displaystyle
z^{*}_{S,\mathrm{int}}<1$
$\displaystyle\iff\frac{\gamma}{\beta_{p}}<1-y^{*}_{\mathrm{int}}\iff
y^{*}_{\mathrm{int}}<1-\frac{\gamma}{\beta_{p}}.$

Thus, $\mathbf{E3}$ exists when $y^{*}_{p}<y^{*}_{\mathrm{int}}<y^{*}_{u}$.
Note that the third row of the Jacobian of the dynamics at $\mathbf{E3}$,
$J_{\mathbf{E3}}$, would be $[0\quad 0\quad c_{\mathtt{IP}}-c_{\mathtt{IU}}]$
as before, and as a result, $c_{\mathtt{IP}}-c_{\mathtt{IU}}<0$ would be an
eigenvalue. Thus, we focus on the $2\times 2$ sub-matrix containing the first
two rows and columns of the Jacobian matrix which simplifies to

$\displaystyle\hat{J}_{\mathbf{E3}}=\left[\begin{smallmatrix}\frac{-\gamma
y^{*}_{\mathrm{int}}}{1-y^{*}_{\mathrm{int}}}&d_{4}\\\
-z^{*}_{S,\mathrm{int}}(1-z^{*}_{S,\mathrm{int}})L(1-\alpha)\beta_{p}&0\\\
\end{smallmatrix}\right],$

where
$d_{4}=y^{*}_{\mathrm{int}}(1-y^{*}_{\mathrm{int}})(1-\alpha)\beta_{p}>0$. For
the above matrix, the sum of the eigenvalues is negative and the determinant
is positive, and as a result, all three eigenvalues of $J_{\mathbf{E3}}$ are
negative. Therefore, $\mathbf{E3}$, when it exists, is a stable equilibrium of
the coupled dynamics.

It can be verified that $\mathbf{E4}$ exists when $y^{*}_{p}\in(0,1)$ or
equivalently, when $\gamma<\alpha\beta_{p}$. The Jacobian matrix at
$\mathbf{E4}$ is

$J_{\mathbf{E4}}=\left[\begin{smallmatrix}(1-2y^{*}_{p})\alpha\beta_{p}-\gamma&d_{5}&d_{6}\\\
0&d_{7}&0\\\ 0&0&c_{\mathtt{IP}}-c_{\mathtt{IU}}\end{smallmatrix}\right],$
(12)

where $d_{5}=y^{*}_{p}(1-y^{*}_{p})(1-\alpha)\beta_{p}>0$,
$d_{6}=y^{*}_{p}(1-y^{*}_{p})\alpha(\beta_{u}-\beta_{p})>0$ and
$d_{7}=c_{P}-L(1-\alpha)\beta_{p}y^{*}_{p}$. Thus, $J_{\mathbf{E4}}$ is a
diagonal matrix with the first diagonal entry

$\displaystyle(1-2y^{*}_{p})\alpha\beta_{p}-\gamma$
$\displaystyle=\gamma-\alpha\beta_{p}<0,$

in the regime where $\mathbf{E4}$ exists. Thus, $\mathbf{E4}$ is stable if

$\displaystyle c_{P}$ $\displaystyle<L(1-\alpha)\beta_{p}y^{*}_{p}\iff
y^{*}_{p}>y^{*}_{\mathrm{int}}.$

This concludes the proof. ∎

###  A.2 Proof of Proposition 2

###### Proof.

It follows from Mei et al. (2017) that under the scalar SIS epidemic dynamics

$\dot{y}(t)=[(1-y(t))\beta-\gamma]y(t),$

the infected proportion $y(t)$ monotonically converges to $0$ if
$\gamma\geq\beta$ and converges to $y^{*}=1-\frac{\gamma}{\beta}$, otherwise.
Furthermore, it is easy to see that in the latter case, $|y(t)-y^{*}|$ is
monotonically decreasing.444For $V(y)=(y-y^{*})^{2}$, we have
$\dot{V}(y)=-2\beta yV(y)\leq 0$.

Note that when $y<y^{*}_{\mathrm{int}}$ (respectively,
$y>y^{*}_{\mathrm{int}}$), (6) is analogous to the above dynamics with
$\beta=\beta_{p}$ (respectively, $\beta=\alpha\beta_{p}$). Recall further that
since $\alpha\in(0,1)$, we have $y^{*}_{p}<y^{*}_{u}$. We now analyze the four
cases stated above.

Case 1: $y^{*}_{u}\leq 0$. In this case, we have
$\alpha\beta_{p}<\beta_{p}\leq\gamma$ and as a result, $y=0$ is the only
equilibrium of (6) both when $y>y^{*}_{\mathrm{int}}$ and
$y<y^{*}_{\mathrm{int}}$. Further, $\dot{y}(t)<0$ for any $y(t)>0$ including
at the neighborhood of $y(t)=y^{*}_{\mathrm{int}}$.

Case 2: $0<y^{*}_{u}<y^{*}_{\mathrm{int}}$. We first show that $y(t)$ is
monotonically decreasing when $y(t)>y^{*}_{\mathrm{int}}$. Note that when
$y(t)>y^{*}_{\mathrm{int}}$, (6) resembles the scalar SIS dynamics with
$\beta=\alpha\beta_{p}$. We have two possibilities.

  * a.

$\gamma\geq\alpha\beta_{p}$: It is easy to see that $\dot{y}(t)<0$ for
$y(t)>y^{*}_{\mathrm{int}}$ in (6), and as a result, $y(t)$ is monotonically
decreasing.

  * b.

$\gamma<\alpha\beta_{p}$: From the above discussion, we note that $y(t)$ under
the scalar SIS dynamics with $\beta=\alpha\beta_{p}$ would converge
monotonically to $1-\frac{\gamma}{\alpha\beta_{p}}=y^{*}_{p}$. Since
$\alpha<1$, $y^{*}_{p}<y^{*}_{u}$ which implies
$y^{*}_{p}<y^{*}_{\mathrm{int}}$. Therefore, $y(t)$ is monotonically
decreasing for $y(t)>y^{*}_{\mathrm{int}}$.

We now focus on the dynamics when $y(t)<y^{*}_{\mathrm{int}}$ which resembles
the scalar SIS dynamics with $\beta=\beta_{p}$. Since $0<y^{*}_{u}$, we have
$\beta_{p}>\gamma$. As a result, for any $y(0)\in(0,y^{*}_{\mathrm{int}})$,
$y(t)$ converges monotonically to $y^{*}_{u}$. Since
$y^{*}_{u}<y^{*}_{\mathrm{int}}$, $\dot{y}(t)<0$ in the neighborhood of
$y^{*}_{\mathrm{int}}$. Thus, the claim follows.

Case 3: $y^{*}_{p}<y^{*}_{\mathrm{int}}<y^{*}_{u}$. Following analogous
arguments as above, we observe that $y(t)$ is monotonically decreasing when
$y(t)>y^{*}_{\mathrm{int}}$ as $\max(0,y^{*}_{p})<y^{*}_{\mathrm{int}}$, and
monotonically increasing when $y(t)<y^{*}_{\mathrm{int}}$ since the
corresponding stable equilibrium $y^{*}_{u}>y^{*}_{\mathrm{int}}$. Thus,
$y^{*}_{\mathrm{int}}$ acts as a sliding surface for the dynamics (6) with
$y(t)=y^{*}_{\mathrm{int}}$ being a Filippov solution of (6) with
$0\in\Big{\\{}\big{[}(1-y^{*}_{\mathrm{int}})\beta_{p}(z_{S}+\alpha(1-z_{S}))-\gamma\big{]}y^{*}_{\mathrm{int}}\;\Big{|}\;z_{S}\in[0,1]\Big{\\}}$.

Case 4: $y^{*}_{p}\geq y^{*}_{\mathrm{int}}$. In this case, $y(t)$ is
monotonically increasing when $y(t)<y^{*}_{\mathrm{int}}$ since
$y^{*}_{u}>y^{*}_{p}>y^{*}_{\mathrm{int}}$. When $y>y^{*}_{\mathrm{int}}$, the
corresponding stable equilibrium is $y^{*}_{p}\in(y^{*}_{\mathrm{int}},1)$. As
a result, $y(t)$ converges monotonically to $y^{*}_{p}$ for
$y(t)>y^{*}_{\mathrm{int}}$. This concludes the proof. ∎

##  Appendix B Proofs omitted from Section 3

###  B.1 Proof of Proposition 4

###### Proof.

We prove each of the cases below.

Case 1: $\gamma>{\beta}_{p}$. Since $\alpha\in(0,1)$ together with
$\gamma>\beta_{p}>\hat{\beta}_{p}$, we have
$1-\frac{\gamma}{\alpha\hat{\beta}_{p}}<1-\frac{\gamma}{\hat{\beta}_{p}}<0$.
In other words, the endemic infection level is negative at $\mathbf{E2}$ and
$\mathbf{E3}$. As a result, $\mathbf{E2}$ and $\mathbf{E3}$ are not
equilibrium points. The local stability of the IFE depends on the dynamics
(LABEL:eq:siri_strong_lowy). Drawing analogy with the SIRI dynamics in (7), we
have $R_{0}=\frac{\beta_{p}}{\gamma}<1$ and
$R_{1}=\frac{\hat{\beta}_{p}}{\gamma}<1$ for the dynamics
(LABEL:eq:siri_strong_lowy). Since $\alpha<1$, we have $R_{0}<1$ and $R_{1}<1$
for the dynamics (LABEL:eq:siri_strong_mediumy) and
(LABEL:eq:siri_strong_highy) as well. As a result, following Case 1 of Theorem
1, all points in the IFE of (10) are locally stable and $y(t)$ decays
monotonically to $0$.

Case 2: $\hat{\beta}_{p}<\gamma<{\beta}_{p}$. $\mathbf{E2}$ and $\mathbf{E3}$
not being equilibrium points follows from the arguments in Case 1 above. We
have $R_{0}=\frac{\beta_{p}}{\gamma}>1$ and
$R_{1}=\frac{\hat{\beta}_{p}}{\gamma}<1$ for the dynamics
(LABEL:eq:siri_strong_lowy). Following Case 3 of Theorem 1, points in the IFE
with $s^{*}<\frac{1-R_{1}}{R_{0}-R_{1}}$ or equivalently with
$r^{*}>1-\frac{1-R_{1}}{R_{0}-R_{1}}=\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$
are locally stable, and unstable otherwise.

We now argue that any trajectory of (10) with $y(0)\in(0,1]$ converges to a
stable IFE. Observe that at any IFE point $(y=0,r=\bar{r})$, the Jacobian
matrix is given by

$J_{IFE}=\left[\begin{smallmatrix}(\hat{\beta}_{p}-\beta_{p})\bar{r}+\beta_{p}-\gamma&0\\\
-\beta_{p}+\gamma&0\\\ \end{smallmatrix}\right],$

i.e., it has $0$ as one eigenvalue and the other eigenvalue
$(\hat{\beta}_{p}-\beta_{p})\bar{r}+\beta_{p}-\gamma>0$ for
$\bar{r}\in[0,\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}})$. Thus, the
vector field at the IFE points towards increasing value of $y(t)$. Let
$y(0)\neq 0$ and let $t_{0},\epsilon_{y}>0$ be suitable constants such that
$y(t)\geq\epsilon_{y}$ for all $t\geq t_{0}$. Consequently, there exists $c>0$
such that under the dynamics (10),

$\dot{s}=-(\dot{y}+\dot{r})\leq-c\beta_{p}\epsilon_{y}s\implies s(t)\leq
s(0)e^{-c\beta_{p}\epsilon_{y}t}.$

As a result, for any $\epsilon_{s}>0$, there exists $T_{0}$ such that
$s(t)\leq\epsilon_{s}$ for $t>T_{0}$. In other words, $s(t)$ acts as a
vanishing input to the dynamics of infected proportion. We now exploit this
property to prove convergence to a stable IFE.

First we show that if $y(t)>\hat{y}^{*}_{\mathrm{int}}$, infected proportion
eventually decreases. Indeed, we have

$\displaystyle\dot{y}$
$\displaystyle=[\alpha\beta_{p}s+\alpha\hat{\beta}_{p}(1-s-y)-\gamma]y$
$\displaystyle=[\alpha(\beta_{p}-\hat{\beta}_{p})s+\alpha\hat{\beta}_{p}(1-y)-\gamma]y$
$\displaystyle\leq[\alpha(\beta_{p}-\hat{\beta}_{p})\epsilon_{s}+\alpha\hat{\beta}_{p}(1-y)-\gamma]y$
$\displaystyle<[\alpha(\beta_{p}-\hat{\beta}_{p})\epsilon_{s}+\hat{\beta}_{p}-\gamma]y;$

the last step follows since $\alpha\in(0,1)$ and $y\leq 1$. Since the second
term is strictly negative, there exists an $\epsilon_{s}$ and $T_{0}$ such
that $\dot{y}(t)<0$ for $t>T_{0}$ when $y>\hat{y}^{*}_{\mathrm{int}}$, and as
a result, the trajectory will eventually remain confined to
(LABEL:eq:siri_strong_lowy) and (LABEL:eq:siri_strong_mediumy).

Next we show that if
$y(t)\in({y}^{*}_{\mathrm{int}},\hat{y}^{*}_{\mathrm{int}})$, infected
proportion eventually decreases as well. In this regime, we have

$\displaystyle\dot{y}$
$\displaystyle=[\alpha\beta_{p}s+\hat{\beta}_{p}(1-s-y)-\gamma]y$
$\displaystyle=[(\alpha\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}(1-y)-\gamma]y$
$\displaystyle<[(\alpha\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}-\gamma]y.$

Note that if $\alpha\beta_{p}<\hat{\beta}_{p}$, we have $\dot{y}<0$ for all
$s(t)\in[0,1]$. If $\alpha\beta_{p}>\hat{\beta}_{p}$, following identical
arguments as above, then there exists an $\epsilon^{\prime}_{s}$ and
$T^{\prime}_{0}$ such that $\dot{y}(t)<0$ for $t>T^{\prime}_{0}$ since
$\hat{\beta}_{p}-\gamma<0$.

It remains to show that when the trajectory is eventually confined to
(LABEL:eq:siri_strong_lowy), i.e., $y(t)\in(0,{y}^{*}_{\mathrm{int}})$,
infected proportion eventually converges to an IFE. In this regime, we have

$\displaystyle\dot{y}$
$\displaystyle=[\beta_{p}s+\hat{\beta}_{p}(1-s-y)-\gamma]y$
$\displaystyle<[(\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}-\gamma]y,$

where the term $\hat{\beta}_{p}-\gamma<0$. Thus, following an identical
argument as above, we have $\dot{y}(t)<0$ for some $t>T^{\prime\prime}_{0}$.
Therefore, $y(t)$ asymptotically converges to an IFE.

Case 3:
$\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}$. For
the dynamics (LABEL:eq:siri_strong_lowy), we have
$R_{0}=\frac{\beta_{p}}{\gamma}>1$ and
$R_{1}=\frac{\hat{\beta}_{p}}{\gamma}>1$ in this regime. Following Theorem 1,
we conclude that all points in the IFE are unstable.

Since $\mathbf{E3}$ is the endemic equilibrium for
(LABEL:eq:siri_strong_highy), its existence as an endemic equilibrium for (10)
requires

$1-\frac{\gamma}{\alpha\hat{\beta}_{p}}>\hat{y}^{*}_{\mathrm{int}}\iff\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]>\gamma,$

which is not satisfied in this parameter regime as $\alpha<1$.

Now observe that $\mathbf{E2}$ is the endemic equilibrium for both
(LABEL:eq:siri_strong_lowy) and (LABEL:eq:siri_strong_mediumy), but not for
the dynamics (LABEL:eq:siri_strong_highy). Thus, for $\mathbf{E2}$ to be an
equilibrium for the hybrid system (10), we must have

$0<1-\frac{\gamma}{\hat{\beta}_{p}}<\hat{y}^{*}_{\mathrm{int}},$

or equivalently, $\hat{\beta}_{p}>\gamma$ and
$\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma$, i.e., precisely the
parameter regime in this case. Further, for both (LABEL:eq:siri_strong_lowy)
and (LABEL:eq:siri_strong_mediumy), $R_{1}=\frac{\hat{\beta}_{p}}{\gamma}>1$,
and as a result from Theorem 1, $\mathbf{E2}$ is locally stable. Let
$y_{E2}=1-\frac{\gamma}{\hat{\beta}_{p}}$.

It remains to show that the infected proportion $y(t)$ converge to
$\mathbf{E2}$. Following analogous arguments as Case 2 above, we can show that
there exists $T_{0}$ such that $\dot{y}(t)<0$ for $t>T_{0}$,
$y(t)>\hat{y}^{*}_{\mathrm{int}}$. In other words, the trajectory will
eventually remain confined to (LABEL:eq:siri_strong_lowy) and
(LABEL:eq:siri_strong_mediumy). We now consider the following two sub-cases.

Case 3A: $y_{E2}\in(0,{y}^{*}_{\mathrm{int}})$. Following analogous arguments
as Case 2 above, we first claim that there exists $T^{\prime}_{0}$ such that
$\dot{y}(t)<0$ for $t>T^{\prime}_{0}$,
$y(t)\in({y}^{*}_{\mathrm{int}},\hat{y}^{*}_{\mathrm{int}})$. In other words,
the trajectory eventually remains confined to (LABEL:eq:siri_strong_lowy).
Observe now that the dynamics of infected proportion in
(LABEL:eq:siri_strong_lowy) can be viewed as a perturbed system with
$s(t)=1-r(t)-y(t)$ playing the role of a vanishing input. We now show that the
equilibrium point $\mathbf{E2}$ is input-to-state stable (ISS) Khalil (2002).
Let $V(y)=(y-y_{E2})^{2}$ be the candidate ISS-Lyapunov function. Under the
dynamics (LABEL:eq:siri_strong_lowy), we have

$\displaystyle\dot{V}(y)$
$\displaystyle=2(y-y_{E2})[\beta_{p}s+\hat{\beta}_{p}r-\gamma]y$
$\displaystyle=2(y-y_{E2})[\beta_{p}s+\hat{\beta}_{p}(1-s-y)-\hat{\beta}_{p}(1-y_{E2})]y$
$\displaystyle=2(y-y_{E2})[(\beta_{p}-\hat{\beta}_{p})s-\hat{\beta}_{p}(y-y_{E2})]y$
$\displaystyle=2(y-y_{E2})(\beta_{p}-\hat{\beta}_{p})ys-2\hat{\beta}_{p}yV(y)$
$\displaystyle\leq-2(\hat{\beta}_{p}-\delta)yV(y),$

when $|y(t)-y_{E2}|\geq\delta^{-1}(\beta_{p}-\hat{\beta}_{p})|s(t)|$. It
follows from Khalil (2002) that $y_{E2}$ is ISS for the dynamics of infected
proportion when $y(t)\in[0,{y}^{*}_{\mathrm{int}})$. Since $s(t)\to 0$ as
$t\to\infty$, we have $y(t)\to y_{E2}$.

Case 3B: $y_{E2}\in({y}^{*}_{\mathrm{int}},\hat{y}^{*}_{\mathrm{int}})$. First
observe that if $y(t)\in(0,{y}^{*}_{\mathrm{int}})$, we have

$\displaystyle\dot{y}$
$\displaystyle=[\beta_{p}s+\hat{\beta}_{p}(1-s-y)-\gamma]y$
$\displaystyle=[(\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}(1-y)-\gamma]y$
$\displaystyle>[(\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}(1-{y}^{*}_{\mathrm{int}})-\hat{\beta}_{p}(1-y_{E2})]y$
$\displaystyle=[(\beta_{p}-\hat{\beta}_{p})s+\hat{\beta}_{p}(y_{E2}-{y}^{*}_{\mathrm{int}})]y>0,$

since $y_{E2}>{y}^{*}_{\mathrm{int}}$ in this sub-case. Consequently, $y(t)$
remains confined to (LABEL:eq:siri_strong_mediumy). Let $V(y)=(y-y_{E2})^{2}$
as before. Under the dynamics (LABEL:eq:siri_strong_mediumy), we have

$\displaystyle\dot{V}(y)$
$\displaystyle=2(y-y_{E2})\\![\alpha\beta_{p}s+\hat{\beta}_{p}(1-\\!s-\\!y)-\\!\hat{\beta}_{p}(1-\\!y_{E2})]y$
$\displaystyle=2(y-y_{E2})[(\alpha\beta_{p}-\hat{\beta}_{p})s-\hat{\beta}_{p}(y-y_{E2})]y$
$\displaystyle=2(y-y_{E2})(\alpha\beta_{p}-\hat{\beta}_{p})ys-2\hat{\beta}_{p}yV(y)$
$\displaystyle\leq-2(\hat{\beta}_{p}-\delta)yV(y)$

when $|y(t)-y_{E2}|\geq\delta^{-1}|\alpha\beta_{p}-\hat{\beta}_{p}|s(t)$ which
implies that $y_{E2}$ is ISS for the dynamics of infected proportion when
$y(t)\in({y}^{*}_{\mathrm{int}},\hat{y}^{*}_{\mathrm{int}})$. Since $s(t)\to
0$ as $t\to\infty$, we have $y(t)\to y_{E2}$.

Case 4:
$\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma<\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$.
It follows from Case 3 above that for $\mathbf{E2}$ to be an equilibrium of
the hybrid dynamics, we must have
$\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]<\gamma$ which is not satisfied
in this regime. Similarly, for $\mathbf{E3}$ to be an equilibrium, we must
have $\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]>\gamma$ which is not
satisfied in this regime. IFE being unstable follows from previous arguments
as $R_{0}>1$ and $R_{1}>1$ for the dynamics (LABEL:eq:siri_strong_lowy).

We now focus on the dynamics at the neighborhood of
$(\hat{y}^{*}_{\mathrm{int}},1-\hat{y}^{*}_{\mathrm{int}})$. For
$y=\hat{y}^{*}_{\mathrm{int}}+\epsilon,r=1-\hat{y}^{*}_{\mathrm{int}}-\epsilon$
with a sufficiently small $\epsilon>0$, we have

$\displaystyle\dot{r}$
$\displaystyle=[-\alpha\hat{\beta}_{p}(1-\hat{y}^{*}_{\mathrm{int}}-\epsilon)+\gamma](\hat{y}^{*}_{\mathrm{int}}+\epsilon)>0,$
$\displaystyle\dot{y}$
$\displaystyle=[\alpha\hat{\beta}_{p}(1-\hat{y}^{*}_{\mathrm{int}}-\epsilon)-\gamma](\hat{y}^{*}_{\mathrm{int}}+\epsilon)<0,$

since $\gamma>\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$. Similarly,
at
$y=\hat{y}^{*}_{\mathrm{int}}-\epsilon,r=1-\hat{y}^{*}_{\mathrm{int}}+\epsilon$
with a sufficiently small $\epsilon>0$, we have

$\displaystyle\dot{r}$
$\displaystyle=[-\hat{\beta}_{p}(1-\hat{y}^{*}_{\mathrm{int}}+\epsilon)+\gamma](\hat{y}^{*}_{\mathrm{int}}-\epsilon)<0,$
$\displaystyle\dot{y}$
$\displaystyle=[\hat{\beta}_{p}(1-\hat{y}^{*}_{\mathrm{int}}+\epsilon)-\gamma](\hat{y}^{*}_{\mathrm{int}}-\epsilon)>0,$

since $\gamma<\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$. Since the
differential inclusion in (LABEL:eq:siri_strong_th2) is a convex combination
of the dynamics in (LABEL:eq:siri_strong_mediumy) and
(LABEL:eq:siri_strong_highy),
$(\hat{y}^{*}_{\mathrm{int}},1-\hat{y}^{*}_{\mathrm{int}})$ acts as a sliding
mode of (10).

Following identical arguments as Case 3B above, it can be shown that when
$y(t)\in(0,\hat{y}^{*}_{\mathrm{int}})$, $\dot{y}>0$ and when
$y(t)>\hat{y}^{*}_{\mathrm{int}}$, $\dot{y}<0$ for some $t>T_{0}$.
Consequently, $y(t)\to{y}^{*}_{\mathrm{int}}$ as $t\to\infty$.

Case 5: $\gamma<\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$. It
follows from Case 3 that for $\mathbf{E3}$ to be an equilibrium of the hybrid
dynamics, we require
$\gamma<\alpha\hat{\beta}_{p}[1-\hat{y}^{*}_{\mathrm{int}}]$ and
$\gamma<\alpha\hat{\beta}_{p}$, both of which are satisfied in this regime.
Further, for (LABEL:eq:siri_strong_highy),
$R_{1}=\frac{\alpha\hat{\beta}_{p}}{\gamma}>1$, and consequently,
$\mathbf{E3}$ is locally stable. The points in IFE being unstable and
$\mathbf{E2}$ not being an equilibrium point follows from previous arguments.

Following identical arguments as Case 3B above, it can be shown that when
$y(t)<\hat{y}^{*}_{\mathrm{int}}$, $\dot{y}>0$ and when
$y(t)>\hat{y}^{*}_{\mathrm{int}}$, the endemic infection level at
$\mathbf{E3}$ is ISS for the infection dynamics with $s(t)$ being the
vanishing input. This concludes the proof. ∎

###  B.2 Proof of Proposition 5

###### Proof.

The proof of statements $1-4$ follows from analogous arguments as the proof of
Proposition 4 and is omitted in the interest of space. We now focus on
establishing the claims regarding the stability of IFE. Note that local
stability of the IFE depends on the dynamics (11a). For statements $2-4$, we
have $\gamma<\hat{\beta}_{p}$, and as a result,
$R_{1}=\frac{\hat{\beta}_{p}}{\gamma}>1$ for (11a).

When $\gamma<\beta_{p}$, we have $R_{0}=\frac{\beta_{p}}{\gamma}>1$ for the
dynamics (11a). Following Case 3 of Theorem 1, all points in the IFE are
unstable.

When $\gamma>\beta_{p}$, we have $R_{0}=\frac{\beta_{p}}{\gamma}<1$ for the
dynamics (11a). It follows from Case 4 of Theorem 1 that IFE with

$\displaystyle
s^{*}>\frac{1-\frac{\hat{\beta}_{p}}{\gamma}}{\frac{\beta_{p}}{\gamma}-\frac{\hat{\beta}_{p}}{\gamma}}=\frac{\gamma-\hat{\beta}_{p}}{\beta_{p}-\hat{\beta}_{p}}\iff
r^{*}<\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$

are locally stable and IFE with
$r^{*}>\frac{\beta_{p}-\gamma}{\beta_{p}-\hat{\beta}_{p}}$ are unstable. ∎

## References

  * Al-Radhawi et al. [2021] M Ali Al-Radhawi, Mahdiar Sadeghi, and Eduardo D Sontag.  Long-term regulation of prolonged epidemic outbreaks in large populations via adaptive control: A singular perturbation approach.  _IEEE Control Systems Letters_ , 2021.
  * Arcak and Martins [2020] Murat Arcak and Nuno C Martins.  Dissipativity tools for convergence to Nash equilibria in population games.  _IEEE Transactions on Control of Network Systems_ , 8(1):39–50, 2020.
  * Berglund and Gentz [2006] Nils Berglund and Barbara Gentz.  _Noise-Induced Phenomena in Slow-Fast Dynamical Systems: A Sample-Paths Approach_.  Springer Science & Business Media, 2006.
  * Blanchini and Miani [2008] Franco Blanchini and Stefano Miani.  _Set-Theoretic Methods in Control_.  Springer, 2008.
  * Chang et al. [2020] Sheryl L Chang, Mahendra Piraveenan, Philippa Pattison, and Mikhail Prokopenko.  Game theoretic modelling of infectious disease dynamics and intervention methods: A review.  _Journal of Biological Dynamics_ , 14(1):57–89, 2020.
  * Cortes [2008] Jorge Cortes.  Discontinuous dynamical systems.  _IEEE Control Systems Magazine_ , 28(3):36–73, 2008.
  * Crellen et al. [2021] Thomas Crellen, Li Pi, Emma L Davis, Timothy M Pollington, Tim CD Lucas, Diepreye Ayabina, Anna Borlase, Jaspreet Toor, Kiesha Prem, Graham F Medley, Petra Klepac, and T. Déirdre Hollingsworth.  Dynamics of SARS-CoV-2 with waning immunity in the UK population.  _Philosophical Transactions of the Royal Society B_ , 376(1829):20200274, 2021.
  * Daskalakis et al. [2009] Constantinos Daskalakis, Paul W Goldberg, and Christos H Papadimitriou.  The complexity of computing a Nash equilibrium.  _SIAM Journal on Computing_ , 39(1):195–259, 2009\.
  * Dhooge et al. [2003] Annick Dhooge, Willy Govaerts, and Yu A Kuznetsov.  MATCONT: A MATLAB package for numerical bifurcation analysis of ODEs.  _ACM Transactions on Mathematical Software_ , 29(2):141–164, 2003.
  * Eksin et al. [2017] Ceyhun Eksin, Jeff S Shamma, and Joshua S Weitz.  Disease dynamics on a network game: A little empathy goes a long way.  _Scientific Reports_ , 7:44122, 2017.
  * Elokda et al. [2021a] Ezzat Elokda, Saverio Bolognani, and Ashish R Hota.  A dynamic population model of strategic interaction and migration under epidemic risk.  In _60th IEEE Conference on Decision and Control (CDC)_ , pages 2085 – 2091, 2021a.
  * Elokda et al. [2021b] Ezzat Elokda, Andrea Censi, and Saverio Bolognani.  Dynamic population games.  _arXiv:2104.14662_ , 2021b.
  * Frieswijk et al. [2022] Kathinka Frieswijk, Lorenzo Zino, Mengbin Ye, Alessandro Rizzo, and Ming Cao.  A mean-field analysis of a network behavioral–epidemic model.  _IEEE Control Systems Letters_ , 6:2533–2538, 2022.
  * Gomes et al. [2004] M Gabriela M Gomes, Lisa J White, and Graham F Medley.  Infection, reinfection, and vaccination under suboptimal immune protection: Epidemiological perspectives.  _Journal of Theoretical Biology_ , 228(4):539–549, 2004.
  * Gong and Cao [2021] Lulu Gong and Ming Cao.  Different environment feedback in fast-slow eco-evolutionary dynamics.  _arXiv preprint arXiv:2105.04659_ , 2021.
  * Hofbauer and Sigmund [2003] Josef Hofbauer and Karl Sigmund.  Evolutionary game dynamics.  _Bulletin of the American Mathematical Society_ , 40(4):479–519, 2003.
  * Hota and Sundaram [2019] Ashish R Hota and Shreyas Sundaram.  Game-theoretic vaccination against networked SIS epidemics and impacts of human decision-making.  _IEEE Transactions on Control of Network Systems_ , 6(4):1461–1472, 2019.
  * Hota et al. [2022] Ashish R Hota, Tanya Sneh, and Kavish Gupta.  Impacts of game-theoretic activation on epidemic spread over dynamical networks.  _SIAM Journal on Control and Optimization_ , pages S92–S118, 2022\.
  * Huang and Zhu [2020] Yunhan Huang and Quanyan Zhu.  A differential game approach to decentralized virus-resistant weight adaptation policy over complex networks.  _IEEE Transactions on Control of Network Systems_ , 7(2):944–955, 2020.
  * Huang and Zhu [2021] Yunhan Huang and Quanyan Zhu.  Game-theoretic frameworks for epidemic spreading and human decision making: A review.  _arXiv preprint arXiv:2106.00214_ , 2021.
  * Khalil [2002] H K Khalil.  _Nonlinear Systems_.  Patience Hall, 2002.
  * Khazaei et al. [2021] Hossein Khazaei, Keith Paarporn, Alfredo Garcia, and Ceyhun Eksin.  Disease spread coupled with evolutionary social distancing dynamics can lead to growing oscillations.  In _60th IEEE Conference on Decision and Control (CDC)_ , pages 4280 – 4286. IEEE, 2021.
  * Marshall [2020] Michael Marshall.  The lasting misery of coronavirus long-haulers.  _Nature_ , 585(7825):339–342, 2020.
  * Martins et al. [2022] Nuno C Martins, Jair Certorio, and Richard J La.  Epidemic population games and evolutionary dynamics.  _arXiv preprint arXiv:2201.10529_ , 2022.
  * Mei et al. [2017] Wenjun Mei, Shadi Mohagheghi, Sandro Zampieri, and Francesco Bullo.  On the dynamics of deterministic epidemic propagation over networks.  _Annual Reviews in Control_ , 44:116–128, 2017.
  * Nowzari et al. [2016] Cameron Nowzari, Victor M Preciado, and George J Pappas.  Analysis and control of epidemics: A survey of spreading processes on complex networks.  _IEEE Control Systems_ , 36(1):26–46, 2016.
  * Pagliara and Leonard [2020] Renato Pagliara and Naomi Ehrich Leonard.  Adaptive susceptibility and heterogeneity in contagion models on networks.  _IEEE Transactions on Automatic Control_ , 66(2):581–594, 2020.
  * Pagliara et al. [2018] Renato Pagliara, Biswadip Dey, and Naomi Ehrich Leonard.  Bistability and resurgent epidemics in reinfection models.  _IEEE Control Systems Letters_ , 2(2):290–295, 2018.
  * Park et al. [2019] Shinkyu Park, Nuno C Martins, and Jeff S Shamma.  Payoff dynamics model and evolutionary dynamics model: Feedback and convergence to equilibria.  _arXiv preprint arXiv:1903.02018_ , 2019.
  * Sandholm [2010] William H Sandholm.  _Population Games and Evolutionary Dynamics_.  MIT Press, 2010.
  * Satapathi et al. [2022] Abhisek Satapathi, Narendra Kumar Dhar, Ashish R. Hota, and Vaibhav Srivastava.  Epidemic propagation under evolutionary behavioral dynamics: Stability and bifurcation analysis.  In _American Control Conference (ACC)_ , pages 3662 – 3667. IEEE, 2022.
  * Theodorakopoulos et al. [2013] George Theodorakopoulos, Jean-Yves Le Boudec, and John S Baras.  Selfish response to epidemic propagation.  _IEEE Transactions on Automatic Control_ , 58(2):363–376, 2013.
  * Trajanovski et al. [2015] Stojan Trajanovski, Yezekael Hayel, Eitan Altman, Huijuan Wang, and Piet Van Mieghem.  Decentralized protection strategies against SIS epidemics in networks.  _IEEE Transactions on Control of Network Systems_ , 2(4):406–419, 2015.
  * Weitz et al. [2016] Joshua S Weitz, Ceyhun Eksin, Keith Paarporn, Sam P Brown, and William C Ratcliff.  An oscillating tragedy of the commons in replicator dynamics with game-environment feedback.  _Proceedings of the National Academy of Sciences_ , 113(47):E7518–E7525, 2016.
  * Zhou et al. [2020] Yunxiu Zhou, Simon A Levin, and Naomi Ehrich Leonard.  Active control and sustained oscillations in actSIS epidemic dynamics.  In _3rd IFAC Workshop on Cyber-Physical & Human Systems_, pages 807–812, Beijing, China, 2020. Elsevier.