\documentclass[reqno]{amsart} \usepackage{hyperref} \usepackage{graphicx} \AtBeginDocument{{\noindent\small \emph{Electronic Journal of Differential Equations}, Vol. 2012 (2012), No. 233, pp. 1--7.\newline ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu \newline ftp ejde.math.txstate.edu} \thanks{\copyright 2012 Texas State University - San Marcos.} \vspace{9mm}} \begin{document} \title[\hfilneg EJDE-2012/233\hfil A Lie algebra approach] {A Lie algebra approach to susceptible-infected-susceptible epidemics} \author[Y. Shang\hfil EJDE-2012/233\hfilneg] {Yilun Shang} % in alphabetical order \address{Yilun Shang \newline Institute for Cyber Security, University of Texas at San Antonio\\ San Antonio, Texas 78249, USA} \email{shylmath@hotmail.com} \thanks{Submitted August 24, 2012. Published December 21, 2012.} \subjclass[2000]{92D30, 17B80, 60J22} \keywords{Epidemic dynamics; Lie algebra; Riccati equation; \hfill\break\indent susceptible-infected-susceptible} \begin{abstract} The susceptible-infected-susceptible (SIS) epidemic model can be represented by a continuous-time Markov chain, which is governed by a set of deterministic differential equations (Kolmogorov forward equations). In this paper, a Lie algebra approach is applied to solve an SIS model where infection rate and recovery rate are time-varying. The method presented here has been used widely in chemical and physical sciences but not in epidemic applications due to insufficient symmetries. \end{abstract} \maketitle \numberwithin{equation}{section} \allowdisplaybreaks \section{Introduction} Analytical description of epidemic spreading has a long history and can be traced back to the seminal work of Kermack and McKendrick \cite{r1,r2}, where only three simple ordinary differential equations are used following the mass action assumption; i.e., the rate of increase in epidemic incidence is proportional to the product of the number of infectious and susceptible individuals. It is also possible to capture the propagation phenomena by mean-field theory \cite{r3,r4,r5,r27} or generating function formalism \cite{r7,r8,r15} especially when the host population is modeled by a network. Such methods, however, are generally more accurate (and valid in essence) when the population size is relatively large. Recently, Keeling and Ross \cite{r10} proposed a time homogeneous Markov chain model to characterize the stochastic nature of epidemic spreading. The complete ensemble of behavior can be predicted by $N+1$ differential equations for susceptible-infected-susceptible (SIS) dynamics \cite{r23} by virtue of the Kolmogorov forward equation \cite{r11}, which governs the rates of transition between states of the disease. The solution of the system can be expressed by the form of matrix exponentials \cite{r10,r12}. Indeed, continuous-time Markovian models are shown to be powerful tools to study stochastic evolutionary processes and have been widely used in other biological and metapopulation models \cite{r14,r13,r16,r25}. In this paper, we further investigate the SIS paradigm represented by a time inhomogeneous Markov chain. In this model there is a fixed population of size $N$, where $S(t)$ and $I(t)$ represent the number of susceptibles and infectives, respectively, in the population at time $t$, $t\ge0$, and $S(t)+I(t)=N$. No immunity is conferred upon recovery from infection, and recovered individuals return immediately to the susceptible state. Members of the population transmit the infection immediately upon becoming infected. Based on the Lie algebraic method developed in \cite{r17}, we generate a low-dimensional Lie algebra and solve the Markovian model by deriving matrix exponential solutions. Different from many physical or chemical systems \cite{r19,r20}, biological or epidemic models often lack of symmetries, which adds difficulty in finding a proper Lie algebra. It is worth noting that Lie algebra solution of some birth-and-death type population models is recently established by House \cite{r18}. The rest of the paper is organized as follows. In Section 2, we briefly review the Lie algebraic methodology for continuous-time Markov chains. we then apply it to an SIS epidemic model in Section 3, and conclude the paper in Section 4. \section{Lie algebra solution of time inhomogeneous Markov chains} In algebra theory, a Lie algebra \cite{r21} is a vector space $V$ over some field $F$ together with a bilinear map $[\cdot,\cdot]:V\times V\rightarrow V$ called the Lie bracket, which satisfies $[X,X]=0$ and the Jacobi identity $$[X,[Y,Z]]+[Y,[Z,X]]+[Z,[X,Y]]=0,\label{e1}$$ for all $X,Y,Z\in V$. For $X\in V$, we define an adjoint operator $\text{ad}X$ by $$(\text{ad}X)Y=[X,Y],\label{e2}$$ for $Y\in V$. In doing so, multiple Lie brackets can be expressed in a succinct way; e.g., $(\text{ad}X)^2Y=[X,[X,Y]]$, etc. Every associate algebra gives rise to a Lie algebra $V$ by defining the Lie bracket as a commutator $$[X,Y]=XY-YX,\label{e3}$$ where $X,Y\in V$. In what follows, we will focus on this Lie product. The classical Baker-Campbell-Hausdorff formula can be written in terms of \eqref{e2} as $$e^XYe^{-X}=(e^{\text{ad}X})Y,\label{e4}$$ where $e^X=\sum_{i=0}^{\infty}X^i/i!$. The type of processes we consider here are continuous-time Markov chains \cite{r10,r11}, taking values in a finite or countably infinite state space $\mathcal{S}$. The dynamical behavior of the Markov chain is specified by a matrix $Q(t)=(q_{ij}(t),i,j\in\mathcal{S})$, where $q_{ij}(t)$ is the rate of transition from state $i$ to state $j$, for $j\not=i$, and $-q_{ii}(t)=q_i(t)=\sum_{j\not=i}q_{ij}(t)$ is the total rate at which we move out of state $i$ at time $t$. In light of the Kolmogorov forward equation (also called the ensemble or master equation), the probability distribution of the process at time $t$, $p(t)=(p_i(t),i\in\mathcal{S})$, is given by $$\frac{{\rm d}p(t)}{{\rm d}t}=H(t)p(t),\label{e5}$$ where $H(t)=Q(t)^T$ ($T$ means transpose), and $p(t)$ is a column probability vector with component $p_i(t)$ representing the probability of finding the system in state $i$ at time $t$. Making use of the Dirac notation for vectors (kets $|\cdot\rangle$), the probability vector can alternatively be written as $$|p(t)\rangle=\sum_{i\in\mathcal{S}}P(i|t)|i\rangle,\label{e5a}$$ where $P(i|t)$ is the probability that the Markov chain in question taking the value of $i$ at time $t$, and $|i\rangle$ is a basis vector, linearly independent of any other basis vector with different value. Note that $H(t)$ in \eqref{e5} is time-dependent implying that the process is time inhomogeneous. The Lie algebraic method introduced in \cite{r17} requires a decomposition of the operator $H(t)$ as $$H(t)=\sum_{i=1}^ma_i(t)H_i,\label{e6}$$ such that $a_i(t)$ are real-valued functions, and $H_i$ are linearly independent constant operators generating a Lie algebra $V=\text{span}\{H_1,\cdots,H_m\}$ by implementing a Lie bracket $$[H_i,H_j]=H_iH_j-H_jH_i=\sum_{k=1}^m\xi_{ijk}H_k\label{e7}$$ for some real $\xi_{ijk}$. The solution of system \eqref{e5} can be uncoupled into a product of exponentials \cite{r17} $$p(t)=e^{g_1(t)H_1}\cdots e^{g_m(t)H_m}p(0)=U(t)p(0),\label{e8}$$ where $g_i(t)$ are real-valued functions and $g_i(0)=0$. Substituting \eqref{e6} and \eqref{e8} into \eqref{e5}, we obtain \begin{aligned} \frac{{\rm d}p(t)}{{\rm d}t} &= \sum_{i=1}^ma_i(t)H_iU(t)p(0) \\ &=\sum_{i=1}^m\dot{g}_i(t)\Big(\prod_{j=1}^{i-1}e^{g_j(t)H_j}\Big) H_i\Big(\prod_{j=i}^{m}e^{g_j(t)H_j}\Big)p(0). \end{aligned}\label{e9} On multiplying $U(t)^{-1}$ on both sides of \eqref{e9}, we have \begin{aligned} &\sum_{i=1}^ma_i(t)H_i\Big(\prod_{j=1}^me^{g_j(t)\text{ad}H_j}\Big)p(0)\\ &= \sum_{i=1}^ma_i(t)H_iU(t)p(0)U(t)^{-1} \\ &= \sum_{i=1}^m\dot{g}_i(t)\Big(\prod_{j=1}^{i-1}e^{g_j(t)H_j}\Big) H_i\Big(\prod_{j=i}^{m}e^{g_j(t)H_j}\Big)p(0)U(t)^{-1} \\ &= \sum_{i=1}^m\dot{g}_i(t)\Big(\prod_{j=1}^{i-1}e^{g_j(t)\text{ad}H_j}\Big) H_i\Big(\prod_{j=1}^me^{g_j(t)\text{ad}H_j}\Big)p(0). \end{aligned}\label{e10} Since $p(0)$ is arbitrary, we conclude that $$\sum_{i=1}^ma_i(t)H_i=\sum_{i=1}^m\dot{g}_i(t) \Big(\prod_{j=1}^{i-1}e^{g_j(t)\text{ad}H_j}\Big)H_i.\label{e12}$$ From \eqref{e12} we derive a linear relation between $a_i(t)$ and $\dot{g}_i(t)$ with initial values $g_i(0)=0$ (involving $\xi_{ijk}$), as the operators $H_i$ are linearly independent. The calculation of $p(t)$ is achieved in $O(1)$ through \eqref{e12} rather than $O(t)$ by means of incremental direct integrations. Therefore, the computation complexity can be dramatically reduced. The matrix exponential form \eqref{e8} would be useful if the derivative of the solution with respect to some model parameter is required in subsequent calculations \cite{r17,r22}. \section{An example: SIS epidemic spreading} The susceptible-infected-susceptible (SIS) epidemiological model \cite{r23} is an accurate yet simple representation of endemic infections. It is often used as a paradigm for many sexually transmitted infections and computer virus propagations \cite{r2,r9,r26}. The model describes the evolution of an infection in a fixed population, where $N$ individuals in the population are divided into two subclasses: the susceptible pool, of size $S$, and the infected (and infectious) class, of size $I$, with $S+I=N$. Susceptible individuals become infected at a rate $\beta(t)$ by contagion from infected individuals, and infected individuals, in turn, recover (and once again become susceptible) at a rate $\gamma(t)$. The above description of the SIS model leads to a Markovian process \cite{r10} whose probability vector can be written as $$|p(t)\rangle=\sum_{S,I}P(S,I|t)|S,I\rangle,\label{e13}$$ where $P(S,I|t)$ denotes the probability that there are $S$ susceptible individuals and $I$ infected ones at time $t$. $|S,I\rangle$ is a basis vector, linearly independent of other basis vectors with different susceptible and infected numbers. The state space $\mathcal{S}$ consists of $N+1$ elements. The Kolmogorov equation governing this process can be written as $$\frac{{\rm d}}{{\rm d}t}|p(t)\rangle=H(t)|p(t)\rangle,\label{e15}$$ with $$H(t)=\gamma(t)(\hat{\rho}-\hat{I})+\beta(t)(\hat{\sigma}-\hat{S}),\label{e16}$$ where $$\begin{gathered} \hat{S}|S,I\rangle=S|S,I\rangle\\ \hat{I}|S,I\rangle=I|S,I\rangle\\ \hat{\rho}|S,I\rangle=I|S+1,I-1\rangle\\ \hat{\sigma}|S,I\rangle=S|S-1,I+1\rangle \end{gathered}\label{e17}$$ and all these operators are linear operators (note that a similar collection is derived for SIR model in \cite{r18}). Table 1 shows the complete set of Lie brackets, under which the algebra $V=\operatorname{span}\{\hat{S},\hat{I},\hat{\rho},\hat{\sigma}\}$ is closed. \begin{table}[ht] \centering \begin{tabular*}{8cm}{@{\extracolsep{\fill}}|ccccc|} \hline $\hat{X}$&$[\hat{X},\hat{S}]$&$[\hat{X},\hat{I}]$&$[\hat{X},\hat{\rho}]$&$[\hat{X},\hat{\sigma}]$ \\\hline $\hat{S}$&0&0&$\hat{\rho}$&$-\hat{\sigma}$\\$\hat{I}$&0&0&$-\hat{\rho}$&$\hat{\sigma}$\\ $\hat{\rho}$&$-\hat{\rho}$&$\hat{\rho}$&0&$\hat{S}-\hat{I}$\\ $\hat{\sigma}$&$\hat{\sigma}$&$-\hat{\sigma}$&$\hat{I}$&0\\ \hline \end{tabular*} \caption{Values of $[\hat{X},\hat{Y}]$ for SIS model.} \end{table} \begin{table}[ht] \centering \begin{tabular*}{10cm}{@{\extracolsep{\fill}}|ccccc|} \hline $\hat{X}$&$e^{g(\text{ad}\hat{X})}\hat{S}$&$e^{g(\text{ad}\hat{X})}\hat{I}$&$e^{g(\text{ad}\hat{X})}\hat{\rho}$ &$e^{g(\text{ad}\hat{X})}\hat{\sigma}$ \\\hline $\hat{S}$&$\hat{S}$&$\hat{I}$&$e^g\hat{\rho}$&$e^{-g}\hat{\sigma}$\\ $\hat{I}$&$\hat{S}$&$\hat{I}$&$e^{-g}\hat{\rho}$&$e^{g}\hat{\sigma}$\\ $\hat{\rho}$&$\hat{S}-g\hat{\rho}$&$\hat{I}+g\hat{\rho}$&$\hat{\rho}$&$\hat{\sigma}+g\hat{S}-g\hat{I}-g^2\hat{\rho}$\\ $\hat{\sigma}$&$\hat{S}+g\hat{\sigma}$&$\hat{I}-g\hat{\sigma}$&$\hat{\rho}+g\hat{I}-\frac{g^2}{2}\hat{\sigma}$&$\hat{\sigma}$\\ \hline \end{tabular*} \caption{Values of $e^{g(\text{ad}\hat{X})}\hat{Y}$ with a scalar $g$ for SIS model.} \end{table} We need to look for a solution of the form $$|p(t)\rangle=e^{g_1(t)\hat{S}}e^{g_2(t)\hat{I}}e^{g_3(t) \hat{\sigma}}e^{g_4(t)\hat{\rho}}|p(0)\rangle.\label{e18}$$ Employing \eqref{e12} and the action of exponential operators shown in Table 2, we obtain $$\begin{split} &\gamma(t)\hat{\rho}-\gamma(t)\hat{I}+\beta(t)\hat{\sigma}-\beta(t)\hat{S}\\ &= \dot{g}_1(t)\hat{S}+\dot{g}_2(t)\hat{I}+\dot{g}_3(t)e^{g_2}e^{-g_1} \hat{\sigma} +\dot{g}_4(t)\Big(e^{-g_2}e^{g_1}\hat{\rho}+g_3\hat{I} -\frac{g_3^2}{2}e^{g_2}e^{-g_1}\hat{\sigma}\Big). \end{split}\label{e19}$$ Equating terms in \eqref{e19} in front of the same basis matrices yields $$\begin{gathered} g_1(t)=-{\rm B}(t),\\ g_4(t)=\int_0^t\gamma(u)e^{g_2(u)+{\rm B}(u)}{\rm d}u,\\ \end{gathered}\label{e20}$$ where ${\rm B}(t)=\int_0^t\beta(u){\rm d}u$; $g_2(t)$ and $g_3(t)$ are determined by the initial value problem $$\begin{gathered} \dot{g}_2(t)=-\gamma(t)-g_3(t)\gamma(t)e^{g_2(t)+{\rm B}(t)},\\ \dot{g}_3(t)=\beta(t)e^{-{\rm B}(t)-g_2(t)} +\frac{\gamma(t)}{2}e^{g_2(t)+{\rm B}(t)}g_3(t)^2,\\ g_2(0)=g_3(0)=0. \end{gathered} \label{e21}$$ The function $g_3$ satisfies a Riccati equation, which may be solved by standard reduction techniques or numerical integration; see e.g. \cite{r24}. Let $|\mathcal{I}(t)\rangle=\sum_{S,I}I|S,I\rangle$, and then $I(t)=\langle\mathcal{I}(t)|p(t)\rangle$ is the number of infected individuals in the population at time $t$. In Fig. 1 we illustrate $I(t)$ with respect to different choices of $\beta(t)$ and $\gamma(t)$ in a population of size $N=100$. \begin{figure}[htb] \includegraphics[width=0.7\textwidth]{fig1} \caption{Dynamics of the SIS model with $N=100$ and $|p(0)\rangle=|0.9\cdot N,0.1\cdot N\rangle$. $I(t)=\langle\mathcal{I}(t)|p(t)\rangle$ are plotted with respect to $\beta(t)=1,\gamma(t)=0.5$ (circles), $\beta(t)=0.5,\gamma(t)=1$ (squares), and $\beta(t)=e^t,\gamma(t)=1$ (triangles), where $|p(t)\rangle$ are obtained from \eqref{e18}.} \end{figure} When $\gamma=0$, the model reduces to a simple susceptible-infected (SI) epidemics, where individuals, once infected, are infected (and infectious) forever. In this case, the solution of \eqref{e19} can be obtained as $$\begin{gathered} g_1(t)=-{\rm B}(t),\\ g_3(t)=\int_0^te^{-{\rm B}(u)}\beta(u){\rm d}u,\\ g_2(t)=g_4(t)=0, \end{gathered}\label{e23}$$ where ${\rm B}(t)=\int_0^t\beta(u){\rm d}u$. Note that this can be derived similarly from the SIR model addressed in \cite{r18}. The consistency confirms that our result is valid. \subsection*{Conclusion} In this paper, we showed that it is possible to solve susceptible-infected-susceptible (SIS) model via a Lie algebra methodology. Lie algebra solution of differential equations has found host of useful applications in physical systems, where wealthy symmetries exist. Due to insufficient symmetry, this method is not widely used in biological or social systems. For future work, more complex and realistic epidemiological mechanisms, such as susceptible-exposed-infectious-recovered (SEIR) model, are worthy of further research. \begin{thebibliography}{00} \bibitem{r14} D. Alonso, A. McKane; \emph{Extinction dynamics in mainland-island metapopulations: an N-patch stochastic model}, Bull. Math. Biol. 64(2002) 913--958. \bibitem{r19} A. Alvermann, H. Fehske; \emph{High-order commutator-free exponential time-propagation of driven quantum systems}, J. Comput. Phys. 230(2011) 5930--5956. \bibitem{r2} R. M. Anderson, R. M. May; \textit{Infectious Diseases of Humans}, Oxford University Press, Oxford, 1991. \bibitem{r13}A. D. Barbour, M. J. Luczak; \emph{A law of large numbers approximation for Markov population processes with countably many types}, Probab. Theory Relat. Fields 153(2012) 727--757. \bibitem{r4}M. Barth\'el\'emy, A. Barrat, R. Pastor-Satorras, A. Vespignani; \emph{Dynamical patterns of epidemic outbreaks in complex heterogeneous networks}, J. Theor. Biol. 235(2005) 275--288. \bibitem{r20}S. Blanes, F. Casas, J. A. Oteo, J. Ros; \emph{The Magnus expansion and some of its applications}, Phys. Rep. 470(2009) 151--238. \bibitem{r5}S. G\'omez, J. G\'omez-Garde{\~n}es, Y. Moreno, A. Arenas; \emph{Nonperturbative heterogeneous mean-field approach to epidemic spreading in complex networks}, Phys. Rev. E 84(2011) 036105. \bibitem{r18}T. House; \emph{Lie algebra solution of population models based on time-inhomogeneous Markov chains}, J. Appl. Probab. 49(2012) 472--481. \bibitem{r21} J. E. Humphreys; \emph{Introduction to Lie Algebras and Representation Theory}, Springer, New York, 1972. \bibitem{r10} M. J. Keeling, J. V. Ross; \emph{On methods for studing stochastic disease dynamics}, J. R. Soc. Interface 5(2008) 171--181. \bibitem{r1} W. O. Kermack, A. G. McKendrick; \emph{A contribution to the mathematical theory of epidemics}, Proc. R. Soc. Lond. A 115(1927) 700--721. \bibitem{r7} M. E. J. Newman; \emph{Spread of epidemic disease on networks}, Phys. Rev. E 66(2002) 016128. \bibitem{r11}J. R. Norris; \textit{Markov Chains}, Cambridge University Press, New York, 1997. \bibitem{r8} M. Marder; \emph{Dynamics of epidemics on random networks}, Phys. Rev. E 75(2007) 066103. \bibitem{r16} O. Ovaskainen; \emph{The quasistationary distribution of the stochastic logistic model}, J. Appl. Probab. 38(2001) 898--907. \bibitem{r3} R. Pastor-Satorras, A. Vespignani; \emph{Epidemic spreading in scale-free networks}, Phys. Rev. Lett. 86(2001) 3200--3203. \bibitem{r24} A. D. Polyanin, V. F. Zaitsev; \emph{Handbook of Exact Solutions for Ordinary Differential Equations}, 2nd Edition, CRC Press, Boca Raton, 2003. \bibitem{r12}J. V. Ross; \emph{Computationally exact methods for stochastic periodic dynamics spatiotemporal dispersal and temporally forced transmission}, J. Theor. Biol. 262(2010) 14--22. \bibitem{r25}Y. Shang; \emph{Likelihood estimation for stochastic epidemics with heterogeneous mixing populations}, Int. J. Comput. Math. Sci. 6(2012) 34--38. \bibitem{r15}Y. Shang; \emph{Distribution dynamics for SIS model on random networks}, J. Biol. Syst. 20(2012) 213--220. \bibitem{r9}Y. Shang; \emph{Optimal attack strategies in a dynamic botnet defense model}, Appl. Math. Inf. Sci. 6(2012) 29--33. \bibitem{r26}Y. Shang; \emph{Optimal control strategies for virus spreading in inhomogeneous epidemic dynamics}, Canad. Math. Bull. doi:10.4153/CMB-2012-007-2. \bibitem{r27}Y. Shang; \emph{Mixed SI(R) epidemic dynamics in random graphs with general degree distributions}, Appl. Math. Comput. doi:10.1016/j.amc.2012.11.026. \bibitem{r17} J. Wei, E. Norman; \emph{Lie algebra solution of linear differential equations}, J. Math. Phys. 4(1963) 575--581. \bibitem{r23} G. H. Weiss, M. Dishon; \emph{On the asymptotic behavior of the stochastic and deterministic models of an epidemic}, Math. Biosci. 11(1971) 261--265. \bibitem{r22} R. M. Wilcox; \emph{Exponential operators and parameter differentiation in quantum physics}, J. Math. Phys. 8(1967) 962--982. \end{thebibliography} \end{document}