Modelling of Rabies Transmission Dynamics Using Optimal Control Analysis

We examine an optimal way of eradicating rabies transmission from dogs into the human population, using preexposure prophylaxis (vaccination) and postexposure prophylaxis (treatment) due to public education. We obtain the disease-free equilibrium, the endemic equilibrium, the stability, and the sensitivity analysis of the optimal control model. Using the Latin hypercube sampling (LHS), the forward-backward sweep scheme and the fourth-order Range-Kutta numerical method predict that the global alliance for rabies control’s aim of working to eliminate deaths from canine rabies by 2030 is attainable throughmass vaccination of susceptible dogs and continuous use of preand postexposure prophylaxis in humans.


Introduction
Rabies is an infection that mostly affects the brain of an infected animal or individual, caused by viruses belonging to the genus Lyssavirus of the family Rhabdoviridae and order Mononegavirales [1,2].This disease has become a global threat and it is also estimated that rabies occurs in more than 150 countries and territories [2].Raccoons, skunks, bats, and foxes are the main animals that transmit the virus in the United States [2].In Asia, Africa, and Latin America, it is known that dogs are the main source of transmission of the rabies virus into the human population [2].When the rabies virus enters the human body or that of an animal, the infection (virus) moves rapidly along the neural pathways to the central nervous system; from there the virus continues to spread to other organs and causes injury by interrupting various nerves [2].The symptoms of rabies are quite similar to those of encephalitis (see [3]).Due to movement of dogs in homes or the surroundings, the risk of not being infected by a rabid dog can never be guaranteed.Rabies is a major health problem in many populations dense with dogs, especially in areas where there are less or no preventive measures (vaccination and treatment) for dogs and humans.Treatment after exposure to the rabies virus is known as postexposure prophylaxis (PEP) and vaccination before exposure to the infection is known as preexposure prophylaxis.
The study of optimal control analysis in maximizing or minimizing a said target was introduced by Pontryagin and his collaborators around 1950.They developed the key idea of introducing the adjoint function to a differential equation, by forming an objective functional [4], and since then there has been a considerable study of infectious disease using optimal control analysis (see [4][5][6][7][8][9][10][11][12]).
Research published by Aubert [13], on the advancement of the expense of wildlife rabies in France, incorporated various variables.They follow immunization of domestic animals, the reinforcement of epidemiological reconnaissance system and the bolster given to indicative research laboratories, the costs connected with outbreaks of rabies, the clinical perception of those mammals which had bitten humans, the preventive immunization, and postexposure treatment of people.A significant percentage (72%) of the cost was the preventive immunization of local animals.In France, as in other European nations in which the red fox (Vulpes) is the species most affected, two primary procedures for controlling rabies were assessed in [13] at the repository level to be specific: fox termination and the oral immunization of foxes.The consolidated costs and advantages of both systems were looked at and included either the expenses of fox separation or the cost of oral immunization.The total yearly costs of both techniques stayed practically identical until the fourth year, after which the oral immunization methodology turned out to be more cost effective.This estimate was made in 1988 and readjusted in 1993 and affirmed by ex-postinvestigation five years later.Accordingly, it was presumed that fox termination brought about a transient diminishment in the event of the infection while oral immunization turned out to be equipped for wiping out rabies even in circumstances in which fox population was growing.Anderson and May [14] formulated a mathematical model based on each time step dynamic which was calculated independently in every cell.Later, Bohrer et al. [15] published a paper on the viability of different rabies spatial immunization designs in a simulated host population.
The research presented by Bohrer [15] stated that, in desert environments, where host population size varies over time, nonuniform spreading of oral rabies vaccination may, under certain circumstances, be more effective than the commonly used uniform spread.The viability of a nonarbitrary spread of the immunization depends, to some extent, on the dispersal behavior of the carriers.The outcomes likewise exhibit that, in a warm domain in a few high-density regions encompassed by populations with densities below the critical threshold for the spread of the disease, the rabies infection can persist.
Levin et al. [16] also presented a model for the immune responses to rabies virus in bats.Coyne et al. [17] proposed an SEIR model, which was also used in a study predicting the local dynamics of rabies among raccoons in the United States.Childs et al. [18] also researched rabies epidemics in raccoons with a seasonal birth pulse, using optimal control of an SEIRS model which describes the population dynamics.Hampson et al. [19] also noted that rabies epidemic cycles have a period of 3-6 years in dog populations in Africa, so they built a susceptible, exposed, infectious, and vaccinate model with an intervention response variable, which showed significant synchrony.
Carroll et al. [20] also used compartmental models to describe rabies epidemiology in dog populations and explored three control methods: vaccination, vaccination pulse fertility control, and culling.An ordinary differential equation model was used to characterize the transmission dynamics of rabies between humans and dogs by [21,22].The work by Zinsstag et al. [23] further extended the existing models on rabies transmission between dogs to include dogto-human transmission and concluded that human postexposure prophylaxis (PEP) with a dog vaccination campaign was the more cost effective in controlling the disease in the long run.Furthermore, Ding et al. [24] formulated an epidemic model for rabies in raccoons with discrete time and spatial features.Their goal was to analyze the strategies for optimal distribution of vaccine baits to minimize the spread of the disease and the cost of carrying out the control.[25] show that culling could be more effective than vaccination, given the same efficacy of control, but Tchuenche and Bauch suggest that culling could be counterproductive, for some parameter values (see [26]).

Smith and Cheeseman
The work in [27,28] also presented a mathematical model of rabies transmission in dogs and from the dog population to the human population in China.Their study did not consider the use optimal control analysis to the study of the rabies virus in dogs and from the dog population to the human population.Furthermore, the insightful work of Wiraningsih et al. [29] studied the stability analysis of a rabies model with vaccination effect and culling in dogs, where they introduced postexposure prophylaxis to a rabies transmission model, but the paper did not consider the noneffectiveness of the pre-and postprophylaxis on the susceptible humans and exposed humans and that of the dog population and the use of optimal control analysis.Therefore, motivated by the research predictions of the global alliance of rabies control [30] and the work mention above, we seek to adjust the model presented in [27][28][29], by formulating an optimal control model, so as to ascertain an optimal way of controlling rabies transmission in dogs and from the dog population to the human population taking into account the noneffectiveness (failure) of vaccination and treatment.
The paper is petition as follows.Section 2 contains the model formulation, mathematical assumptions, the mathematical flowchart, and the model equations.Section 3 contains the model analysis, invariant region, equilibrium points, basic reproduction number R 0 , and the stability analysis of the equilibria.In Section 4 we present the parameter values leading to numerical values of the basic reproduction number R 0 , the herd immunity threshold and sensitivity analysis using Latin hypercube sampling (LHS), and some numerical plots.Section 5 contains the objective functional and the optimality system of the model.Finally, Sections 6 and 7 contain discussion and conclusion, respectively.

Model Formulation
We present two subpopulation transmission models of rabies virus in dogs and that of the human population (see Figure 1), based on the work presented in [27][28][29].The dog population has a total of four compartments.The compartments represent the susceptible dogs,   (), exposed dogs,   , infected dogs,   (), and partially immune dogs,   ().Thus, the total dog population is   () =   () +   () +   () +   ().The human population also has four compartments representing susceptible humans,   (), exposed humans,   (), infected humans,   (), and partially immune humans,   ().Thus, the total human population is   () =   () +   () +   () +   ().It is assumed that there is no human to human transmission of the rabies virus in the human submodel (see [29]).In the dog submodel, it is assumed that there is a direct transmission of the rabies virus from one dog to the other and from the infected dog compartment to the susceptible human population.It is further assumed that the susceptible dog population,   (), is increased by recruitment at a rate   and   is the birth or immigration rate into the susceptible human population,   ().It is assumed that the transmission and contact rate of the rabid dog into the dog compartment is   .Suppose that ]  represents the control strategy due to public education and vaccination in the dogs compartment; then the transmission dynamics become (1 − ]  )      , where (1 − ]  ) is the noneffectiveness (failure) of the vaccine.

Dog population dynamics
It is also assumed that the contact rate of infectious dogs to the human population is   .Similarly, administrating vaccination to the susceptible humans the progression rate of the susceptible humans to the exposed stage becomes , where ]  is the preexposure prophylaxis (vaccination), (1 − ]  ) represents the failure of the preexposure prophylaxis in the human compartment.Furthermore, administrating postexposure prophylaxis (treatment) to affected humans at the rate   decreases the progression rate of the rabies virus, at the exposed class to the infectious class as (1 −   )      , where (1 −   ) is the failure rate of the postexposure prophylaxis and     represents the rate at which exposed humans progress to the infected compartment [27].The rate of losing immunity in both compartments is represented by   and   , respectively.
The exposed humans without clinical rabies that move back to the susceptible population are denoted by the rate     .The natural death rate of dogs is   , and   denotes the mortality rate of humans (natural death rate),   represents the death rate associated with rabies infection in dogs, and   represents the disease induce death in humans.The rate at which exposed dogs die due to culling is   , and   represents the rate at which exposed dogs without clinical rabies move back to the susceptible dog compartment.Subsequently, using the idea presented in [29], we assumed that the exposed dogs are treated or quarantined by their owners at the rate   ; this implies that (1 −   )    is the progression rate of the exposed dogs to the infectious compartment, where (1 −   ) is the failure of the treatment or quarantined strategy, and     denotes those exposed dogs that develop clinical rabies [27].Figure 1 shows the mathematical dynamics of the rabies virus in both compartments.

Model Analysis
Model system (1) will be studied in a biological feasible region as outlined below.+ for all  > 0. We want to show that the region Ω is positively invariant, so that it becomes sufficient to look at the dynamics of model system (1), given that () =   () +   () +   () +   () , where   () is the total population of dogs at any time () and   () is total population of humans at any time ().Equation (2) gives which yields Similarly (3) gives Now, assuming that there are no disease induced death rate and culling effect in the dogs' compartment, it implies that ( 5) and ( 6) become Suppose   / ≤ 0,   / ≤ 0,   ≤   /  , and   ≤   /  , and then imposing the theorem proposed in [32] on differential inequality results in 0 ≤  D ≤   /  and 0 ≤   ≤   /  .Therefore (7) becomes Solve ( 8) and ( 9) using the integrating factor (IF) method.Thus / + () = ,  =  ∫ () .After some algebraic manipulation the feasible solution of the dogs' population in model system (1) is in the region Similarly the human population follows suit, and from (9) this implies that the feasible solution of the human population of model system (1) is in the region Therefore, the feasible solutions are contained in Ω.Thus Ω = Ω  × Ω  .From the standard comparison theorem used on differential inequality in [33], it implies that Hence, the total dog population size   () →   /  as  → ∞.Similarly, the total human population size   () →   /  as  → ∞.This means that the infected state variables (  ,   ,   ,   ) of the two populations tend to zero as time goes to infinity.Therefore, the region Ω is pulling (attracting) all the solutions in R 8 + .This gives the feasible solution set of model system (1) as Hence, (1) is mathematically well posed and epidemiologically meaningful.

Basic Reproduction
Here, the basic reproduction number (R 0 ) measures the average number of new infections produced by one infected dog in a completely susceptible (dog and human) population (see also [34]).Now taking   ,   ,   , and   as our infected compartments gives where  1 =   /,  2 =   /,  3 =   /, and  4 =   /.Now, using the next generation matrix operator  =  −1 and the Jacobian matrix as described in [34], results in Using the fact that  =  −  gives  and  evaluated at E 0 as where the element in matrix  constitutes the new infection terms, while that of matrix  constitutes the new transfer of infection terms from one compartment to another.Now, splitting matrix  into four 2 × 2 submatrices and finding its corresponding inverses result in  =  −1 , given by Letting Finding the matrix determinant of ( 22 This gives a characteristic equation of the form  3 ( − ) = 0; solving the characteristic polynomial results in the following eigenvalues:   = [0, 0, 0, ].The basic reproduction number R 0 is the spectral radius (largest eigenvalue) ( −1 ), also defined as the dominant eigenvalue of  −1 . Therefore, Remark 2. R 0 contains the secondary infection produced by the infectious compartment of dogs (in the presence of preexposure prophylaxis (vaccination), postexposure prophylaxis (treatment/quarantine), and culling of exposed dogs).When R 0 < 1, the infection gradually leaves the dog compartment, but when R 0 > 1, the rabies virus remains in the dog compartments for a longer time, thereby increasing the rate at which the susceptible dogs and humans get infected by a rabid dog.

3.3.
Endemic Equilibrium E 1 .The endemic equilibrium is given as Note that if R 0 = 1, it results in the disease-free equilibrium; if R 0 > 1, then there exists a unique endemic equilibrium; if R 0 < 1, then there exist two endemic equilibriums.

Stability
Analysis of E 0 .Linearizing (1) at E 0 and subtracting eigenvalue  along the main diagonal yield where Simplifying matrix J(E 0 ) gives where From (28) the four characteristic factors that are negative are where The other four characteristic factors can be obtained using the Routh-Hurwitz criterion.Routh-Hurwitz stability criterion is a test to ascertain the nature of the eigenvalues.If the roots of the polynomial are all positive, then the polynomial has a negative real part [35,36].The remaining four characteristic eigenvalues are obtained as follows: Hence, simplifying the coefficient of the above characteristic polynomial in (31) yields Therefore, from the Routh-Hurwitz criterion of order four, it implies that the conditions,  11 > 0,  12 > 0,  13 > 0,  14 > 0, and  11  12  13 >  2  13 +  2 11  14 , are satisfied if R 0 < 1.Hence, the disease-free equilibrium E 0 is locally asymptotically stable when R 0 < 1 (see [37]).

Global Stability of E 0
Theorem 3. The disease-free equilibrium E 0 of model ( 1) is globally asymptotically stable if R 0 ≤ 1 and unstable if R 0 > 1.

Global Stability of Endemic Equilibrium E 1
Theorem 4. The endemic equilibrium E 1 of model ( 1) is globally asymptotically stable whenever R 0 > 1.
Proof.Suppose R 0 > 1; then the existence of the endemic equilibrium point is assured.Using the common quadratic Lyapunov function as illustrated in [40], we consider a Lyapunov function with the following candidate: Now, differentiating (42) along the solution curve of ( This also implies that This shows that V/ is negative and V/ = 0, if and only if Additionally every solution of (1) with the initial conditions approaches E 1 as  → ∞ (see [38,39]); therefore, the largest compact invariant set in {(  ,   ,   ,   ,   ,   ,   ,   ) ∈ Ω : V/ ≤ 0} is the singleton set {E 1 }.Therefore, from Lasalle's invariant principle [41], it implies that the endemic equilibrium E 1 is globally asymptotically stable in Ω whenever R 0 > 1.

Numerical Analysis
Considering the parameter values in Table 1, we will ascertain the numerical importance of our analysis.

Different Scenarios of the Basic Reproduction Number
R 0 .We shall denote R 0 without pre-and postexposure prophylaxis (treatment) as R * 0 and R 0 without preexposure prophylaxis and culling as R * * 0 and the R 0 without postexposure prophylaxis (treatment) and culling as R * * * 0 .Therefore, using the parameter values in Table 1, R * 0 , R * * 0 , and R * * * 0 are given as follows: Therefore, from the above calculations it indicates that the best way in reducing or minimizing the rabies virus in the dogs compartment is to use more of preexposure prophylaxis (vaccination).4.1.1.Herd Immunity Threshold  1 .Therefore, from the above numerical values, we are motivated to know the number of humans or dogs that should be vaccinated when R * 0 = 3.027.
This shows that if R * 0 = 3.027, then 66% of individuals and dogs should receive vaccination.

Sensitivity Analysis.
To determine parameters that contribute most to the rabies transmission, we used two sensitivity analysis approach: the normalised forward sensitivity index as presented in [37] and the Latin hypercube sampling as described in [42].To determine the dependence of parameters in R 0 , using a sampling size,  = 1000, the partial rank correction coefficients (PRCC) value of the ten parameters in R 0 are shown in Figure 2(a).The longer the bar in Figure 2(a) suggests that the statistical influence of those parameters to changes in R 0 is high.Also, using the normalised forward sensitivity index gives the following values and the nature of their signs in Table 2, based on the parameter value given in Table 1.The plus sign or minus sign signifies that the influence is positive or negative, respectively [42], Therefore, from Table 2 it shows that an addition or a reduction in the values of   ,   ,   , and   will have an increase or a decrease in the spread of the rabies virus.For

example, Γ
R 0 = 1 indicates that increasing or reducing the transmission rate by 5% may increase or reduce the number of secondary infection by 5%.The negative sign in Table 2 will have a reduction in the basic reproduction number, R 0 , when the values of those parameters are increased, and a reduction in the values of   , ]  ,   ,   , and   will lead to an increase in the number of secondary infections.
The Latin hypercube sampling (LHS) in Figure 2(a) shows that   ,   ,   , and   have a minimal influence on the rate at which the rabies virus is spread.The Latin hypercube sampling (LHS) plots for the ten parameters in R 0 show that culling of exposed dogs does not actually minimize the spread of rabies as compared to vaccination of susceptible dogs.Figure 2(a) also shows that the most influential parameter in spreading the infection is   followed by   .Figure 2(c) shows that an increase in the basic reproduction number will contribute to a high level of secondary infection in the human population.Similarly, Figure 2(a) shows that vaccination of dogs ]  is the most effective way of controlling the rabies virus in the dog population as compared to the treatment/quarantine of exposed dogs,   .Figure 3(a) gives the contour nature of ]  and   , which shows a more saturated effect on the basic reproduction number.Figure 3(b) shows that   and   have a positive relation with the basic   reproduction number R 0 .Therefore, an increase in   and   will have a direct increase in the spread of the rabies virus.Figure 2(b) indicates that with a high number of recruitment of dogs into the susceptible dog's compartment will have a corresponding high increase in the number of infected humans.Figure 2(d) demonstrates that a high number of infected dogs in the compartment will lead to an increase in the number of infected humans.Figure 3(c) shows that a high increase in the number of disease induce death rate and natural death rate will have a negative reflection on R 0 ; biologically, we would not recommend this approach in minimizing the spread of the disease, since an increase in both   and   may result in a high rate of the disease in the human population, even though   and   naturally reduce the number of susceptible and infected dogs in the population.Finally, Figure 3(d) shows the 3D plot of Figure 3(a).

Objective Functional
Given that () ∈  ∈ R  is a state variable of model system (1) and () ∈  ∈ R  are the control variables at any time () with  (0) ≤  ≤  () , then an optimal control problem consists of finding a piecewise continuous control () and its corresponding state ().This optimizes the cost functional [(), ()] using Pontryagin's maximum principle [43].Therefore we set the following likelihood control strategies: (1)  1 = ]  is the control effort aimed at increasing the immunity of susceptible dogs (preexposed prophylaxis).
(3)  3 = ]  is the control effort aimed at increasing the immunity of susceptible humans (preexposure prophylaxis).
Our goal is to seek optimal controls such as ] *  ,  *  , ] *  , and  *  that minimize the objective functional: Therefore, (51) is subject to From (51) the quantities  1 and  2 denote the weight constants of the exposed classes and  3 and  4 are the weight of the infectious classes, respectively.1, 2, 3, 4 are the weight constants for the dog and human controls.1] 2  , 2 (60)

Numerical Simulations of the Optimality System.
To determine the control strategies ]  ,   , ]  , and   , as given in the objective functional, we began an iteration of the model until convergence is achieved.The results of the simulation of the control strategies are displayed below.We consider equal weights of ( 1 = 1,  2 = 1,  3 = 1,  4 = 1) for both exposed and infected classes.We varied the cost associated with the objective functional, which indicate that, with low cost of vaccination, the rate at which individuals will seek for vaccination of their susceptible dogs will increase, and this could result in low transmission of rabies in a heterogeneous population.We consider the various cost of preexposure prophylaxis and postexposure prophylaxis to be (1 = 1, Exposed dogs Exposed humans 2 = 4, 3 = 1, 4 = 4).We found that the optimal time in controlling the infection using preexposure prophylaxis in dogs is much better than using postexposure prophylaxis in dogs, as shown by the trajectories of the red line and blue line in Figure 4, respectively.The blue line in Figure 4 indicates that applying postexposure prophylaxis will considerably take a longer time in controlling of rabies in dogs.The green line in Figure 4 signifies that preexposure prophylaxis in humans increases the immunity levels of humans and hence reduces the rate at which individuals move to the infected stage.prophylaxis in the dog population will result in a high of the rabies infection in the human population.Figure 6(a) also shows that combining pre-and postexposure prophylaxis (vaccination and treatment/quarantine) in the dog compartment will reduce the spread of the rabies virus, thereby reducing the using of pre-and postexposure prophylaxis (vaccination and treatment) in humans.Figure 6(b) indicates that a rapid use of pre-and postexposure prophylaxis in the human population will reduce the number of rabies deaths in the human population.Figure 7 shows the simulation effects of applying both controls on the model.Figure 7(a) shows that, with the use of the optimal control strategies, the rate of the infection in the susceptible dogs will reduce significantly.Figures 7(b) and 7(c) show that there is a proportional decrease in the number of exposed and infected dogs when the control measures are applied.Similarly, Figures 7(e) and 7(f) show a significant decrease in the number of infected and exposed humans when the control measures are applied.Figure 7(d) shows that there is a proportional increase in the number of recovered dogs when the control measures are applied.Finally, Figures 8(a)-8(h) show the simulation effect of corresponding adjoint functions.

Discussion
The numerical simulations of the resulting optimality system show that, during the case where it is more expensive to vaccinate than treatment, more resources should be invested in treating affected individuals until the disease prevalence begins to fall.This option, however, does not reduce the number of individuals expose to the disease quickly enough, thus resulting in an overall increase in the infected human population.On the other hand, if it is more expensive to treat than to vaccinate, then more susceptible dogs should be vaccinated, so as to lower the rate at which newborn dogs get infected.Nevertheless, in the case where both measures are equally expensive, the simulation shows that the optimal way to drive the epidemic towards eradication within any specified period is to use more preexposure prophylaxis in both compartments.

Conclusion
We studied an optimal control model of rabies transmission dynamics in dogs and the best way of reducing death rate of rabies in humans.The stability analysis shows that the disease-free equilibrium is locally and globally asymptotically stable.We also obtained an optimal control solution for the model which predicts that the optimal way of eliminating deaths from canine rabies as projected by the global alliance for rabies control [30] is using more of preexposure prophylaxis in both dogs and humans and public education; however, the results show that the effective and optimal consideration of preexposure prophylaxis and postexposure prophylaxis in humans without an optimal use of vaccination in the dog population is not beneficial if total elimination of the disease is desirable in Africa and Asia.Any combination strategy which involves vaccination in the dogs' population gives a better result and hence it may be beneficial in eliminating the disease in Asia, Africa, and Latin America.

Figure 1 :
Figure 1: Optimal control model of rabies transmission dynamics.
) and denoting it by  give the expression  = |−I|, where I is the identity matrix of a 4 × 4 matrix; thus  =                            −   0                           = 0.

A
LHS plot for the parameters in R 0 Effect of varying recruitment rate on the infected humans Effect of an increase in R 0 on the infected humans Effect of varying the initial infected dog population size on the infected humans

Figure 2 :
Figure 2: The graphical representation of some parameters in R 0 and the effect of varying some initial state values on the model.
The contour plot of ]  and   to R 0 The 3D plot of R 0 to   and   The 3D plot of R 0 to   and   The 3D plot of R 0 to   and ]

Figure 3 :
Figure 3: The graphical representation of some parameters in R 0 .

Figure 4 :
Figure 4: The simulation effect of the controls.

Figure 5 :
Figure 5: The trajectories of the model with and without pre-and postexposure prophylaxis on exposed humans and that of the exposed dogs.

]]Figure 6 :
Figure 6: The trajectories of the model with and without pre-and postexposure prophylaxis on infected humans and that of the infected dogs.
With optimal control(e) Infected humans with and without control With optimal control (f) Exposed humans with and without control

Figure 7 : 8 ( 6 (
Figure 7: The trajectories of the model with and without optimal control on individual compartments.

Figure 8 :
Figure 8: The trajectories of the model with and without optimal control on individual compartments and corresponding adjoint function.