\documentclass[reqno]{amsart} \usepackage{hyperref} \usepackage{graphicx} \AtBeginDocument{{\noindent\small \emph{Electronic Journal of Differential Equations}, Vol. 2015 (2015), No. 50, pp. 1--14.\newline ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu \newline ftp ejde.math.txstate.edu} \thanks{\copyright 2015 Texas State University - San Marcos.} \vspace{9mm}} \begin{document} \title[\hfilneg EJDE-2015/50\hfil Optimal control of an SIR model] {Optimal control of an SIR model with changing behavior through an education campaign} \author[H. R. Joshi, S. Lenhart, S. Hota, F. Agusto \hfil EJDE-2015/50\hfilneg] {Hem Raj Joshi, Suzanne Lenhart, Sanjukta Hota, Folashade Agusto} \address{Hem Raj Joshi \newline Dept. of Mathematics and Computer Science, Xavier University, Cincinnati, OH 45207-4441, USA} \email{joshi@xavier.edu} \address{Suzanne Lenhart \newline Dept. of Mathematics, University of Tennessee Knoxville, TN 37996-1320, USA} \email{lenhart@math.utk.edu} \address{Sanjukta Hota \newline Dept. of Mathematics and Computer Science, Fisk University, Nashville, TN 37208, USA} \email{sanjuktahota@gmail.com} \address{Folashade Agusto \newline Dept. of Mathematics and Statistics, Austin Peay State University, Clarksville, TN 37044, USA} \email{fbagusto@gmail.com} \thanks{Submitted February 4, 2014. Published February 19, 2015.} \subjclass[2000]{35K55, 49K20} \keywords{Effect of education on SIR models; optimal control} \begin{abstract} An SIR type model is expanded to include the use of education or information given to the public as a control to manage a disease outbreak when effective treatments or vaccines are not readily available or too costly to be widely used. The information causes a change in behavior resulting in three susceptible classes. We study stability analysis and use optimal control theory on the system of differential equations to achieve the goal of minimizing the infected population (while minimizing the cost). We illustrate our results with some numerical simulations. \end{abstract} \maketitle \numberwithin{equation}{section} \newtheorem{theorem}{Theorem}[section] \newtheorem{lemma}[theorem]{Lemma} \allowdisplaybreaks \section{Introduction}\label{S:1} The effects of changing behavior is important in epidemic outcomes, and now such effects are beginning to be included in models \cite{ Fenichel, Funk}. Management strategies of how to motivate people to make such behavior changes will become increasingly important. Before the HIV drugs become readily available, the decrease in Uganda HIV rates, in contrast to other countries in the region, was an interesting phenomenon. Some studies indicate that the abstinence, be faithful (AB) and condoms (C) campaigns started by the Ugandan government in 1992 had changed people's behaviors and attitudes, and thus reversed a troubling pattern of increase in HIV/AIDS \cite{Green, Okware}. Other government and nongovernmental agencies began campaigns to distribute information and educational materials about the disease; some organizations emphasized the AB behavior and others the C behavior. Some argue that the Ugandan government initially emphasized only AB strategies, and the effects of C did not appear until later \cite{Singh}. Once the emphasis shifted to the C-type campaigns, there was a dramatic drop in new HIV infections. See \cite{Albright, deWalque, Del Valle, Kelly} about the advances in HIV education and prevention and their effects on the epidemic. Using data about the numbers of organizations giving information and the percentages of the types of behavior recommended and HIV epidemic data in Uganda, Joshi et al \cite{Joshi} developed a SIR model of differential equations, that divided the susceptible class into subclasses based on the AB and C behavior and the resulting different infectivity rates. Parameters in this work \cite{Joshi} were estimated to fit with data about numbers of deaths and infected cases, and the level and type of information given by organizations in Uganda from 1997-2005. When a new virus strain surfaces, vaccines, the usual first line of defense, may not be available. Education or public health campaigns may encourage individuals to change behavior. Recent studies showed the effects of face masks on the transmission of H1N1 influenza during the 2009 pandemic \cite{Aiello}. An economic analysis with SEIR models of differential equations involving face masks for that pandemic showed the reduction of infected cases and of financial losses \cite{Tracht2, Tracht}. Persons could change their behavior to wear surgical face masks or N95 respirators and obtain some level of protection against influenza \cite{Killingley, Loeb}. We are investigating the level of education or information given to the public as a control to manage a disease outbreak when the effective treatments or vaccines are not readily available or too costly to be widely used. We adapt the model from \cite{Joshi} to have three susceptible classes depending on behavior and having different transmission rates and with time-varying education campaign level. With limited resources, the balance between benefits of lower numbers of infecteds and the cost of the education campaign is investigated using optimal control theory on this system of differential equations; the level of education is taken as the control. In \cite{Beh,Cas}, simpler models with differential equations and controls on the level of information represented the influence of the education campaign on the infectivity coefficients and the death rate due to the disease. See \cite{Lun} for stability analysis of a model with two susceptible classes with an education feature. In the next section, we formulate our model and discuss briefly its stability analysis. The optimal control problem, with our goal expressed as an objective functional is given in section 3. Our numerical illustrations and some concluding comments are presented in the following section. \section{Mathematical Model}\label{S:2} We develop an optimal control model of Susceptibles, Infected and Recovered- an SIR type model. In the system of differential equations of the model, the control is the education (or information) level, which helps to change the behavior of some individuals in the susceptible class. Here by ``education", we mean information campaigns or educational outreach materials that give needed information about the disease. This change in behavior leads to subdividing susceptibles into three subclasses, namely $S, S_1$ and $S_2$. A proportion of the susceptible populations, $S$, decide to change their behavior due to an effect of a successful education campaign and thus enter in the $S_1$ or $S_2$ class. These two classes, $S_1$ and $S_2$, have lower transmission rates than the $S$ class and will contribute to lower the number of new infections and thus also lower the recovered/removed population. We consider optimal control of an ordinary differential equation model, which describes the interaction of education with Susceptibles as following: \begin{equation} \label{eq1} \begin{aligned} S'(t) &= -(\alpha_1 + \alpha_2) E(t) S(t) - \beta_1 S(t) I(t) + \Lambda - d S(t) \\ S'_1 (t)&= \alpha_1 E(t) S(t) - \beta_2 S_1(t) I(t) - d S_1(t) \\ S'_2 (t)&= \alpha_2 E(t) S(t) - \beta_3 S_2(t) I(t) - d S_2 (t)\\ I' (t)&= \beta_1 S(t) I(t) + \beta_2 S_1(t) I (t)+ \beta_3 S_2(t) I(t) - d I(t) - \gamma I(t) \\ R' (t)& = \gamma I(t) \end{aligned} \end{equation} with initial conditions $S(0), S_1(0), S_2(0), I(0)$, and $R(0)$. The rate of entering into the $S$ class is $\Lambda$ and the natural death rate is $d$. Individuals only enter the general Susceptible class $S$. Now that we have three susceptible classes, we need three infection rates $\beta_1, \beta_2, \beta_3 $ for $S, S_1$, and $S_2$ respectively for their interactions with the Infected class $I$. Notice that, as a result of interactions of individuals in class $S$ with the control, education $E$, a proportion of the susceptibles leave the general susceptible class $S$ and move to $S_1$ and $S_2$. The rate of moving into class $S_i$ for $i=1,2$ is $\alpha_i E S$. Also, as a result of each susceptible class interacting with the infected class we have individuals leaving at their respective rates and moving to the infected class. The rate $\gamma$ is the transition rate where individuals leave the infected class $I$ and move to the removed class $R$. The removed class $R$ could represent recovered, infected or removed individuals due to disease related deaths (like in Uganda as in \cite{Joshi}). Since model \eqref{eq1} represents human populations, all parameters in the model are non-negative and one can show that the solutions of the system are non-negative, given non-negative initial values. The model \eqref{eq1} will be analyzed in a biologically-feasible region, $\Gamma \subset \mathbb{R}^5_+$ with \[ \Gamma = \big\{(S(t),S_1(t),S_2(t),I(t),R(t))\in \mathbb{R}^{5}_+: 0\leq N(t)\leq \frac{\Lambda }{d} \big\}, \] where $N = S + S_1 +S_2 + I +R$. The following steps are followed to establish the positive invariance of $\Gamma$ (i.e., solutions in $\Gamma$ remain in $\Gamma$ for all $t>0$). The rate of change of the total populations is obtained by adding the equations of the model \eqref{eq1} to give \begin{equation} N'(t)=\Lambda - d N(t). \label{efn} \end{equation} Solving the differential equation \eqref{efn}, we find that $$ N(t)= N(0)e^{-dt}+\frac{\Lambda}{d}(1-e^{-dt}). $$ In particular, $N(t)=\frac{\Lambda}{d}$, if $N(0)=\frac{\Lambda}{d}$. Thus, the region $\Gamma$ is positively-invariant. Hence, it is sufficient to consider the dynamics of the flow generated by \eqref{eq1} in $\Gamma$. In this region, the model is epidemiologically and mathematically well-posed \cite{Het, Lak}. Thus, every solution of the basic model \eqref{eq1} with initial conditions in $\Gamma$ remains in $\Gamma$ for all $t > 0$. Therefore, the $\omega$-limit sets of the system \eqref{eq1} are contained in $\Gamma$. This result is summarized below. \begin{lemma} \label{lem1} The region $\Gamma \subset \mathbb{R}^5_+$ is positively-invariant for the basic model \eqref{eq1} with non-negative initial conditions in $\mathbb{R}^{5}_+$. \end{lemma} To consider the stability of the model, we temporarily assume that the control $E$ is just a constant parameter. Under this assumption, $E(t) =e$, where $e$ is a constant and the model \eqref{eq1} has a disease free equilibrium (DFE), obtained by setting the right-hand sides of the equations in the model to zero, given by \begin{align*} \mathcal{E}_{0} & = (S^*,S_1^*,S_2^*,I^*,R^*) \\ &= \Big(\frac{\Lambda }{(\alpha_1+\alpha_2)e+d},\frac{\Lambda } {(\alpha_1+\alpha_2)e+d}\frac{ \alpha_1e}{d},\frac{\Lambda } {(\alpha_1+\alpha_2)e+d}\frac{ \alpha_2e}{d},0,0\Big). \end{align*} The stability of $\mathcal{E}_{0}$ can be established using the next generation operator method on the system \eqref{eq1}. We take, $I$, as our infected compartment, then using the notation in \cite{van}, the Jacobian matrices $F$ and $V$ for the new infection terms and the remaining transfer terms are respectively given by, $$ F=[\beta_1 S^*+\beta_2 S^*_1+\beta_3 S^*_2] ~~{\rm and}~~V = [d+\gamma]. $$ It follows that the basic reproduction number of the system \eqref{eq1}, denoted by $\mathcal{R}_0$, is given by \begin{equation} \label{Rf0} \mathcal{R}_0 = \rho(FV^{-1}) = \frac{\beta_1 S^*+\beta_2 S^*_1+\beta_3 S^*_2}{d+\gamma}, \end{equation} where $\rho$ is the spectral radius. Further, using \cite[Theorem 2]{van}, the following result is established. \begin{lemma} \label{lem2} The DFE of model \eqref{eq1} (with $E(t)=e$), given by $\mathcal{E}_{0}$, is locally asymptotically stable (LAS) if $\mathcal{R}_0 < 1$, and unstable if $\mathcal{R}_0 > 1$. \end{lemma} The basic reproduction number ($\mathcal{R}_0$) measures the average number of new infections generated by a single infected individual in a completely susceptible population \cite{And,Dic,Het,van}. Thus, Lemma \ref{lem2} implies that the infection can be eliminated from the population (when $\mathcal{R}_0 < 1)$ if the initial sizes of the sub-populations are in the basin of attraction of the DFE, $\mathcal{E}_{0}$. We do not consider an endemic equilibrium since we are considering the case when a disease outbreak has just started. \section{Formulation and analysis of the optimal control problem}\label{S:3} Now we turn our focus to using a time-varying control function $E(t)$, which represents the level of the educational campaign that causes susceptible individuals to change their behavior. The control set $\mathbf{E}$ is \[ \mathbf{E}=\{E(t) : 0\le a \le E(t)\le b < 1, \, 0\le t\le T, \, E(t)\ \text{ is Lebesgue measurable}\}. \] Our goal is to find the control $E(t)$ and associated state variables $S(t)$, $S_1(t)$, $S_2(t)$, $I(t)$, and $R(t)$ to minimize the following objective functional: \[ J[E] = \int_0^T (I(t) - A(S+ S_1 + S_2) + B E) \, dt. \] By choosing appropriate positive balancing constants $A$ and $B$, our goal is to minimize the infected population, and maximize the susceptible population while minimizing the cost of the control. If one only wants to minimize the infected population and not be concerned with the level of the $S, S_1$ and $S_2$ populations, one would take $A=0,$ The structure of this model gives bounded solutions for finite final $T$. This objective functional and the differential equations are linear in the control with bounded states, and one can show by standard results that an optimal control and corresponding optimal states exist \cite{Fleming}. By using Pontryagin's Maximum Principle \cite{Fleming, Lenhart, Pontryagin} we derive necessary conditions for our optimal control and corresponding states. The Hamiltonian is \begin{equation}\label{Hamiltonian} \begin{split} H &= I(t) - A(S(t) + S_1(t) + S_2(t) ) + B E(t) \\ &\quad + \lambda_1 (-\alpha_1 E(t) S(t) - \alpha_2 E(t) S(t) - \beta_1 S(t) I(t) + \Lambda - d S(t)) \\ &\quad + \lambda_2 (\alpha_1 E(t) S(t) - \beta_2 S_1(t) I(t) - d S_1(t)) \\ &\quad + \lambda_3 (\alpha_2 E(t) S(t) - \beta_3 S_2(t) I(t) - d S_2(t)) \\ &\quad + \lambda_4 (\beta_1 S(t) I(t) + \beta_2 S_1(t) I(t) + \beta_3 S_2(t) I(t) - d I(t) - \gamma I(t)) \\ &\quad + \lambda_5 (\gamma I(t)) \end{split} \end{equation} Given an optimal control $E^\ast $, there exist adjoint functions, $\lambda_1, \lambda_2, \lambda_3, \lambda_4, \lambda_5$, corresponding to the states $S, S_1, S_2, I$, and $R$ such that: \begin{equation}\label{eq1.2} \begin{aligned} \lambda_1'& = -\frac{\partial H}{\partial S}\\ &= -[ -A + \lambda_1 (-\alpha_1 E - \alpha_2 E - \beta_1 I -d) + \alpha_1\lambda_2 E + \alpha_2\lambda_3 E + \beta_1 \lambda_4 I ], \\ \lambda_2' &= -\frac{\partial H}{\partial S_1} = -[ -A + \lambda_2 (-\beta_2 I - d) + \beta_2 \lambda_4 I ], \\ \lambda_3'& = -\frac{\partial H}{\partial S_2} =-[ -A + \lambda_3 (-\beta_3 I - d) + \beta_3 \lambda_4 I ], \\ \lambda_4' &= -\frac{\partial H}{\partial I}\\ &= -[ 1 + \lambda_1(-\beta_1 S ) + \lambda_2 (- \beta_2 S_1) -\beta_3 \lambda_3 S_2 + \lambda_4( \beta_1 S + \beta_2 S_1\\ &\quad + \beta_3 S_2 - d - \gamma) + \gamma \lambda_5 ], \\ \lambda_5' &= -\frac{\partial H}{\partial R} =0. \end{aligned} \end{equation} where $\lambda_1(T)=0, \lambda_2(T)=0, \lambda_3(T)=0, \lambda_4(T)=0$, and $\lambda_5(T)=0$ are the transversality conditions. The Hamiltonian is minimized with respect to the control variable at $E^*$. Since the Hamiltonian is linear in the control, we must consider if the optimal control is bang-bang (at its lower or upper bound), singular or a combination. The singular case could occur if the slope or the switching function, \begin{equation} \label{Hamil} \frac{\partial H}{\partial E} = B +[- (\alpha_1 + \alpha_2)\lambda_1 + \alpha_1 \lambda_2 + \alpha_2\lambda_3]S , \end{equation} is zero on non-trivial interval of time. Note that the optimal control would be at its upper bound or its lower bound according to: \begin{equation*} \frac{\partial H}{\partial E} < 0 \quad \text{or}\quad > 0. \end{equation*} To investigate the singular case, let us suppose ${\frac{\partial H}{\partial E} =0}$ on some non-trivial interval. In this case, we calculate $$ {\frac{d}{dt} \big(\frac{\partial H}{\partial E}\big)=0} $$ and then we will show that control is not present in that equation. To solve for the value of the singular control, we will further calculate $$ { \frac{d^2} {dt^2} \big(\frac{\partial H}{\partial E}\big)=0 }. $$ We simplify the time derivative of $ \frac{\partial H}{\partial E}$, \begin{equation} \begin{aligned} 0& = \frac{d}{dt} \big(\frac{\partial H}{\partial E}\big) = \frac{d}{dt}\{B +[- (\alpha_1 + \alpha_2)\lambda_1 + \alpha_1 \lambda_2 + \alpha_2\lambda_3]S \} \\ & = [-(\alpha_1 + \alpha_2)\lambda_1 + \alpha_1 \lambda_2 + \alpha_2\lambda_3]S' + [-(\alpha_1 + \alpha_2)\lambda_1' + \alpha_1 \lambda_2' + \alpha_2\lambda_3']S \end{aligned} \label{eq4} \end{equation} We calculate both sums separately and add them together. The first sum can be written as: \begin{align*} &[-(\alpha_1 + \alpha_2)\lambda_1 + \alpha_1 \lambda_2 + \alpha_2\lambda_3]S' \\ &= [-(\alpha_1 + \alpha_2)\lambda_1 + \alpha_1 \lambda_2 + \alpha_2\lambda_3] [-(\alpha_1 + \alpha_2) E S - \beta_1 S I + \Lambda-dS] \\ &= (\alpha_1 + \alpha_2)^2 \lambda_1 ES + \beta_1 (\alpha_1 + \alpha_2) \lambda_1S I - (\Lambda-d S)(\alpha_1 + \alpha_2) \lambda_1 \\ &\quad - \alpha_1(\alpha_1 + \alpha_2) \lambda_2 E S - \alpha_1 \beta_1 \lambda_2 S I + (\Lambda-d S) \alpha_1 \lambda_2\\ &\quad - \alpha_2(\alpha_1 + \alpha_2) \lambda_3 E S - \alpha_2 \beta_1 \lambda_3 S I + (\Lambda-d S) \alpha_2 \lambda_3 \end{align*} The second sum can be written as: \begin{align*} & (\alpha_1 + \alpha_2)\{-A + \lambda_1 [-(\alpha_1 + \alpha_2) E - \beta_1 I -d] + \alpha_1 \lambda_2 E + \alpha_2 \lambda_3 E + \beta_1 \lambda_4 I \}S \\ & - \alpha_1 [ -A + \lambda_2 (-\beta_2 I - d) + \beta_2 \lambda_4 I ]S - \alpha_2 [ -A + \lambda_3 (-\beta_3 I - d) + \beta_3 \lambda_4 I ] S \\ &= -(\alpha_1 + \alpha_2)^2 \lambda_1 ES - \beta_1 (\alpha_1 + \alpha_2) \lambda_1 IS -d(\alpha_1 + \alpha_2) \lambda_1 S + \alpha_1 (\alpha_1 + \alpha_2) \lambda_2 ES \\ &\quad + \alpha_2(\alpha_1 + \alpha_2) \lambda_3 ES + \beta_1 (\alpha_1 + \alpha_2)\lambda_4 S I +\alpha_1 (\beta_2 I + d) \lambda_2 S -\alpha_1 \beta_2 \lambda_4 S I \\ &\quad + \alpha_2 (\beta_3 I + d)\lambda_3 S - \alpha_2 \beta_3 \lambda_4 S I \end{align*} Thus combining, we have \begin{align*} 0& = \frac{d}{dt} \big(\frac{\partial H}{\partial E}\big) = - \Lambda(\alpha_1 + \alpha_2) \lambda_1 - \alpha_1 \beta_1 \lambda_2 S I + \Lambda \alpha_1 \lambda_2 \\ &\quad - \alpha_2 \beta_1 \lambda_3 S I + \Lambda \alpha_2 \lambda_3 + \beta_1 (\alpha_1 + \alpha_2)\lambda_4 SI \\ &\quad + \alpha_1 \beta_2 \lambda_2 S I -\alpha_1 \beta_2 \lambda_4 S I + \alpha_2 \beta_3 \lambda_3 S I - \alpha_2 \beta_3 \lambda_4 S I \\ &= [ - \Lambda(\alpha_1 + \alpha_2) \lambda_1 + \Lambda \alpha_1 \lambda_2 + \Lambda \alpha_2\lambda_3 ] + (\alpha_1 \beta_2 - \alpha_1 \beta_1) \lambda_2 S I \\ &\quad + (\alpha_2 \beta_3 - \alpha_2 \beta_1) \lambda_3 S I + [\beta_1 (\alpha_1 + \alpha_2) -\alpha_1 \beta_2 - \alpha_2 \beta_3 ] \lambda_4 S I \\ &= \Lambda[ \alpha_1 (\lambda_2 -\lambda_1) + \alpha_2 ( \lambda_3- \lambda_1) ] \\ &\quad + \{ \alpha_1( \beta_2 - \beta_1) \lambda_2 + \alpha_2( \beta_3 - \beta_1) \lambda_3 + [\alpha_1(\beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3)] \lambda_4 \} S I. \end{align*} We see that the control does not explicitly show in this expression, so next we calculate the second derivative with respect to time. \begin{equation} \begin{aligned} 0 &= \frac{d^2}{dt^2} \big(\frac{\partial H}{\partial E}\big) \\ &= \Lambda [ \alpha_1 (\lambda_2' -\lambda_1') + \alpha_2 ( \lambda_3'-\lambda_1') ] + \Big\{ \alpha_1( \beta_2-\beta_1) \lambda_2'+\alpha_2(\beta_3-\beta_1)\lambda_3'\\ &\quad +[\alpha_1(\beta_1-\beta_2) + \alpha_2( \beta_1 - \beta_3)] \lambda_4'\Big\} S I + \Big\{ \alpha_1( \beta_2 - \beta_1) \lambda_2 \\ &\quad + \alpha_2( \beta_3 - \beta_1) \lambda_3 + [\alpha_1(\beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3) ] \lambda_4 \Big\} (S I'+ S' I ) \end{aligned} \label{eq8} \end{equation} Using systems \eqref{eq1} and \eqref{eq1.2}, we simplify \eqref{eq8} as follows \begin{align*} 0 &= \frac{d^2}{dt^2} \big(\frac{\partial H}{\partial E}\big) \\ &= - \Lambda \Big\{ [\beta_1(\alpha_1 + \alpha_2) \lambda_1 - \alpha_1 \beta_2 \lambda_2 - \alpha_2 \beta_3 \lambda_3 + (\alpha_1(\beta_2- \beta_1) + \alpha_2( \beta_3 - \beta_1)) \lambda_4 ]I \\ & \quad + d[( \alpha_1 + \alpha_2) \lambda_1 - \alpha_1 \lambda_2 - \alpha_2 \lambda_3] + \big[(\alpha_1+ \alpha_2)^2 \lambda_1 \\ &\quad -(\alpha_1+ \alpha_2) (\alpha_1 \lambda_2 +\alpha_2 \lambda_3)\big] E \Big\} + \Big\{ \alpha_1( \beta_2 - \beta_1) (\beta_2 I + d) \lambda_2 \\ &\quad + \alpha_2( \beta_3 - \beta_1) (\beta_3 I + d) \lambda_3 + (\alpha_1 \beta_2( \beta_1 - \beta_2) + \alpha_2 \beta_3( \beta_1 - \beta_3)) \lambda_4 I \\ &\quad - (\alpha_1( \beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3))(( \beta_1 S + \beta_2 S_1 + \beta_3 S_2 - d - \gamma )\lambda_4 \\ & \quad +1 + A - \beta_1 \lambda_1 S - \beta_2 \lambda_2 S_1 - \beta_3 \lambda_3 S_2 +\gamma \lambda_5) \Big\} SI \\ &\quad + [ \alpha_1( \beta_2 - \beta_1) \lambda_2 + \alpha_2( \beta_3 - \beta_1) \lambda_3 + (\alpha_1(\beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3) ) \lambda_4 ] \\ &\quad\times \Big\{ (\beta_1 S + \beta_2 S_1 + \beta_3 S_2 - (d + \gamma )) S I + (-(\alpha_1 + \alpha_2) E S - \beta_1 S I + \Lambda - d S) I \Big\}. \end{align*} The above equation can be written in the form $$ \frac{d^2}{dt^2} \big(\frac{\partial H}{\partial E}\big)=\Phi_1(t) E(t) + \Phi_2(t) =0 $$ and we can solve for the singular control as $$ {E_{\rm singular} (t) = - \frac{\Phi_2(t)}{\Phi_1(t)}}, $$ if $$ \Phi_1(t) \neq 0 \quad\text{and}\quad a \leq - \frac{\Phi_2(t)}{\Phi_1(t)} \leq b $$ with \begin{align*} &\Phi_1(t) \\ & = -\Lambda [(\alpha_1+ \alpha_2)^2 \lambda_1 - (\alpha_1+ \alpha_2) ( \alpha_1 \lambda_2 + \alpha_2 \lambda_3)] - \big[ \alpha_1( \beta_2 - \beta_1) \lambda_2\\ &\quad + \alpha_2( \beta_3 - \beta_1) \lambda_3 + (\alpha_1(\beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3) ) \lambda_4 \big] (\alpha_1 + \alpha_2) SI \\ & = -\Lambda [(\alpha_1+ \alpha_2)^2 \lambda_1 - (\alpha_1+ \alpha_2) ( \alpha_1 \lambda_2 + \alpha_2 \lambda_3)] \\ & \quad - [ \alpha_1( \beta_1 - \beta_2) (\lambda_4- \lambda_2) + \alpha_2( \beta_1 - \beta_3) (\lambda_4 -\lambda_3) ] (\alpha_1 + \alpha_2) SI \\ & = - \Lambda (\alpha_1 + \alpha_2) \frac{B}{S} - [ \alpha_1( \beta_1 - \beta_2) (\lambda_4- \lambda_2) + \alpha_2( \beta_1 - \beta_3) (\lambda_4 -\lambda_3) ] (\alpha_1 + \alpha_2) SI \end{align*} and \begin{align*} &\Phi_2 (t)\\ &= -\Lambda \Big\{ [\beta_1(\alpha_1 + \alpha_2) \lambda_1 - \alpha_1 \beta_2 \lambda_2 - \alpha_2 \beta_3 \lambda_3 + (\alpha_1(\beta_2- \beta_1) + \alpha_2( \beta_3 - \beta_1)) \lambda_4 ]I \\ &\quad + d \frac{B}{S} \Big\} + \Big\{ \alpha_1( \beta_2 - \beta_1) (\beta_2 I + d) \lambda_2 + \alpha_2( \beta_3 - \beta_1) (\beta_3 I + d) \lambda_3 + (\alpha_1 \beta_2( \beta_1 - \beta_2) \\ &\quad + \alpha_2 \beta_3( \beta_1 - \beta_3)) \lambda_4 I - (\alpha_1( \beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3))(( \beta_1 S + \beta_2 S_1+ \beta_3 S_2 \\ &\quad - d - \gamma )\lambda_4 +1 + A - \beta_1 \lambda_1 S - \beta_2 \lambda_2 S_1 - \beta_3 \lambda_3 S_2 +\gamma \lambda_5) \Big\} SI + \big[ \alpha_1( \beta_2 - \beta_1) \lambda_2\\ &\quad + \alpha_2( \beta_3 - \beta_1) \lambda_3 + (\alpha_1(\beta_1 - \beta_2) + \alpha_2( \beta_1 - \beta_3) ) \lambda_4 ] \Big\{ (\beta_1 S + \beta_2 S_1 + \beta_3 S_2 \\ & \quad - (d + \gamma )) S I + ( - \beta_1 S I + \Lambda - d S) I \Big\} \end{align*} To check the generalized Legendre-Clebsch condition for the singular control to be optimal, we require ${ \frac{d}{dE} \frac{d^2}{dt^2} \left(\frac{\partial H}{\partial E}\right)=\Phi_1(t) }$ to be negative \cite{Krener}. To summarize, our control characterization is: On a nontrivial interval, \begin{gather*} \text{if $\frac{\partial H}{\partial E} < 0$ at $t$, then $E^*(t) =b$}, \\ \text{if $\frac{\partial H}{\partial E} > 0$ at $t$, then $E^*(t) =a$}, \\ \text{if $\frac{\partial H}{\partial E} = 0$, then $E_{\rm singular} (t) = - \frac{\Phi_2}{\Phi_1}$}. \end{gather*} Hence, our control is optimal at $t$ provided $\Phi_1(t) < 0$ and $ a \leq - \frac{\Phi_2(t) }{\Phi_1(t) } \leq b$. \section{Numerical results and conclusions}\label{S:4} \begin{table}[htb] \caption{Description of the variables and parameters for model \eqref{eq1}} \label{t1} \begin{center} \footnotesize \begin{tabular}{lllll}\hline \hline Variable&Description && &\\ \hline $S(t)$ &Susceptible humans& & & \\ $S_1(t),~S_2(t)$ &Susceptible humans who change their & &\\ & behavior due to education campaign&&& \\ $I(t)$ &Infected humans && \\ $R(t)$ &Removed humans& & & \\ \hline \hline Parameter& Description & Baseline value & \\ \hline $\Lambda$ &Recruitment rate & $0.005$ & \\ $d$ & Death rate & $0.0015$ & \\ $\alpha_1~\alpha_2$ &Transfer rate to the educated susceptible classes & $0.0019, ~0.0152$ & \\ $\beta_1,~\beta_2,~\beta_3$ &Infection rate & $0.0040, ~0.0002, ~0.0016$ & \\ $\gamma$ & Removal rate & $0.005$ & \\ $a,~ b$ & Control lower and upper bound & $0,~~0.85$ & \\ $A,~ B$ & Balancing constant & $0,~~5\times 10^{-2}$ & \\ \hline \hline \end{tabular} \end{center} \end{table} \begin{figure}[htb] \begin{center} \includegraphics[width=0.98\textwidth]{fig1.pdf} % lenhartfeb.pdf} \end{center} \caption{Simulation results for \eqref{eq1}, using the parameter values in Table \ref{t1}. Dashed lines are for the \lq\lq without control\rq\rq case, and solid lines for the \lq\lq with control\rq\rq case.} \label{fig1} % \label{f0} \label{f01}} \end{figure} \begin{figure}[htb] \begin{center} \includegraphics[width=0.98\textwidth]{fig2.pdf} % JoshichangeAlpha2.pdf \end{center} \caption{Simulation results for \eqref{eq1} varying $\alpha_2$, using the parameter values in Table \ref{t1}, solid lines for $\alpha_2 = 0.0152$, and dashed lines for $\alpha_2 = 0.152$.} \label{fig2} % \label{f1} }\label{f11}} \end{figure} \begin{figure}[htb] \begin{center} \includegraphics[width=0.98\textwidth, trim = 40mm 85mm 40mm 85mm, clip ]{fig3.pdf} % JoshiCompareB.pdf \end{center} \caption{Simulation results for \eqref{eq1} varying $\Lambda$, using the parameter values in Table \ref{t1}, solid lines for $\Lambda = 0.005$, $\mathcal{R}_0 = 2.0513$, and dashed lines for $\Lambda = 0.05$ $\mathcal{R}_0 = 20.5128$.} \label{fig3} % \label{f2} \label{f21} \end{figure} \begin{figure}[htb] \begin{center} \includegraphics[width=0.98\textwidth, trim = 40mm 85mm 40mm 85mm, clip ]{fig4.pdf} % JoshiCompareA.pdf \end{center} \caption{Simulation results for \eqref{eq1} varying the balancing constant $A$, using parameter the values in Table \ref{t1}, solid lines for $A = 0$, and dashed lines for $A = 10$.} \label{fig4} % \label{f3} \label{f31}} \end{figure} \begin{figure}[htb] \begin{center} \includegraphics[width=0.98\textwidth, trim = 40mm 85mm 40mm 85mm, clip ]{fig5.pdf} % JoshiCompareNEWUpperBound.pdf \end{center} \caption{Simulation results for \eqref{eq1} varying upper bound, using the parameter values in Table \ref{t1}, solid lines for upper bound, and dashed lines for upper bound $= 2$.} \label{fig5} %\label{f4} }\label{f41}} \end{figure} The optimality system is the state and adjoint systems coupled with the optimal control characterization. We solved optimality system numerically using the forward-backward sweep method \cite{ Hackbush, Lenhart}. Starting with an initial guess for the control, the state system is solved forward in time. Using those new state values, the adjoint system is solved backward in time. The control is updated using a convex combination of the old control values and the new control values from the characterization. The iterative method is repeated until convergence. We note that the uniqueness of the optimal control can be proven for the final time $T$ sufficiently small. But in our numerical simulations, we did not find an indication of non-uniqueness and did not encounter any occurrence of the singular case. We explore the transmission model \eqref{eq1} to study the effects of time dependent control measures using parameter values in Table \ref{t1} and initial conditions, $S(0) = 5$, $S_1(0)=0$, $S_2(0)=0$, $I(0) = 1.2$, $R(0) = 0.08$, $E(0) = 0.5$, except when otherwise stated. With no control, the basic reproductive number $\mathcal{R}_0 $ is $2.0513$, thus, indicating the disease free equilibrium is unstable. Here $S(0)$, $S_1(0)$, $S_2(0)$, $I(0) $ and $R(0)$, as well as the corresponding states in the figures, are in millions of individuals. Figure \ref{fig1} shows a higher number of susceptible individuals in the absence of educational campaign (without control) compared to the presence of educational campaigns (with control). This is due to the fact that susceptible individuals in the community are not changing their behavior which causes them to move to either of the two other susceptible classes $S_1$ and $S_2$. In Figure \ref{fig2}, we varied the transfer rate into the susceptible class $S_2$ by increasing $\alpha_2$ ten fold (i.e. $\alpha_2=0.152$) and then compare the two educational campaign cases. We observed that the increase naturally lead to an increase in the $S_2$ class which results in a subsequent reduction in the total number of infected individuals in the community and the length of the time with positive control from about $31$ days to $26$ days. Less control effort is needed due to an increasing rate of behavior rate (change of transition rate) for $S_2$. In varying the recruitment or birth rate $\Lambda$ from $\Lambda=0.005$ to $\Lambda=0.05$ (a ten fold increase), we observed by comparing the two educational campaign cases in Figure \ref{fig3}, an increase in the various classes which resulted in an increase in the control time from about $31$ days to about $35$ days. We equally observed a ten fold increase in $\mathcal{R}_0 $ from $2.0513$ to $20.513.$ Next we consider the balancing constant $A$. We increase the constant $A$ from $A=0$ to $A=10$, this indicates that it is important to maximize the various susceptible populations. We observed from Figure \ref{fig4}, an increase in the various classes which resulted in an increase in the control time from about $31$ days to about $45$ days. Lastly, we varied the control upper bound from $0.85$ to $2$, and we observed from Figure \ref{fig5}, a decrease in the susceptible class $S$ and an increase in classes $S_1$ and $S_2$ leading to a reduction in the total number of infected. With this increase in the upper bound, there is a decrease in the control time from about $31.5$ days to about $26.5$ days. In conclusion, we illustrated optimal controls for several scenarios with a model with three susceptible classes due to changing behavior. The behavior changes result from information distributed to susceptibles. This work demonstrates an optimal control tool in making decisions about allocating efforts to slow down an epidemic with an information educational campaign. In future work, we will investigate models in which education campaigns and treatment are both important options for the disease management. \subsection*{Acknowledgments} HRJ was partially supported by Xavier Universities Jesuit Fellowship. The authors acknowledge partial support for short-term visits at National Institute for Mathematical and Biological Synthesis (NIMBioS). NIMBioS is an institute sponsored by the National Science Foundation, the U.S. Department of Homeland Security, and the U.S. Department of Agriculture through NSF Award \#EF-0832858, with additional support from The University of Tennessee, Knoxville. Lenhart is also partially supported by the Center for Business and Economic Research at the University of Tennessee. \begin{thebibliography}{99} \bibitem{Aiello} Aiello, A.; Murray, G.; Coulborn, R.; Davis, B. M.; Uddin, M.; Shay, D. K.; Watermal, S. H.; Moto, A. S.\; Mask Use, Hand Hygiene, and Seasonal Influenza-Like Illness Among Young DULTS: A Randomized Intervention Trial. \emph{Journal of Infectious Diseases} \textbf{201} (2010), 491-498. \bibitem{Albright} Albright, K.; Kawooya, D.; The Role of Information in Uganda's Reduction of HIV/AIDS Prevalence: Individual Perceptions of HIV/AIDS Information. \emph{Information Development}, \textbf{21}(2) (2005), 106-112. \bibitem{And} Anderson, R. M.; May, R.; Infectious Diseases of Humans, Oxford University Press, New York (1991). \bibitem{Beh} Behncke, H., (2000) Optimal control of deterministic epidemics. \emph{Optimal Control Applications and Methods} \textbf{21}., 269-275. \bibitem{Cas} Castilho, C.; Optimal control of an epidemic through education campaigns. \emph{Electronic Journal of Differential Equations } No. ?? (2006) 1-11. \bibitem{deWalque} de Walque, D.; How does the Impact of an HIV/AIDS Information Campaign vary with Educational Attainment? Edidence from Uganda. \emph{Journal of Development Economics}, \textbf{84} (2007), 686-714. \bibitem{Del Valle} Del Valle, S.; Evangelista, A. M.; Velasco, M. C.; Kribs-Zaleta, C. M.; Hsu Schmitz, S-F.; Efforts of Education, Vaccination, and Treatment on HIV Transmission in Homosexuals with Genetic Heterogeneity, \emph{Math. Biosciences} \textbf{187} (2004), 111-122. \bibitem{Dic} Diekmann, O.; Heesterbeek, J. A. P.; Metz, J. A. P.; On the definition and computation of the basic reproduction ratio $R_0$ in models for infectious diseases in heterogeneous populations. \emph{J. Math. Biol}, \textbf{28} (1990), 503-522. \bibitem{Fleming} Fleming, W. H.; Rishel, R. W.; Deterministic and Stochastic Optimal Control. Springer Verlag, New York, (1975). \bibitem{Fenichel} Fenichel, E.; Castillo-Chavez, C.; Ceddia, M. G.; Chowell, G.; Gonzalez Parra, P. A; Hickling, G. J.; Holloway, G.; Horan, R.; Morin, B.; Perrings, C.; Springborn, M.; Velazquez, L.; Villalobos, C,; Adaptive human behavior in epidemiological models, \emph{Proceedings Of The National Academy Of Sciences Of The United States Of America}, \textbf{108}(15) (2011), 6306-6311. \bibitem{Funk} Funk, S.; Salath, M.; Jansen, V.; Modelling the influence of human behaviour on the spread of infectious diseases: A review. \emph{J R Soc Interface}, \textbf{(7)} (2010), 1247-1256. \bibitem{Green} Green, E. C.; Halperin, D. T.; Nantulya; Hogle, J. A.; Uganda's HIV PrevensionSuccess: The Role of Sexual Behavior Change and The National Response. \emph{AIDS Behavior}, \textbf{10}(4) (2006), 335-346. \bibitem{Hackbush} Hackbush, W.; A Numerical Method for Solving Parabolic Equations with Opposite Orientation, {\em Computing}, \textbf{20}(3) (1978), 229-240. \bibitem{Het} Hethcote, H. W.; The mathematics of infectious diseases, \emph{SIAM Rev}, \textbf{42}(4) (2000), 599-653. \bibitem{joshi} Joshi, H.; Lenhart, S.; Albright, K.; Gipson, K.; Modeling the Effect of Information Campaigns On the HIV Epidemic In Uganda, \emph{Mathematical Biosciences and Engineering}, \textbf{5}(4), (2008), 757-770. \bibitem{Joshi} Margevicius, R.; Joshi, H. R.; The Influence of Education in Reducing the HIV Epidemic, \emph{Involve},\textbf{6}(2) (2013), 127-135. \bibitem{Kelly} Kelly, J. A.; Advances in HIV/AIDS Education and Prevention, \emph{Family Relations}, \textbf{44} (1995), 345-352. \bibitem{Killingley} Killingley, B.; Respirators Versus Medical Masks: Evidence Accumulates but the Jury Remains Out, \emph{Influenza Other Respiratory Virus}, \textbf{5} (2011), 143-45. \bibitem{Krener} Krener, A. J.; The High Order Maximum Principle and its Application to Singular Exteremals, \emph{SIAM Journal on Control and Optimization}, \textbf{15} (1997), 256-293. \bibitem{Lak} Lakshmikantham, V.; Leela, S.; Martynyuk, A. A.; Stability Analysis of Nonlinear Systems. Marcel Dekker, Inc., New York and Basel, 1989. \bibitem{Lenhart} Lenhart, S.; Workman, J. T.; \emph{Optimal Control Applied to Biological Models}, Chapman and Hall, 2007. \bibitem{Lun} Lungu, E. M.; Kgosimore, M.; Nyababza; Models for the spread of HIV/AIDS: Trends in Southern Africa, \emph{Mathematical Studies on Human Dynamics: Emerging Paradigms and Challenges}, editors, A. B. Gumel, C. Castillo-Chavez, R. E. Mickens, and D. P. Clemence, American Mathematical Society , 2007. \bibitem{Okware} Okware, S.; Kinsman, J.; Onyango, S.; Opio, A.; Kaggwa, P.; Revisiting the ABC strategy: HIV prevention in Uganda in the era of antiretroviral therapy, \emph{Postgraduate Medical Journal}, \textbf{81}(960) (2005), 625-628. \bibitem{Loeb} Loeb, M.; Dafoe, N.; Mahony, J.; John, M.; Sarabia, A.; Glavin, V.; Webby, R.; Smieja, M.; Earn, D. J. D.; Chong, S.; Webb, A.; Walter, S. D.; Surgical Mask vs N95 Respirartor for Preventing Influenza Amnong Health Care Workers, \emph{J. Am. Med. Assoc.}, \textbf{302} (2009), 1865-1871. \bibitem{Pontryagin} Pontryagin, L. S.; Boltyanskii, V. G.; Gamkrelidze, R. V.; Mishchenko, E. F.; \emph{The mathematical theory of optimal processes}, Wiley, New York, 1962. \bibitem{Singh} Singh, S.; Darroch, J. E.; Bankole, A.; The roles of abstinence, monogamy and condom use in HIV decline, \emph{Reproductive Health Matters}, \textbf{12} (2004), 129--135. \bibitem{van} van den Driessche, P.; Watmough, J.; Reproduction Numbers and sub-threshold Endemic Equilibria for Compartmental Models of Disease Transmission, \emph{Math. Biosci.} \textbf{180} (2002), 29-48. \bibitem{Tracht2} Tracht, S. M.; Del Valle, S. Y.; Edwards, B. K.; Economic Analysisof the Use of Facemasks During Pandemic (H1N1) 2009, \emph{Journal of Theoretical Biology} \textbf{300} (2012), 161-172. \bibitem{Tracht} Tracht, S. M.; Del Valle, S. M.; Hyman, J. M.; Mathematical Modelling of the Effectiveness of Facemasks in Reducing the Spread of Novel Inuenza A (H1N1), \emph{PloS One} \textbf{5}(2) (2010), 1-12. \end{thebibliography} \end{document}