# Mathematical Modelling of Bacterial Meningitis Transmission Dynamics with Control Measures

Authors: Joshua Kiddy K. Asamoah; Farai Nyabadza; Baba Seidu; Mehar Chand; Hemen Dutta
Journal: Computational and Mathematical Methods in Medicine
 (2018)

---

## Abstract

Vaccination and treatment are the most effective ways of controlling the transmission of most infectious diseases. While vaccination helps susceptible individuals to build either a long-term immunity or short-term immunity, treatment reduces the number of disease-induced deaths and the number of infectious individuals in a community/nation. In this paper, a nonlinear deterministic model with time-dependent controls has been proposed to describe the dynamics of bacterial meningitis in a population. The model is shown to exhibit a unique globally asymptotically stable disease-free equilibriumE0, when the effective reproduction number RVT≤1, and a globally asymptotically stable endemic equilibrium E1, when RVT>1; and it exhibits a transcritical bifurcation at RVT=1. Carriers have been shown (by Tornado plot) to have a higher chance of spreading the infection than those with clinical symptoms who will sometimes be bound to bed during the acute phase of the infection. In order to find the best strategy for minimizing the number of carriers and ill individuals and the cost of control implementation, an optimal control problem is set up by defining a Lagrangian function L to be minimized subject to the proposed model. Numerical simulation of the optimal problem demonstrates that the best strategy to control bacterial meningitis is to combine vaccination with other interventions (such as treatment and public health education). Additionally, this research suggests that stakeholders should press hard for the production of existing/new vaccines and antibiotics and their disbursement to areas that are most affected by bacterial meningitis, especially Sub-Saharan Africa; furthermore, individuals who live in communities where the environment is relatively warm (hot/moisture) are advised to go for vaccination against bacterial meningitis.


---

## Body

## 1. Introduction

Meningitis is an inflammation of the meninges which are membranes that surround the spinal cord and the brain [1]. It is often caused by viruses, bacteria, and protozoa. Bacterial meningitis is common in children and young adults. This disease mostly spreads in communities/societies that live in close quarters (e.g., police staff, police cells, college students, military staff, and prisons) [2]. Bacterial meningitis is generally caused by germs such asListeria monocytogenes,Streptococcus pneumoniae, Group BStreptococcus,Neisseria meningitidis, andHaemophilus influenzae, which spreads from one person to another [3]. This infection varies by age groups: Group BStreptococcus,Streptococcus pneumoniae,Listeria monocytogenes, andEscherichia coli are mostly found in newborn babies;Streptococcus pneumoniae,Neisseria meningitidis,Haemophilus influenzae type b (Hib), and Group BStreptococcus are common in babies and children;Neisseria meningitidis andStreptococcus pneumoniae are predominant in teens and young adults; andStreptococcus pneumoniae,Neisseria meningitidis,Haemophilus influenzae type b (Hib), Group BStreptococcus, andListeria monocytogenes are commonly found in older adults [3]. Bacterial meningitis is characterized by intense headache and fever, vomiting, sensitivity to light, and stiff neck, which result in convulsion, delirium, and death.It is estimated that meningococcal meningitis causes over 10,000 deaths annually in Sub-Saharan Africa [4]. About 4,100 cases of bacterial meningitis occurred between 2003 and 2007 in the United States [3, 5]. Between 5% to 40% of children and 20% to 50% of adults with this condition die [6]. Infections from bacterial meningitis can cause permanent disabilities such as brain damage, hearing loss, and learning disabilities [3]. The illness of bacterial meningitis becomes worse when symptoms are not detected early enough; even with proper treatment, the individual could die [2].Prevention of bacterial meningitis can be achieved through vaccination and/or preventing contact with infectious individuals. Vaccination is the most effective way of protecting children against certain types of bacterial meningitis [3]. Vaccines that can prevent meningitis includeHaemophilus influenza type B (Hib), pneumococcal conjugate, and meningococcal vaccine [6, 7]. The conjugate meningitis A vaccine, MenAfrivac, is recommended to protect people in Sub-Saharan Africa against the most common type, serotype A [8]. In the United States, the primary means of preventing meningococcal meningitis is antimicrobial chemoprophylaxis [9]. Empirical therapy includes ceftriaxone or cefotaxime and vancomycin forStreptococcus pneumoniae [2]. There is a vaccine against meningococcal disease which is 85%–100% effective in preventing four kinds of bacteria (serogroups A, C, Y, and W-135) that cause about 70% of the disease in the United States [2].Trotter and Ramsay [10] outlined some recommendations on the use of conjugate vaccines in Europe based on the experience with meningococcal C conjugate (MCC) vaccines. In areas with limited health infrastructure and resources, there are a number of antibiotics including penicillin, ampicillin, and chloramphenicol that can be used to treat the infection meningitis.Mathematical models have been shown to help increase the understanding of the spread and control of infectious diseases. Martínez et al. [2] studied the spread of meningococcal meningitis with the use of a discrete mathematical model, based on cellular automata where the population was divided into five classes: susceptible, asymptomatic infected, infected with symptoms, carriers, recovered, and died classes. Broutin et al. [11] studied the dynamics of meningococcal meningitis in nine African countries by adopting some mathematical tools to time series analysis and wavelet method, the results of their studies suggest that “international cooperation in Public Health and cross disciplines studies are highly recommended to help in controlling this infectious disease.” Miller and Shahab [12] studied the cost effectiveness of immunisation strategies for the control of epidemic meningococcal meningitis. The research work in [13] gives a detailed description of the use of antibiotics for the prevention and treatment of meningitis infection. Irving et al. [14] used deterministic compartmental models to investigate how well simple model structures with seasonal forcing were able to qualitatively capture the patterns of meningitis infection. They demonstrated that the complex and irregular timing of epidemics could be caused by the interaction of temporary immunity conferred by carriage of the bacteria together with seasonal changes in the transmissibility of infection. Actually, there have been a significant number of studies of various types of Meningitis in Africa and Europe without the use of optimal control analysis (see [15–28]).It is obvious that mathematical modelling has become crucial in investigating the epidemiological behaviour of meningitis. Furthermore, mathematical modelling helps to identify the risk factors for diseases, so as to find out why everyone does not have the same infection uniformly [29].The application of optimal control in disease modelling gives valuable information on how to apply control measures. Through vaccination, treatment, public education, and so forth, many infectious diseases have been controlled [29]. Since the introduction of optimal control theory in disease modelling, there have been a considerable number of studies of infectious diseases using optimal control analysis (see [30–41]). With the significant influence of optimal control theory in disease modelling, this paper presents an optimal control model for bacterial meningitis in the presence of vaccination and treatment due to public health education. The model is qualitatively analyzed and numerically simulated in order to help give policy direction on how to control the spread of the disease.The rest of the paper is organized as follows. Section2 presents the model formulation and analysis. Section 3 presents the analysis of the optimal control problem, leading to the existence and characterization of the control measures. Section 4 contains the numerical simulations and discussion. Section 5 presents the conclusion of the study.


## 2. Model Formulation and Analysis

Adopting the epidemiological studies of a meningitis model as presented in [4], we consider four mutually exclusive compartments to indicate individuals with unique natures (i.e., susceptibles, S(t), carriers, C(t), ill individuals, I(t), and recovered individuals, R(t)) in relation to the disease. It is assumed that the susceptible compartment, S(t), is populated through recruitment at the rate, π (thus migration and/or birth rate), and β is the rate of effective contact of carriers and/or infected (ill) individuals in the susceptible population. The carrier compartment consists of individuals that have the infection and do not show any clinical symptoms but contribute to the spread of the disease. When a susceptible individual is exposed to this infection, that individual can harbour the bacterium for weeks or even months [42]; but, in a normal circumstance, an individual develops symptoms of the infection within 3 to 7 days after exposure [3]. Carriers are assumed to develop clinical symptoms (i.e., move to ill individuals compartment, I(t)) at rate α. Ill individuals who are seriously infected are assumed to have no natural recovery except when given treatment on time. From an epidemiological perspective, individuals in the removed/recovered compartment, R(t), do not attain permanent immunity. After vaccination, immunity develops within 7–10 days and remains effective for approximately 3–5 years [2]. Therefore, it is assumed that the immunity acquired from developing the diseases or carrying the bacteria or through vaccination is of the same intensity (they all lead to the recovered compartment from which people return to the susceptible compartment at a given unique rate θ). Research indicates that carriers may recover naturally from the infection without treatment, and we denote such natural recovery rate as ω. Infected individuals are assumed to die from disease at rate δ. Since natural death is inevitable, μ is assumed to be the natural death rate of individuals in all the compartments. Vaccination and treatment due to public health education have been shown to be strategies of control of diseases. Therefore we introduced this two control measures in the model as u1(t) and u2(t), respectively, here u1; thus vaccination is comprised of both reactive vaccination and preventative vaccination. The effectiveness of both control measures in minimizing the disease is denoted by σ and γ, respectively. If σ=γ=0, it signifies that vaccination and treatment have no effect on the model; if σ=γ=1, it also signifies that vaccination and treatment are perfectly effective (i.e., 100% effectiveness) (see [43]). If 0<σandγ<1, it signifies that both vaccine and treatment are imperfect [44]. In view of this, it is assumed that administering treatment to the ill individuals leads to recovery at rate γu2. It is also assumed that the vaccinated individuals develop partial immunity at rate σu1. From an epidemiological perspective treatment is not given to carriers in real life (since we do not know who carries the bacteria or not), but it is assumed that carriers of meningitis are just like people with the HIV virus who do not know their status unless they go for medical test; hence this paper seeks to encourage individuals to go for regular test of this bacterial disease; therefore it is assumed that a certain portion of carriers could be treated before any symptoms of the infection show up, which results in γu2C as shown in Figure 1. The vaccine is assumed to be imperfect and thus has a failure rate of (1-u1σ). Therefore the force of new infections is given by(1)λ=β1-u1σηC+η1IN.Equation (1) is often referred to as standard incidence rate of new infections, which is normalized by the total population N=S+C+I+R. Table 1 gives a full description of parameters used in the model.Table 1
Description of parameter values used in the model. “Due to lack of relevant data, most of the parameter values are assumed within realistic ranges for a typical scenario in a rural community/society for the purpose of illustration” [44], with most of the parameter values where appreciable per year.
Parameter Description Value/Range Reference π Recruitment rate (migration and/or birth rate) 100–100000 Assumed β Effective contact rate 0.88 [28] θ Loss of immunity 0.04–2 [14, 47] μ Natural death rate 0.02 [14] δ Disease-induced mortality 0.05–0.5 [48] α Rate of progression fromC to I 0.1–0.52 [14] ω Natural recovery rate 0.06–0.2 Assumed η1 Per capita infection rate by ill individuals 0.2–0.95 Assumed η Per capita infection rate by Carriers 0.2–0.85 Assumed u1 Vaccination rate of susceptible individuals against Meningitis 0,1 Assumed u2 Treatment rate for carriers without natural recovery and ill individuals 0,1 Assumed σ Effectiveness of vaccination 0.85–1 [2] γ Effectiveness of Treatment 0.1–0.9 AssumedFigure 1
The flowchart diagram describing bacterial meningitis transmission dynamics within the population. The four circles represent the four compartments of individuals, the movement between the compartments is indicated by the continuous arrows,u1 is a control measure (vaccination), and u2 is the second control measure (treatment), with the consideration that both control measures lies in 0<u1andu2≤1.The set of differential equations and flow diagram corresponding to the bacterial meningitis dynamics and disease pathway with control terms is given in system (2) and Figure 1.(2)dSdt=π-1-u1σβSηC+η1IN-μ+u1σS+θR,dCdt=1-u1σβSηC+η1IN-α+ω+μ+u2γC,dIdt=αC-δ+μ+u2γI,dRdt=u2γI+ω+u2γC+u1σS-θ+μR,S0>0,C0≥0,I0≥0,R0≥0.Adding the equations in model (2) gives the rate of change of total population as(3)dNdt=π-μN-δI≤π-μN.Let(4)Ω=S,C,I,R∈R+4:N≤πμ. System (2) is well-posed with all solutions in Ω remaining in Ω if initial conditions are positive. It can easily be shown that if the initial conditions start outside Ω the solutions tend to Ω.


### 2.1. Equilibrium Points

To obtain the equilibrium points of system (2), we assume that the control measures are time-independent (see [45] for similar analysis).


#### 2.1.1. Disease-Free EquilibriumE0

To obtain the disease-free equilibrium,C(t),  I(t), and the right-hand-side of system (2) are set to zero. If susceptible individuals are assumed to receive vaccination against the disease at a constant rate, then the disease-free equilibrium will be given by(5)E0S0,C0,I0,R0=πθ+μμu1σ+θ+μ,0,0,πu1σμu1σ+θ+μ.


#### 2.1.2. Effective Reproduction NumberRVT

Using the next generation matrix method [46], the effective reproduction number of the bacterial meningitis model with vaccination and treatment is obtained as(6)RVT=RCVT+RIVT,where(7)RCVT=βη1-u1σθ+μα+ω+μ+u2γθ+μ+u1σ,RIVT=βαη11-u1σθ+μα+ω+μ+u2γδ+μ+u2γθ+μ+u1σ.To determine how the two control measures impact on the reproduction number, we make a plot ofRVT, on the u1-u2 plane for arbitrary constant values of model parameters in Figure 2. From the figure, it is shown that u2 decreases with increasing u1. So an increase in vaccination levels decreases the need for treatment. Figure 2 also shows the region in which the vaccination and treatment values should lie for the disease control.Figure 2
Region whereRVT<1 and RVT>1 in the u1-u2 parameter space, with the parameter values β=0.88,η=0.2,ω=0.06,δ=0.03,θ=0.0839,α=0.05,γ=0.2andσ=0.7.


#### 2.1.3. Endemic EquilibriumE1

System (2) can be shown to have a unique endemic equilibrium of the form S∗,C∗,I∗,R∗, where(8)S∗=απμ+θ+ω+u2γδ+μ+u2γθ+θαu2γ-α+ω+μ+u2γδ+μ+u2γθ+μπαRVT-1αμθ+μ+u1σARVT-1+B,C∗=δ+μ+u2γπαRVT-1αARVT-1+B,R∗=απu1σ+δ+μ+u2γω+u2γμ+u1σ+αu2γμ+u1σ-α+ω+μ+u2γδ+μ+u2γu1σπαRVT-1αμθ+μ+u1σARVT-1+B,I∗=παRVT-1ARVT-1+B, with A=(δ+μ+u2γ)(α+ω+μ+u2γ) and B=(δ+μ+u2γ)ω+u2γ-u1σ)+α(u2γ+μ).Remark 1.
(i) IfRVT<1, then system (2) will have only one equilibrium: the disease-free equilibrium.
(ii) IfRVT>1, then system (2) will have two equilibria: the disease-free equilibrium, E0, and the endemic equilibrium, E1.
(iii) The caseRVT=1 is a critical threshold point where the disease-free equilibrium E0 loses its local asymptotic stability. Thus RVT=1 gives the idea of transcritical bifurcation where the stability of system (2) moves between E0 and E1 [37].


### 2.2. Stability Analysis

In analyzing the local stability of the disease-free equilibrium, the Routh-Hurwitz criteria are used, and, for the global stability of the two equilibria, the direct Lyapunov technique is employed.


#### 2.2.1. Local Stability Analysis ofE0

Theorem 2.
The disease-free equilibrium of system (2) is locally asymptotically stable if RVT<1 and unstable if RVT>1.Proof.
Evaluating the Jacobian matrix of system (2) at the disease-free equilibrium gives(9)JE0=-μ+u1σ-1-u1σβηθ+μu1σ+θ+μ-1-u1σβη1θ+μu1σ+θ+μθ01-u1σβηθ+μu1σ+θ+μ-α+ω+μ+u2γ1-u1σβη1θ+μu1σ+θ+μ00α-δ+μ+u2γ0u1σω+u2γu2γ-θ+μ.
The characteristic polynomial of the Jacobian matrixJ(E0) is given by(10)Pλ=λ4+DTλ3+DT1λ2+DT2λ+DT3,where(11)T=A+B+C-RCVT,T1=A+B-RCVT2μ+u1σ+θ+μ+u1σθ+μ+μθ1-RVT,T2=θ+μ1-RVT+A-RCVT-RCVT+μ+u1σ1-RIVT+δ+μ+u2γθ+μ+A+B-RCVTμθ,T3=θμRIVT-1+μ+u1σθ+μRIVT-1,A=1δ+μ+u2γ,B=1α+ω+μ+u2γ,C=μ+u1σ+θ+μα+ω+μ+u2γδ+μ+u2γ,D=α+ω+μ+u2γδ+μ+u2γ.
The Routh-Hurwitz conditions [46] that guarantee that the eigenvalues of the characteristic polynomial in (10) have negative real parts are given by(12)DT>0,DT1>0,DT2>0,DT3>0,D3TT1T2>D2T22+D3T2T3.
These conditions are easily seen to be satisfied whenRVT<1. Thus, the disease-free equilibrium of system (2) is locally asymptotically stable when RVT<1 and unstable when RVT>1. This completes the proof.


#### 2.2.2. Global Stability ofE0

Theorem 3.
The disease-free equilibriumE0 of system (2) is globally asymptotically stable if RVT≤1 and unstable if RVT>1.Proof.
LetV(S,C,I,R), with positive constants, K1 and K2, be a Lyapunov function defined as(13)VS,C,I,R=S-S0-S0ln⁡SS0+K1C+K2I+R-R0-R0ln⁡RR0.
Taking the time derivative of the Lyapunov function we obtain(14)dVdt=1-S0SdSdt+K1dCdt+K2dIdt+1-R0RdRdt,whereK1≥0,K2≥0.
SubstitutingdS/dt,  dC/dt,  dI/dt, and dR/dt in (2) into (14) gives(15)dVdt=1-SS0π-1-u1σβSηC+η1IN-μ+σu1S+θR+K11-u1σβSηC+η1IN-α+ω+μ+u2γC+K2αC-δ+μ+u2γI+1-RR0u2γI+ω+u2γC+σu1S-θ+μR≤K11-u1σβηθ+μu1σ+θ+μC+1-u1σβη1θ+μu1σ+θ+μI-α+ω+μ+u2γC+K2αC-δ+μ+u2γI.
Since(16)S≤S0=πθ+μμu1σ+θ+μ,R≤R0=πu1σμu1σ+θ+μ,N≤πμonΩ,
this implies that(17)dVdt≤1-u1σβηθ+μu1σ+θ+μK1+αK1-α+ω+μ+u2γK2C+1-u1σβη1θ+μu1σ+θ+μK1-δ+μ+u2γK2I.
Equating the coefficient ofI in (17) to zero gives(18)1-u1σβθ+μK1=δ+μ+u2γu1σ+θ+μK2.
ChoosingK1=(δ+μ+u2γ)(u1σ+θ+μ)  and  K2=(1-u1σ)β(θ+μ) and plugging K1 and K2 into (17), we have (19)dVdt≤α+ω+μ+u2γδ+μ+u2γu1σ+θ+μRVT-1C≤0,ifRVT≤1.
AdditionallydV/dt=0 if and only if C=0. Hence, the largest compact invariant set in (S,C,I,R)∈Ω:dV/dt≤0 is the singleton set E0. Therefore, from LaSalle’s invariance principle, we conclude that E0 is globally asymptotically stable in Ω if RVT≤1 [37, 49].


#### 2.2.3. Global Stability ofE1

Theorem 4.
The endemic equilibriumE1 of system (2) is globally asymptotically stable whenever RVT>1.Proof.
SupposeRVT>1, and then the existence of the endemic equilibrium point is assured. Using the common quadratic Lyapunov function(20)Vx1,x2,…,xn=∑i=1nci2xi-xi∗2,
as illustrated in [50], we consider the following candidate Lyapunov function:(21)VS,C,I,R=12S-S∗+C-C∗+I-I∗+R-R∗2.
The time derivative ofV(S,C,I,R) in (21) is given by(22)dVdtS,C,I,R=S-S∗+C-C∗+I-I∗+R-R∗dS+C+I+Rdt.
Plugging the equations in system (2) into (22) yields(23)dVdt=S-S∗+C-C∗+I-I∗+R-R∗π-μS+C+I+R-δI.
Now setting(24)π=μS∗+C∗+I∗+R∗+δI∗,
we have(25)dVdt=S-S∗+C-C∗+I-I∗+R-R∗μS∗+C∗+I∗+R∗+δI∗-μS+C+I+R-δI=S-S∗+C-C∗+I-I∗+R-R∗-μS-S∗-μC-C∗-μI-I∗-μR-R∗-δI-I∗.
Further simplification gives(26)dVdt=-μS-S∗2-μS-S∗C-C∗+I-I∗+R-R∗+δI-I∗-μC-C∗2-μC-C∗S-S∗+I-I∗+R-R∗+δI-I∗-μ+δI-I∗2-μI-I∗C-C∗+R-R∗-μR-R∗2-μR-R∗S-S∗+C-C∗+I-I∗+δI-I∗.
It has therefore been shown thatdV/dt is negative, and additionally at E1 (i.e., if S=S∗, C=C∗, I=I∗, and R=R∗), dV/dt=0. It follows from LaSalle’s invariant principle [51] that all solutions of system (2) approach E1 as t→∞ if RVT>1. Therefore, the endemic equilibrium E1 is globally asymptotically stable in Ω whenever RVT>1 [37, 49]. This completes the proof.


### 2.3. Sensitivity Analysis

Sensitivity analysis is used to determine the response of a model to variations in its parameter values. In the present case, the focus is given to determining how changes in the model parameters impact the effective reproduction number. This is done through the normalized forward-sensitivity index. We also use the Latin hypercube sampling and the partial rank correlation coefficients (PRCC) to plot scatter diagrams and Tornado plots to determine the relative importance of the parameters inRVT for the disease transmission and prevalence (see also [52]).Definition 5.
The normalized forward-sensitivity index ofRVT to any parameter, say ρ, as given in [46] can be defined as(27)ΓRVTρ=∂RVT∂ρρRVT.The sensitivity indexes ofRVT with respect to its parameters are computed as follows:(28)ΓRVTβ=∂RVT∂ββRVT=1,ΓRVTη=∂RVT∂ηηRVT=δ+μ+u2γηηδ+μ+u2γ+η1α,ΓRVTη1=∂RVT∂η1η1RVT=αη1αη1+δ+μ+u2γ,ΓRVTω=∂RVT∂ωωRVT=-ωα+ω+μ+u2γ,ΓRVTσ=∂RVT∂σσRVT=-u1σu1σ+θ+μ.Similarly, we can compute the sensitivity indexes ofRVT with respect to the remaining parameters in RVT, in the same manner. Using the parameter values γ=0.4, σ=1, η=0.35, ω=0.2, δ=0.1, β=0.88, α=0.2, θ=0.0839, and μ=0.02, with u1=0.5 and u2=0.5, the sensitivity indexes of RVT are shown in Table 2.Table 2
Sensitivity indexes ofRVT to the parameters in (6).
Parameter Description Sensitivity index β Effective contact rate 1 η Per capita infection rate by Carriers 0.6154 η1 Per capita infection rate by ill individuals 0.3846 θ Loss of immunity 0.6686 α Rate of falling ill 0.0620 μ Natural death rate 0.1031 δ Disease-induced death rate −0.1202 ω Natural recovery rate −0.3226 u1 Vaccination rate −1.8280 σ Effectiveness of vaccination −0.8280 u2 Treatment rate −0.5630 γ Effectiveness of treatment −0.5230The corresponding Tornado plots based on a random sample of1000 points for the twelve parameters in RVT are shown in Figure 3(a). The positive values in Table 2 show a promotion of the propagation of the disease. Therefore an increase in the values of β, η, η1, θ, and α will have an increase in the spread of the disease. For example, ΓRVTβ=1 indicates that increasing the effective contact rate by 10% increases the number of secondary infections by 10%. The negative values in Table 2 indicate a reduction in the effective reproduction number RVT if the values of the corresponding parameters are increased. Thus, a reduction in the values of vaccination u1, treatment u2, and natural recovery ω will lead to an increase in the number of secondary infections in the population.Figure 3
Sensitivity plots.
(a)
The Tornado plots for the elven parameters inRVT (b)
The effect of effective contact rate onRVT (c)
The effect of vaccination onRVT (d)
The effect of treatment onRVTFigure3(a) shows the Tornado plots for the twelve parameters in RVT. It can be seen that, in controlling the spread of bacterial meningitis in a population, more susceptible individuals should be given vaccination. Figure 3(a) also suggests that carriers are likely to have more contacts with the susceptible population than the ill individuals who will typically be bound to their beds during the acute phase of the disease. Therefore, the probability of ill individuals transmitting the infections to susceptibles may be lower than that of carriers who are able to mix well with others within the population. Figures 3(b), 3(c), and 3(d) show the regression plots of effective contact rate (β), vaccination rate (u1), and treatment rate (u2), respectively. Figure 3(b) shows that transmission rate has a positive correlation in the spread of bacterial meningitis. Figure 3(c) shows that vaccination has a negative correlation in the spread of bacterial meningitis and hence vaccination increases the immunity of individuals against the meningitis infection, thereby reducing the spread of the infection. Figure 3(d) shows that individuals who receive treatment after being infected with bacterial meningitis have a higher chance of recovery and that reduces the spread of the infection and deaths due to bacterial meningitis.


## 2.1. Equilibrium Points

To obtain the equilibrium points of system (2), we assume that the control measures are time-independent (see [45] for similar analysis).


### 2.1.1. Disease-Free EquilibriumE0

To obtain the disease-free equilibrium,C(t),  I(t), and the right-hand-side of system (2) are set to zero. If susceptible individuals are assumed to receive vaccination against the disease at a constant rate, then the disease-free equilibrium will be given by(5)E0S0,C0,I0,R0=πθ+μμu1σ+θ+μ,0,0,πu1σμu1σ+θ+μ.


### 2.1.2. Effective Reproduction NumberRVT

Using the next generation matrix method [46], the effective reproduction number of the bacterial meningitis model with vaccination and treatment is obtained as(6)RVT=RCVT+RIVT,where(7)RCVT=βη1-u1σθ+μα+ω+μ+u2γθ+μ+u1σ,RIVT=βαη11-u1σθ+μα+ω+μ+u2γδ+μ+u2γθ+μ+u1σ.To determine how the two control measures impact on the reproduction number, we make a plot ofRVT, on the u1-u2 plane for arbitrary constant values of model parameters in Figure 2. From the figure, it is shown that u2 decreases with increasing u1. So an increase in vaccination levels decreases the need for treatment. Figure 2 also shows the region in which the vaccination and treatment values should lie for the disease control.Figure 2
Region whereRVT<1 and RVT>1 in the u1-u2 parameter space, with the parameter values β=0.88,η=0.2,ω=0.06,δ=0.03,θ=0.0839,α=0.05,γ=0.2andσ=0.7.


### 2.1.3. Endemic EquilibriumE1

System (2) can be shown to have a unique endemic equilibrium of the form S∗,C∗,I∗,R∗, where(8)S∗=απμ+θ+ω+u2γδ+μ+u2γθ+θαu2γ-α+ω+μ+u2γδ+μ+u2γθ+μπαRVT-1αμθ+μ+u1σARVT-1+B,C∗=δ+μ+u2γπαRVT-1αARVT-1+B,R∗=απu1σ+δ+μ+u2γω+u2γμ+u1σ+αu2γμ+u1σ-α+ω+μ+u2γδ+μ+u2γu1σπαRVT-1αμθ+μ+u1σARVT-1+B,I∗=παRVT-1ARVT-1+B, with A=(δ+μ+u2γ)(α+ω+μ+u2γ) and B=(δ+μ+u2γ)ω+u2γ-u1σ)+α(u2γ+μ).Remark 1.
(i) IfRVT<1, then system (2) will have only one equilibrium: the disease-free equilibrium.
(ii) IfRVT>1, then system (2) will have two equilibria: the disease-free equilibrium, E0, and the endemic equilibrium, E1.
(iii) The caseRVT=1 is a critical threshold point where the disease-free equilibrium E0 loses its local asymptotic stability. Thus RVT=1 gives the idea of transcritical bifurcation where the stability of system (2) moves between E0 and E1 [37].


## 2.1.1. Disease-Free EquilibriumE0

To obtain the disease-free equilibrium,C(t),  I(t), and the right-hand-side of system (2) are set to zero. If susceptible individuals are assumed to receive vaccination against the disease at a constant rate, then the disease-free equilibrium will be given by(5)E0S0,C0,I0,R0=πθ+μμu1σ+θ+μ,0,0,πu1σμu1σ+θ+μ.


## 2.1.2. Effective Reproduction NumberRVT

Using the next generation matrix method [46], the effective reproduction number of the bacterial meningitis model with vaccination and treatment is obtained as(6)RVT=RCVT+RIVT,where(7)RCVT=βη1-u1σθ+μα+ω+μ+u2γθ+μ+u1σ,RIVT=βαη11-u1σθ+μα+ω+μ+u2γδ+μ+u2γθ+μ+u1σ.To determine how the two control measures impact on the reproduction number, we make a plot ofRVT, on the u1-u2 plane for arbitrary constant values of model parameters in Figure 2. From the figure, it is shown that u2 decreases with increasing u1. So an increase in vaccination levels decreases the need for treatment. Figure 2 also shows the region in which the vaccination and treatment values should lie for the disease control.Figure 2
Region whereRVT<1 and RVT>1 in the u1-u2 parameter space, with the parameter values β=0.88,η=0.2,ω=0.06,δ=0.03,θ=0.0839,α=0.05,γ=0.2andσ=0.7.


## 2.1.3. Endemic EquilibriumE1

System (2) can be shown to have a unique endemic equilibrium of the form S∗,C∗,I∗,R∗, where(8)S∗=απμ+θ+ω+u2γδ+μ+u2γθ+θαu2γ-α+ω+μ+u2γδ+μ+u2γθ+μπαRVT-1αμθ+μ+u1σARVT-1+B,C∗=δ+μ+u2γπαRVT-1αARVT-1+B,R∗=απu1σ+δ+μ+u2γω+u2γμ+u1σ+αu2γμ+u1σ-α+ω+μ+u2γδ+μ+u2γu1σπαRVT-1αμθ+μ+u1σARVT-1+B,I∗=παRVT-1ARVT-1+B, with A=(δ+μ+u2γ)(α+ω+μ+u2γ) and B=(δ+μ+u2γ)ω+u2γ-u1σ)+α(u2γ+μ).Remark 1.
(i) IfRVT<1, then system (2) will have only one equilibrium: the disease-free equilibrium.
(ii) IfRVT>1, then system (2) will have two equilibria: the disease-free equilibrium, E0, and the endemic equilibrium, E1.
(iii) The caseRVT=1 is a critical threshold point where the disease-free equilibrium E0 loses its local asymptotic stability. Thus RVT=1 gives the idea of transcritical bifurcation where the stability of system (2) moves between E0 and E1 [37].


## 2.2. Stability Analysis

In analyzing the local stability of the disease-free equilibrium, the Routh-Hurwitz criteria are used, and, for the global stability of the two equilibria, the direct Lyapunov technique is employed.


### 2.2.1. Local Stability Analysis ofE0

Theorem 2.
The disease-free equilibrium of system (2) is locally asymptotically stable if RVT<1 and unstable if RVT>1.Proof.
Evaluating the Jacobian matrix of system (2) at the disease-free equilibrium gives(9)JE0=-μ+u1σ-1-u1σβηθ+μu1σ+θ+μ-1-u1σβη1θ+μu1σ+θ+μθ01-u1σβηθ+μu1σ+θ+μ-α+ω+μ+u2γ1-u1σβη1θ+μu1σ+θ+μ00α-δ+μ+u2γ0u1σω+u2γu2γ-θ+μ.
The characteristic polynomial of the Jacobian matrixJ(E0) is given by(10)Pλ=λ4+DTλ3+DT1λ2+DT2λ+DT3,where(11)T=A+B+C-RCVT,T1=A+B-RCVT2μ+u1σ+θ+μ+u1σθ+μ+μθ1-RVT,T2=θ+μ1-RVT+A-RCVT-RCVT+μ+u1σ1-RIVT+δ+μ+u2γθ+μ+A+B-RCVTμθ,T3=θμRIVT-1+μ+u1σθ+μRIVT-1,A=1δ+μ+u2γ,B=1α+ω+μ+u2γ,C=μ+u1σ+θ+μα+ω+μ+u2γδ+μ+u2γ,D=α+ω+μ+u2γδ+μ+u2γ.
The Routh-Hurwitz conditions [46] that guarantee that the eigenvalues of the characteristic polynomial in (10) have negative real parts are given by(12)DT>0,DT1>0,DT2>0,DT3>0,D3TT1T2>D2T22+D3T2T3.
These conditions are easily seen to be satisfied whenRVT<1. Thus, the disease-free equilibrium of system (2) is locally asymptotically stable when RVT<1 and unstable when RVT>1. This completes the proof.


### 2.2.2. Global Stability ofE0

Theorem 3.
The disease-free equilibriumE0 of system (2) is globally asymptotically stable if RVT≤1 and unstable if RVT>1.Proof.
LetV(S,C,I,R), with positive constants, K1 and K2, be a Lyapunov function defined as(13)VS,C,I,R=S-S0-S0ln⁡SS0+K1C+K2I+R-R0-R0ln⁡RR0.
Taking the time derivative of the Lyapunov function we obtain(14)dVdt=1-S0SdSdt+K1dCdt+K2dIdt+1-R0RdRdt,whereK1≥0,K2≥0.
SubstitutingdS/dt,  dC/dt,  dI/dt, and dR/dt in (2) into (14) gives(15)dVdt=1-SS0π-1-u1σβSηC+η1IN-μ+σu1S+θR+K11-u1σβSηC+η1IN-α+ω+μ+u2γC+K2αC-δ+μ+u2γI+1-RR0u2γI+ω+u2γC+σu1S-θ+μR≤K11-u1σβηθ+μu1σ+θ+μC+1-u1σβη1θ+μu1σ+θ+μI-α+ω+μ+u2γC+K2αC-δ+μ+u2γI.
Since(16)S≤S0=πθ+μμu1σ+θ+μ,R≤R0=πu1σμu1σ+θ+μ,N≤πμonΩ,
this implies that(17)dVdt≤1-u1σβηθ+μu1σ+θ+μK1+αK1-α+ω+μ+u2γK2C+1-u1σβη1θ+μu1σ+θ+μK1-δ+μ+u2γK2I.
Equating the coefficient ofI in (17) to zero gives(18)1-u1σβθ+μK1=δ+μ+u2γu1σ+θ+μK2.
ChoosingK1=(δ+μ+u2γ)(u1σ+θ+μ)  and  K2=(1-u1σ)β(θ+μ) and plugging K1 and K2 into (17), we have (19)dVdt≤α+ω+μ+u2γδ+μ+u2γu1σ+θ+μRVT-1C≤0,ifRVT≤1.
AdditionallydV/dt=0 if and only if C=0. Hence, the largest compact invariant set in (S,C,I,R)∈Ω:dV/dt≤0 is the singleton set E0. Therefore, from LaSalle’s invariance principle, we conclude that E0 is globally asymptotically stable in Ω if RVT≤1 [37, 49].


### 2.2.3. Global Stability ofE1

Theorem 4.
The endemic equilibriumE1 of system (2) is globally asymptotically stable whenever RVT>1.Proof.
SupposeRVT>1, and then the existence of the endemic equilibrium point is assured. Using the common quadratic Lyapunov function(20)Vx1,x2,…,xn=∑i=1nci2xi-xi∗2,
as illustrated in [50], we consider the following candidate Lyapunov function:(21)VS,C,I,R=12S-S∗+C-C∗+I-I∗+R-R∗2.
The time derivative ofV(S,C,I,R) in (21) is given by(22)dVdtS,C,I,R=S-S∗+C-C∗+I-I∗+R-R∗dS+C+I+Rdt.
Plugging the equations in system (2) into (22) yields(23)dVdt=S-S∗+C-C∗+I-I∗+R-R∗π-μS+C+I+R-δI.
Now setting(24)π=μS∗+C∗+I∗+R∗+δI∗,
we have(25)dVdt=S-S∗+C-C∗+I-I∗+R-R∗μS∗+C∗+I∗+R∗+δI∗-μS+C+I+R-δI=S-S∗+C-C∗+I-I∗+R-R∗-μS-S∗-μC-C∗-μI-I∗-μR-R∗-δI-I∗.
Further simplification gives(26)dVdt=-μS-S∗2-μS-S∗C-C∗+I-I∗+R-R∗+δI-I∗-μC-C∗2-μC-C∗S-S∗+I-I∗+R-R∗+δI-I∗-μ+δI-I∗2-μI-I∗C-C∗+R-R∗-μR-R∗2-μR-R∗S-S∗+C-C∗+I-I∗+δI-I∗.
It has therefore been shown thatdV/dt is negative, and additionally at E1 (i.e., if S=S∗, C=C∗, I=I∗, and R=R∗), dV/dt=0. It follows from LaSalle’s invariant principle [51] that all solutions of system (2) approach E1 as t→∞ if RVT>1. Therefore, the endemic equilibrium E1 is globally asymptotically stable in Ω whenever RVT>1 [37, 49]. This completes the proof.


## 2.2.1. Local Stability Analysis ofE0

Theorem 2.
The disease-free equilibrium of system (2) is locally asymptotically stable if RVT<1 and unstable if RVT>1.Proof.
Evaluating the Jacobian matrix of system (2) at the disease-free equilibrium gives(9)JE0=-μ+u1σ-1-u1σβηθ+μu1σ+θ+μ-1-u1σβη1θ+μu1σ+θ+μθ01-u1σβηθ+μu1σ+θ+μ-α+ω+μ+u2γ1-u1σβη1θ+μu1σ+θ+μ00α-δ+μ+u2γ0u1σω+u2γu2γ-θ+μ.
The characteristic polynomial of the Jacobian matrixJ(E0) is given by(10)Pλ=λ4+DTλ3+DT1λ2+DT2λ+DT3,where(11)T=A+B+C-RCVT,T1=A+B-RCVT2μ+u1σ+θ+μ+u1σθ+μ+μθ1-RVT,T2=θ+μ1-RVT+A-RCVT-RCVT+μ+u1σ1-RIVT+δ+μ+u2γθ+μ+A+B-RCVTμθ,T3=θμRIVT-1+μ+u1σθ+μRIVT-1,A=1δ+μ+u2γ,B=1α+ω+μ+u2γ,C=μ+u1σ+θ+μα+ω+μ+u2γδ+μ+u2γ,D=α+ω+μ+u2γδ+μ+u2γ.
The Routh-Hurwitz conditions [46] that guarantee that the eigenvalues of the characteristic polynomial in (10) have negative real parts are given by(12)DT>0,DT1>0,DT2>0,DT3>0,D3TT1T2>D2T22+D3T2T3.
These conditions are easily seen to be satisfied whenRVT<1. Thus, the disease-free equilibrium of system (2) is locally asymptotically stable when RVT<1 and unstable when RVT>1. This completes the proof.


## 2.2.2. Global Stability ofE0

Theorem 3.
The disease-free equilibriumE0 of system (2) is globally asymptotically stable if RVT≤1 and unstable if RVT>1.Proof.
LetV(S,C,I,R), with positive constants, K1 and K2, be a Lyapunov function defined as(13)VS,C,I,R=S-S0-S0ln⁡SS0+K1C+K2I+R-R0-R0ln⁡RR0.
Taking the time derivative of the Lyapunov function we obtain(14)dVdt=1-S0SdSdt+K1dCdt+K2dIdt+1-R0RdRdt,whereK1≥0,K2≥0.
SubstitutingdS/dt,  dC/dt,  dI/dt, and dR/dt in (2) into (14) gives(15)dVdt=1-SS0π-1-u1σβSηC+η1IN-μ+σu1S+θR+K11-u1σβSηC+η1IN-α+ω+μ+u2γC+K2αC-δ+μ+u2γI+1-RR0u2γI+ω+u2γC+σu1S-θ+μR≤K11-u1σβηθ+μu1σ+θ+μC+1-u1σβη1θ+μu1σ+θ+μI-α+ω+μ+u2γC+K2αC-δ+μ+u2γI.
Since(16)S≤S0=πθ+μμu1σ+θ+μ,R≤R0=πu1σμu1σ+θ+μ,N≤πμonΩ,
this implies that(17)dVdt≤1-u1σβηθ+μu1σ+θ+μK1+αK1-α+ω+μ+u2γK2C+1-u1σβη1θ+μu1σ+θ+μK1-δ+μ+u2γK2I.
Equating the coefficient ofI in (17) to zero gives(18)1-u1σβθ+μK1=δ+μ+u2γu1σ+θ+μK2.
ChoosingK1=(δ+μ+u2γ)(u1σ+θ+μ)  and  K2=(1-u1σ)β(θ+μ) and plugging K1 and K2 into (17), we have (19)dVdt≤α+ω+μ+u2γδ+μ+u2γu1σ+θ+μRVT-1C≤0,ifRVT≤1.
AdditionallydV/dt=0 if and only if C=0. Hence, the largest compact invariant set in (S,C,I,R)∈Ω:dV/dt≤0 is the singleton set E0. Therefore, from LaSalle’s invariance principle, we conclude that E0 is globally asymptotically stable in Ω if RVT≤1 [37, 49].


## 2.2.3. Global Stability ofE1

Theorem 4.
The endemic equilibriumE1 of system (2) is globally asymptotically stable whenever RVT>1.Proof.
SupposeRVT>1, and then the existence of the endemic equilibrium point is assured. Using the common quadratic Lyapunov function(20)Vx1,x2,…,xn=∑i=1nci2xi-xi∗2,
as illustrated in [50], we consider the following candidate Lyapunov function:(21)VS,C,I,R=12S-S∗+C-C∗+I-I∗+R-R∗2.
The time derivative ofV(S,C,I,R) in (21) is given by(22)dVdtS,C,I,R=S-S∗+C-C∗+I-I∗+R-R∗dS+C+I+Rdt.
Plugging the equations in system (2) into (22) yields(23)dVdt=S-S∗+C-C∗+I-I∗+R-R∗π-μS+C+I+R-δI.
Now setting(24)π=μS∗+C∗+I∗+R∗+δI∗,
we have(25)dVdt=S-S∗+C-C∗+I-I∗+R-R∗μS∗+C∗+I∗+R∗+δI∗-μS+C+I+R-δI=S-S∗+C-C∗+I-I∗+R-R∗-μS-S∗-μC-C∗-μI-I∗-μR-R∗-δI-I∗.
Further simplification gives(26)dVdt=-μS-S∗2-μS-S∗C-C∗+I-I∗+R-R∗+δI-I∗-μC-C∗2-μC-C∗S-S∗+I-I∗+R-R∗+δI-I∗-μ+δI-I∗2-μI-I∗C-C∗+R-R∗-μR-R∗2-μR-R∗S-S∗+C-C∗+I-I∗+δI-I∗.
It has therefore been shown thatdV/dt is negative, and additionally at E1 (i.e., if S=S∗, C=C∗, I=I∗, and R=R∗), dV/dt=0. It follows from LaSalle’s invariant principle [51] that all solutions of system (2) approach E1 as t→∞ if RVT>1. Therefore, the endemic equilibrium E1 is globally asymptotically stable in Ω whenever RVT>1 [37, 49]. This completes the proof.


## 2.3. Sensitivity Analysis

Sensitivity analysis is used to determine the response of a model to variations in its parameter values. In the present case, the focus is given to determining how changes in the model parameters impact the effective reproduction number. This is done through the normalized forward-sensitivity index. We also use the Latin hypercube sampling and the partial rank correlation coefficients (PRCC) to plot scatter diagrams and Tornado plots to determine the relative importance of the parameters inRVT for the disease transmission and prevalence (see also [52]).Definition 5.
The normalized forward-sensitivity index ofRVT to any parameter, say ρ, as given in [46] can be defined as(27)ΓRVTρ=∂RVT∂ρρRVT.The sensitivity indexes ofRVT with respect to its parameters are computed as follows:(28)ΓRVTβ=∂RVT∂ββRVT=1,ΓRVTη=∂RVT∂ηηRVT=δ+μ+u2γηηδ+μ+u2γ+η1α,ΓRVTη1=∂RVT∂η1η1RVT=αη1αη1+δ+μ+u2γ,ΓRVTω=∂RVT∂ωωRVT=-ωα+ω+μ+u2γ,ΓRVTσ=∂RVT∂σσRVT=-u1σu1σ+θ+μ.Similarly, we can compute the sensitivity indexes ofRVT with respect to the remaining parameters in RVT, in the same manner. Using the parameter values γ=0.4, σ=1, η=0.35, ω=0.2, δ=0.1, β=0.88, α=0.2, θ=0.0839, and μ=0.02, with u1=0.5 and u2=0.5, the sensitivity indexes of RVT are shown in Table 2.Table 2
Sensitivity indexes ofRVT to the parameters in (6).
Parameter Description Sensitivity index β Effective contact rate 1 η Per capita infection rate by Carriers 0.6154 η1 Per capita infection rate by ill individuals 0.3846 θ Loss of immunity 0.6686 α Rate of falling ill 0.0620 μ Natural death rate 0.1031 δ Disease-induced death rate −0.1202 ω Natural recovery rate −0.3226 u1 Vaccination rate −1.8280 σ Effectiveness of vaccination −0.8280 u2 Treatment rate −0.5630 γ Effectiveness of treatment −0.5230The corresponding Tornado plots based on a random sample of1000 points for the twelve parameters in RVT are shown in Figure 3(a). The positive values in Table 2 show a promotion of the propagation of the disease. Therefore an increase in the values of β, η, η1, θ, and α will have an increase in the spread of the disease. For example, ΓRVTβ=1 indicates that increasing the effective contact rate by 10% increases the number of secondary infections by 10%. The negative values in Table 2 indicate a reduction in the effective reproduction number RVT if the values of the corresponding parameters are increased. Thus, a reduction in the values of vaccination u1, treatment u2, and natural recovery ω will lead to an increase in the number of secondary infections in the population.Figure 3
Sensitivity plots.
(a)
The Tornado plots for the elven parameters inRVT (b)
The effect of effective contact rate onRVT (c)
The effect of vaccination onRVT (d)
The effect of treatment onRVTFigure3(a) shows the Tornado plots for the twelve parameters in RVT. It can be seen that, in controlling the spread of bacterial meningitis in a population, more susceptible individuals should be given vaccination. Figure 3(a) also suggests that carriers are likely to have more contacts with the susceptible population than the ill individuals who will typically be bound to their beds during the acute phase of the disease. Therefore, the probability of ill individuals transmitting the infections to susceptibles may be lower than that of carriers who are able to mix well with others within the population. Figures 3(b), 3(c), and 3(d) show the regression plots of effective contact rate (β), vaccination rate (u1), and treatment rate (u2), respectively. Figure 3(b) shows that transmission rate has a positive correlation in the spread of bacterial meningitis. Figure 3(c) shows that vaccination has a negative correlation in the spread of bacterial meningitis and hence vaccination increases the immunity of individuals against the meningitis infection, thereby reducing the spread of the infection. Figure 3(d) shows that individuals who receive treatment after being infected with bacterial meningitis have a higher chance of recovery and that reduces the spread of the infection and deaths due to bacterial meningitis.


## 3. Optimal Control Problem

Since the goal of this paper is to find the best ways to control the spread of meningitis, we define the following optimal control problem:(29)J=minu1σ,u2⁡∫0TA1Ct+A2It+B12u12t+B22u22tdt,subject to model (2).The admissible control setU is Lebesgue measurable, which is defined by(30)U=u1t,u2t∣0≤u1≤u1max≤1,0≤u2≤u2max≤1,t∈0,T.Our objective is to find(u1∗,u2∗)∈U which minimizes the associated cost of the vaccination and the associated cost of the treatment over the specified time interval, as well as minimizing the number of infections at a terminal time (see also [37]). The coefficients A1>0 and A2>0 are constants that are introduced to maintain a balance in the size of C(t) and I(t), respectively. B1>0 and B2>0 are the corresponding weights associated with the cost of vaccination (u1) and treatment (u2), respectively. The higher bounds (maximum) attainable for the control measures u1σ and u2 are u1max and u2max, respectively. We fix the control measures u1 and u2 to lie between 0 and 1 so that u1max=1 and u2max=1. Therefore the attainment of u1max and u2max depends on the number of resources available [37]. These resources may include the human effort, material resources, cost of producing vaccine and disbursement, infrastructural resources, the number of health facilities in the community, and the number of hospital beds at the health facilities. The cost of hospitalization, medical test, diagnosis, drug cost, and so forth (see [38–40]) can be associated with treatment. The cost of vaccination may include the cost of the vaccine, the cost of production, the cost of disbursement, the vaccine storage cost, and other related overheads [37]. The severity of the side effects and overdoses of the vaccination and treatment is taken care of by squaring the control measures, and T is the final time during the optimal simulation.


### 3.1. Existence of the Optimal Control

Model (2) can be written as(31)Gt=KG+FG,where (32)G=StCtItRt,K=-μ+u1σ00θ0-α+ω+μ+u2γ000α-δ+μ+u1σ0u1σω+u2γu2γ-θ+μ,FG=π-1-u1σβStηCt+η1ItNt1-u1σβStηCt+η1ItNt00,andGt is the times derivative of G(t). System (31) is nonlinear with a bounded coefficient.SettingG1=S1(t),C1(t),I1(t),R1(t) and G2=S2(t),C2(t),I2(t),R2(t) gives(33)FG1-FG2=-1-u1σβηS1tC1tN1t-S2tC2tN2t-1-u1σβS1tI1tN1t-S2tI2tN2t1-u1σβηS1tC1tN1t-S2tC2tN2t-1-u1σβS1tI1tN1t-S2tI2tN2t00.Therefore(34)
F
G
1
-
F
G
2
=
-
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
-
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
+
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
-
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
≤
2
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
2
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
≤
2
1
-
u
1
σ
β
η
S
1
C
1
-
S
2
C
2
+
2
1
-
u
1
σ
β
η
1
S
1
I
1
-
S
2
I
2
≤
2
1
-
u
1
σ
β
η
C
1
S
1
-
S
2
+
S
2
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
I
1
S
1
-
S
2
+
S
2
I
1
-
I
2
≤
S
1
-
S
2
2
1
-
u
1
σ
β
η
C
1
+
2
1
-
u
1
σ
β
η
1
I
1
+
2
1
-
u
1
σ
β
η
S
2
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
S
2
I
1
-
I
2
≤
2
1
-
u
1
σ
β
π
μ
η
+
1
S
1
-
S
2
+
2
1
-
u
1
σ
β
η
π
μ
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
π
μ
I
1
-
I
2
≤
M
S
1
-
S
2
+
C
1
-
C
2
+
I
1
-
I
2
,where(35)M=max21-u1σβπμη+1,21-u1σβηπμ,21-u1σβη1πμ,so that(36)DG1-DG2≤VG1-G2,whereV=max⁡M,K<∞.The functionD is therefore uniformly Lipschitz continuous. From the definition of the control measures u1(t) and u2(t) and the constraint on the state variables, such that S(t)>0, C(t)≥0, I(t)≥0, and R(t)≥0, we observe that a solution of system (31) exists [40, 53, 54]. From the objective functional and its associated constraints in model (2), we can find the optimal solution for our model. Firstly, we find the Lagrangian (L) and Hamiltonian (H) for the control problem [55]. The Lagrangian of the optimal problem is given by(37)LS,C,I,R,u1,u2=A1Ct+A2It+B12u12t+B22u22t.Our focus is to find the minimal value of the Lagrangian function, which is done by a pointwise minimization of the Hamiltonian (H) defined as follows (using Pontryagin’s maximum principle):(38)HS,C,I,R,u1,u2,t=A1C+A2I+B12u12+B22u22+λ1π-1-u1σβSηC+η1IN-μ+u1σS+θR+λ21-u1σβSηC+η1IN-α+ω+μ+u2γC+λ3αC-δ+μ+u2γI+λ4u2γI+ω+u2γC+u1σS-θ+μR,where λi,i=1,2,3,4, are the adjoint variables associated with S(t), C(t), I(t), and R(t), defined by(39)dλ1dt=-∂H∂S,dλ2dt=-∂H∂C,dλ3dt=-∂H∂I,dλ4dt=-∂H∂R.Theorem 6.
There exists an optimal control pairu1∗(t),u2∗(t) such that(40)Ju1∗t,u2∗t=minu1,u2∈U⁡Ju1t,u2t.Proof.
We start our proof by considering the properties of the existence of the optimal control (see [56]). Following [57], the set of control measures with corresponding state variables are positive. The set U is convex and closed by definition. Therefore, our optimal system is closed and bounded which ascertains the compactness required for the existence of the optimal control. Additionally, the integrand in the objective functional (29), A1C(t)+A2I(t)+B1/2u12(t)+B2/2u22(t), is convex on the control set U. Furthermore, we can state that there exists a positive constant ρ>1 [58], and nonnegative numbers ν1 and ν2 such that the objective functional has a lower bound of ν1u12+u22ρ/2-ν2 so that(41)Ju1t,u2t≥ν1u12+u22ρ/2-ν2,since the control measures and the state variables are bounded, this leads us to a compact proof of existence of the optimal control.


### 3.2. Characterization of the Optimal Control

We will apply Pontryagin’s maximum principle to the Hamiltonian function above to derive the necessary condition of optimality for our control problem.Theorem 7.
LetS, C, I, and R be optimal state solutions with corresponding optimal control variables u1∗ and u2∗ for the objective functional and its constraints in model (2) with N=S+C+I+R. Then, there exist four adjoint variables λ1, λ2, λ3, and λ4 that satisfy(42)dλ1dt=λ11-u1σβηC+η1IN-1-u1σβSηC+η1IN2+μ+u1σ+λ21-u1σβSηC+η1IN2-1-u1σβηC+η1IN-λ4u1σ,dλ2dt=-A1+λ11-u1σβηSN-1-u1σβSηC+η1IN2+λ2α+ω+μ+u2γ+λ21-u1σβSηC+η1IN2-1-u1σβηSN-λ3α-λ4ω+u2γ,dλ3dt=-A2+λ11-u1σβη1SN-1-u1σβSηC+η1IN2+λ21-u1σβSηC+η1IN2-1-u1σβη1SN+λ3δ+μ+u2γ-λ4u2γ,dλ4dt=-λ1θ+λ4θ+μ,with transversality conditions(43)λiT=0,i=1,2,3,4for  the  control  setU=u1,u2,  such  that∂H∂ui=0,wherei=1,2.Therefore, the optimal control pair (u1∗,u2∗) is given by(44)u1∗=min⁡max⁡0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1max,u2∗=min⁡max⁡0,γλ2-λ4C+γλ3-λ4IB2,u2max.Proof.
We use the Hamiltonian function in (38) in order to obtain the adjoint relations and the transversality conditions. We set the state variables in the Hamiltonian function to S,C,I, and R, and differentiating the Hamiltonian (H) with respect to S,C,I, and R, respectively, yields (42). Also, differentiating the Hamiltonian (H) with respect to the control measures u1σ and u2 in the interior of U, we obtain the optimality conditions below:(45)∂H∂u1=B1u1-λ1βSηC+η1IS+C+I+R-λ1σS+λ2βSηC+η1IS+C+I+R+λ4σS=0,∂H∂u2=B2u2-λ2γC-λ3γI+λ4γC+λ4γI=0.Plugging u1=u1∗ and u2=u2∗ into (45) and solving the optimal control pair (u1∗,u2∗), we have(46)u1∗=1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u2∗=γλ2-λ4C+γλ3-λ4IB2.The two control measures which are bounded with lower bounds zero and upper bounds uimax=1, where i=1,2, give(47)u1∗∈U⟹u1∗t=0if1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS≤0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σSif0<1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS<u1max,u1maxif1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS≥u1max,u2∗∈U⟹u2∗t=0ifγλ2-λ4C+γλ3-λ4IB2≤0,γλ2-λ4C+γλ3-λ4IB2if0<γλ2-λ4C+γλ3-λ4IB2<u2max,u2maxifγλ2-λ4C+γλ3-λ4IB2≥u2max.Using (47), the optimal control measures are characterized as (44), completing the proof.Therefore, our optimality system is given by(48)dSdt=π-1-minmax0,1B1βSη1I+ηCS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN-μ+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσS+θR,dCdt=1-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN-α+ω+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγC,dIdt=αC-δ+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγI,dRdt=minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγI+ω+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγC-θ+μR+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσS,(49)dλ1dt=λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2+λ1μ+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσ+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN-λ4minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσ,(50)dλ2dt=-A1+λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηSN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2+λ2α+ω+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηSN-λ3α-λ4ω+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ,dλ3dt=-A2+λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβη1SN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN2+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβη1SN+λ3δ+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ-λ4minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ,dλ4dt=-λ1θ+λ4θ+μ,withλiT=0,i=1,2,3,4


## 3.1. Existence of the Optimal Control

Model (2) can be written as(31)Gt=KG+FG,where (32)G=StCtItRt,K=-μ+u1σ00θ0-α+ω+μ+u2γ000α-δ+μ+u1σ0u1σω+u2γu2γ-θ+μ,FG=π-1-u1σβStηCt+η1ItNt1-u1σβStηCt+η1ItNt00,andGt is the times derivative of G(t). System (31) is nonlinear with a bounded coefficient.SettingG1=S1(t),C1(t),I1(t),R1(t) and G2=S2(t),C2(t),I2(t),R2(t) gives(33)FG1-FG2=-1-u1σβηS1tC1tN1t-S2tC2tN2t-1-u1σβS1tI1tN1t-S2tI2tN2t1-u1σβηS1tC1tN1t-S2tC2tN2t-1-u1σβS1tI1tN1t-S2tI2tN2t00.Therefore(34)
F
G
1
-
F
G
2
=
-
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
-
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
+
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
-
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
≤
2
1
-
u
1
σ
β
η
S
1
t
C
1
t
N
1
t
-
S
2
t
C
2
t
N
2
t
+
2
1
-
u
1
σ
β
η
1
S
1
t
I
1
t
N
1
t
-
S
2
t
I
2
t
N
2
t
≤
2
1
-
u
1
σ
β
η
S
1
C
1
-
S
2
C
2
+
2
1
-
u
1
σ
β
η
1
S
1
I
1
-
S
2
I
2
≤
2
1
-
u
1
σ
β
η
C
1
S
1
-
S
2
+
S
2
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
I
1
S
1
-
S
2
+
S
2
I
1
-
I
2
≤
S
1
-
S
2
2
1
-
u
1
σ
β
η
C
1
+
2
1
-
u
1
σ
β
η
1
I
1
+
2
1
-
u
1
σ
β
η
S
2
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
S
2
I
1
-
I
2
≤
2
1
-
u
1
σ
β
π
μ
η
+
1
S
1
-
S
2
+
2
1
-
u
1
σ
β
η
π
μ
C
1
-
C
2
+
2
1
-
u
1
σ
β
η
1
π
μ
I
1
-
I
2
≤
M
S
1
-
S
2
+
C
1
-
C
2
+
I
1
-
I
2
,where(35)M=max21-u1σβπμη+1,21-u1σβηπμ,21-u1σβη1πμ,so that(36)DG1-DG2≤VG1-G2,whereV=max⁡M,K<∞.The functionD is therefore uniformly Lipschitz continuous. From the definition of the control measures u1(t) and u2(t) and the constraint on the state variables, such that S(t)>0, C(t)≥0, I(t)≥0, and R(t)≥0, we observe that a solution of system (31) exists [40, 53, 54]. From the objective functional and its associated constraints in model (2), we can find the optimal solution for our model. Firstly, we find the Lagrangian (L) and Hamiltonian (H) for the control problem [55]. The Lagrangian of the optimal problem is given by(37)LS,C,I,R,u1,u2=A1Ct+A2It+B12u12t+B22u22t.Our focus is to find the minimal value of the Lagrangian function, which is done by a pointwise minimization of the Hamiltonian (H) defined as follows (using Pontryagin’s maximum principle):(38)HS,C,I,R,u1,u2,t=A1C+A2I+B12u12+B22u22+λ1π-1-u1σβSηC+η1IN-μ+u1σS+θR+λ21-u1σβSηC+η1IN-α+ω+μ+u2γC+λ3αC-δ+μ+u2γI+λ4u2γI+ω+u2γC+u1σS-θ+μR,where λi,i=1,2,3,4, are the adjoint variables associated with S(t), C(t), I(t), and R(t), defined by(39)dλ1dt=-∂H∂S,dλ2dt=-∂H∂C,dλ3dt=-∂H∂I,dλ4dt=-∂H∂R.Theorem 6.
There exists an optimal control pairu1∗(t),u2∗(t) such that(40)Ju1∗t,u2∗t=minu1,u2∈U⁡Ju1t,u2t.Proof.
We start our proof by considering the properties of the existence of the optimal control (see [56]). Following [57], the set of control measures with corresponding state variables are positive. The set U is convex and closed by definition. Therefore, our optimal system is closed and bounded which ascertains the compactness required for the existence of the optimal control. Additionally, the integrand in the objective functional (29), A1C(t)+A2I(t)+B1/2u12(t)+B2/2u22(t), is convex on the control set U. Furthermore, we can state that there exists a positive constant ρ>1 [58], and nonnegative numbers ν1 and ν2 such that the objective functional has a lower bound of ν1u12+u22ρ/2-ν2 so that(41)Ju1t,u2t≥ν1u12+u22ρ/2-ν2,since the control measures and the state variables are bounded, this leads us to a compact proof of existence of the optimal control.


## 3.2. Characterization of the Optimal Control

We will apply Pontryagin’s maximum principle to the Hamiltonian function above to derive the necessary condition of optimality for our control problem.Theorem 7.
LetS, C, I, and R be optimal state solutions with corresponding optimal control variables u1∗ and u2∗ for the objective functional and its constraints in model (2) with N=S+C+I+R. Then, there exist four adjoint variables λ1, λ2, λ3, and λ4 that satisfy(42)dλ1dt=λ11-u1σβηC+η1IN-1-u1σβSηC+η1IN2+μ+u1σ+λ21-u1σβSηC+η1IN2-1-u1σβηC+η1IN-λ4u1σ,dλ2dt=-A1+λ11-u1σβηSN-1-u1σβSηC+η1IN2+λ2α+ω+μ+u2γ+λ21-u1σβSηC+η1IN2-1-u1σβηSN-λ3α-λ4ω+u2γ,dλ3dt=-A2+λ11-u1σβη1SN-1-u1σβSηC+η1IN2+λ21-u1σβSηC+η1IN2-1-u1σβη1SN+λ3δ+μ+u2γ-λ4u2γ,dλ4dt=-λ1θ+λ4θ+μ,with transversality conditions(43)λiT=0,i=1,2,3,4for  the  control  setU=u1,u2,  such  that∂H∂ui=0,wherei=1,2.Therefore, the optimal control pair (u1∗,u2∗) is given by(44)u1∗=min⁡max⁡0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1max,u2∗=min⁡max⁡0,γλ2-λ4C+γλ3-λ4IB2,u2max.Proof.
We use the Hamiltonian function in (38) in order to obtain the adjoint relations and the transversality conditions. We set the state variables in the Hamiltonian function to S,C,I, and R, and differentiating the Hamiltonian (H) with respect to S,C,I, and R, respectively, yields (42). Also, differentiating the Hamiltonian (H) with respect to the control measures u1σ and u2 in the interior of U, we obtain the optimality conditions below:(45)∂H∂u1=B1u1-λ1βSηC+η1IS+C+I+R-λ1σS+λ2βSηC+η1IS+C+I+R+λ4σS=0,∂H∂u2=B2u2-λ2γC-λ3γI+λ4γC+λ4γI=0.Plugging u1=u1∗ and u2=u2∗ into (45) and solving the optimal control pair (u1∗,u2∗), we have(46)u1∗=1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u2∗=γλ2-λ4C+γλ3-λ4IB2.The two control measures which are bounded with lower bounds zero and upper bounds uimax=1, where i=1,2, give(47)u1∗∈U⟹u1∗t=0if1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS≤0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σSif0<1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS<u1max,u1maxif1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS≥u1max,u2∗∈U⟹u2∗t=0ifγλ2-λ4C+γλ3-λ4IB2≤0,γλ2-λ4C+γλ3-λ4IB2if0<γλ2-λ4C+γλ3-λ4IB2<u2max,u2maxifγλ2-λ4C+γλ3-λ4IB2≥u2max.Using (47), the optimal control measures are characterized as (44), completing the proof.Therefore, our optimality system is given by(48)dSdt=π-1-minmax0,1B1βSη1I+ηCS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN-μ+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσS+θR,dCdt=1-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN-α+ω+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγC,dIdt=αC-δ+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγI,dRdt=minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγI+ω+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγC-θ+μR+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσS,(49)dλ1dt=λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2+λ1μ+minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσ+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN-λ4minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxσ,(50)dλ2dt=-A1+λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηSN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2+λ2α+ω+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβSηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηSN-λ3α-λ4ω+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ,dλ3dt=-A2+λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβη1SN-λ11-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN2+λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβηC+η1IN2-λ21-minmax0,1B1βSηC+η1IS+C+I+Rλ1-λ2+λ1-λ4σS,u1maxβη1SN+λ3δ+μ+minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ-λ4minmax0,γλ2-λ4C+γλ3-λ4IB2,u2maxγ,dλ4dt=-λ1θ+λ4θ+μ,withλiT=0,i=1,2,3,4


## 4. Numerical Simulations of the Optimal Control

The optimality system consisting of the state equations (2) and the adjoint equations (45) is solved using the forward-backward sweep scheme. The solution is started with an initial guess for the control measures and the final time set to T=30 days and later varied to T=60 days and T=80 days. The state system is solved forward in time while the costate system is solved backward in time. The current solutions of the state system together with the initial guess for the control measures u1 and u2 are used to solve the costate system. The controls u1 and u2 are updated using the characterizations in (47). With the current cosystem solution and updated controls, the state system is solved and the whole solution process continues until convergence is achieved. Since bacterial meningitis is an endemic disease especially in Sub-Sahara Africa, parameter values that make effective reproduction number RVT more than unity are considered. Taking the associated costs on carriers and ill individuals as A1=A2=1, B1=2, and B2=4, the graphical results of the model are as follows. Figures 4 and 5, show the dynamics of the infection in the presences of control measures and without control measures in a more localised form. Figure 5(f) shows the effect of varying the cost associated with the treatment control. It is seen that if the weight of treatment cost B2 is continuously increased, the upper bound time reduces to about 9 days as indicated by the red line in Figure 5(f). This is true because a higher cost of treatment causes the use of the treatment control to be less. Figure 6 shows the spreading rate of meningitis with and without control measures in a highly dense population, thus in colleges, prisons, cities, and so forth, which suggests that when carriers become more in a highly populated settings there is a likelihood of the infection spreading faster, thereby leading to a corresponding increase in the number of ill (sick) individuals; hence the trajectories in Figures 6(a) and 6(b) show that there is a significant difference in the number of carriers and infected individuals with and without control measures. Therefore, applying both vaccination and treatment has a higher rate of controlling bacterial meningitis than depending only on one control method. Figures 6(a) and 6(b) also show that as infections increase in the population, there is a need for introducing vaccination and treatment at the earlier stage so as to minimize the spread of the infection in the population. This further suggests that college students and prisoners should be given regular vaccination against meningitis, since bacterial meningitis is more likely to affect college students and prisoners than other people, due to the close proximate of beds in the dormitories and prison cells, especially in Sub-Saharan Africa. Figure 6(c) shows that the trajectory for treatment moves to an upper bound at about t=35 days and slowly decreases to a lower bound, which implies that, during the period before 35 days, substantial amount of the control u2 should be applied while administering the vaccination control in the susceptible population so as to reduce the number of carriers and ill individuals until the 50th day. Afterwards, less amount of treatment can be used since the number of carriers and ill individuals will be considerably reduced by the earlier investment in treatment within the period before 35 days.Figure 4
(a), (b), (c), and (d) show that administering vaccination on susceptible individuals reduces the number of secondary infections in the carrier and the ill population. (e) and (f) show the optimal control profile of vaccination. Initial conditionsS(0)=700,C(0)=250,I(0)=40, and R(0)=10 and the parameter values π=100,η=0.35,ω=0.2,δ=0.1,σ=1,α=0.2andθ=0.0839 were used in this simulation.
(a)
The optimal effect onC(t) with vaccination control (b)
The optimal effect onI(t) with vaccination control (c)
The optimal effect onC(t) with an increase in time (d)
The optimal effect onI(t) with an increase in time (e)
Vaccination control profile (f)
Optimal control profile with an increase in timeFigure 5
(a) and (b) show that giving treatment to carriers and the ill individuals reduces the number of secondary infections in the carrier and the ill population but the disease may not be eradicated if the use of only treatment is stuck to. Furthermore, we increased time so as to see whether the disease could be eradicated when using treatment only. We see in (c) and (d) that there is a possibility of recurrence of the disease using only treatment control strategy as indicated by the curved green line in (c) and the linear flow of the green line in (d). (e) shows that, to minimize bacterial meningitis outbreak in30 days, the treatment control should be held intensively for 30 days at a constant rate. (f) shows the optimal control plot of treatment with different cost weights.
(a)
The optimal effect onC(t) with treatment control (b)
The optimal effect onI(t) with treatment control (c)
The optimal effect onC(t) with an increase in time (d)
The optimal effect onI(t) with an increase in time (e)
Treatment control effect (f)
Optimal control plot of treatment with different cost weightsFigure 6
Further simulation trajectories. (d) shows the immunity levels when the two control measures are applied with the following parameter values:π=100000,η=0.35,ω=0.2,  δ=0.1,σ=1,α=0.2, and θ=0.0839.
(a)
The optimal effect onC(t) with both control measures (b)
The optimal effect onI(t) with both control measures (c)
Simulation diagram for control measures (d)
Immunity of recovered individualsFigures7(a), 7(b), 8(c), and 8(d) show the impact of varying α and η. It can be seen that an increase in α and η without control measures has a corresponding decrease in the number of carriers and a considerable increase in the number of infected individuals. Therefore, applying the control measures gives a drastic decrease in the number of infections in both carriers and the ill compartment at a higher rate of α=0.5 and η=0.65. Furthermore, in Figures 7(a) and 7(b), authors kept the natural recovery rate (ω) fixed at 0.2 while assessing the effect of varying the rate of falling ill on carriers and on ill individuals, respectively. Hence, considering an increased rate from 0.2 to 0.7 of the rate of falling ill as shown in Figures 8(a) and 8(b) demonstrates how serious the infection is in the absence of vaccination or treatment, which indicates that ill individuals have a lower means of recovering from the disease, not even through a higher rate of natural defences of ω=0.5 as indicated in Figure 8(b) without vaccination and treatment, which clearly shows that vaccination and treatment are essential in controlling meningitis irrespective of the age group; therefore we suggest that health authorities should increase the rate at which susceptible individuals get vaccination through public education and also encourage individuals with and without symptoms of the infection to visit health centers for a quick check-up and if the infection is detected an immediate treatment/vaccination should be given on time to avert the spreading of the infection to community members, relatives, or close friends.Figure 7
Simulation results of the optimal control model withα = 0.2-0.5, on carriers and the ill compartment.
(a)
Optimal simulation effect in varyingα on C(t) (b)
Optimal simulation effect in varyingα on I(t)Figure 8
Simulation results of the optimal control model withη = 0.35-0.65. The solutions for the carriers and ill individuals with α=0.2 and α=0.5 give RVT=4.225 and RVT=5.250, respectively. The case where η=0.35 and η=0.65 gives RVT=4.225 and RVT=4.854, respectively, which indicates the existence of an unstable disease-free equilibrium and a stable endemic equilibrium. In (e) and (f), vaccination and treatment rates are set to 0.5 and the effectiveness of the treatment control γ varied. (e) and (f) show that as the effectiveness of the treatment control increases, the number of ill individuals reduces faster as compared to the number of carriers.
(a)
Numerical simulation without control measures onR(t) (b)
Numerical simulation without control measures onR(t) (c)
Optimal simulation effect in varyingη on C(t) (d)
Optimal simulation effect in varyingη on I(t) (e)
Effectiveness of treatmentγ on I(t) (f)
Effectiveness of treatmentγ on C(t)


### 4.1. Vaccination Control Only

Figures4(a)–4(f) show simulation results of implementing only control strategy u1 with u2 set to zero.


### 4.2. Treatment Control Only

Considering treatment control,u2, by setting the vaccination control, u1, to zero, the simulation results are as shown in Figures 5(a)–5(f).


### 4.3. Applying Both Control Measures

Finally, both control measures are implemented using the following initial population size ofS(0)=700000, C(0)=250000, I(0)=40000, and R(0)=10000 and the results are plotted in Figures 6(a)–6(d).Figures9(a), 9(b), 10(a), and 10(b) show the interepidemic nature of some selected parameters in the SCIRS model presented above. Figures 9(a) and 9(b) indicate that as individuals lose immunity and become susceptible again, there is a corresponding increase in the rate of infection, and this suggests that an individual who receives vaccination and/or treatment and remains immune for 3–5 years should immediately go for another vaccination against the infection so as to avert transmission of the infection to that individual. Figures 10(a) and 10(b) depict that the per capita infection rate by ill individuals, η1, and the per capita infection rate by carriers, η, are interlinked and so have a corresponding influence in the spread of the infection. Figure 10(b) also indicates that as the infectivity rate of carriers and ill individuals goes beyond 0.6, there could be a high spread of the infection in the community/society should there be an outbreak, since the effective reproduction number will be getting closer to one (1) with a possibility of going beyond a unit.Figure 9
The interepidemic period of the SCIRS model depending on parametersθ and β.
(a)
The contour plot ofθ and β to RVT (b)
The 3D plot ofRVT to θ and βFigure 10
The interepidemic period of theSCIRS model depending on parameters η and η1.
(a)
The contour plot ofη1 and η to RVT (b)
The 3D plot ofRVT to η1 and η


## 4.1. Vaccination Control Only

Figures4(a)–4(f) show simulation results of implementing only control strategy u1 with u2 set to zero.


## 4.2. Treatment Control Only

Considering treatment control,u2, by setting the vaccination control, u1, to zero, the simulation results are as shown in Figures 5(a)–5(f).


## 4.3. Applying Both Control Measures

Finally, both control measures are implemented using the following initial population size ofS(0)=700000, C(0)=250000, I(0)=40000, and R(0)=10000 and the results are plotted in Figures 6(a)–6(d).Figures9(a), 9(b), 10(a), and 10(b) show the interepidemic nature of some selected parameters in the SCIRS model presented above. Figures 9(a) and 9(b) indicate that as individuals lose immunity and become susceptible again, there is a corresponding increase in the rate of infection, and this suggests that an individual who receives vaccination and/or treatment and remains immune for 3–5 years should immediately go for another vaccination against the infection so as to avert transmission of the infection to that individual. Figures 10(a) and 10(b) depict that the per capita infection rate by ill individuals, η1, and the per capita infection rate by carriers, η, are interlinked and so have a corresponding influence in the spread of the infection. Figure 10(b) also indicates that as the infectivity rate of carriers and ill individuals goes beyond 0.6, there could be a high spread of the infection in the community/society should there be an outbreak, since the effective reproduction number will be getting closer to one (1) with a possibility of going beyond a unit.Figure 9
The interepidemic period of the SCIRS model depending on parametersθ and β.
(a)
The contour plot ofθ and β to RVT (b)
The 3D plot ofRVT to θ and βFigure 10
The interepidemic period of theSCIRS model depending on parameters η and η1.
(a)
The contour plot ofη1 and η to RVT (b)
The 3D plot ofRVT to η1 and η


## 5. Conclusion

We presented a mathematical framework of vaccination and treatment onSCIRS bacterial meningitis model. The investigation of the stability of the model shows that the disease-free equilibrium is locally and globally asymptotically stable and the endemic equilibrium depicts a global stability. Scatter plots and the Tornado plots of the twelve parameters in RVT show that the effective contact rate β has a major impact in transmitting the disease, followed by the infectivity potential of carriers η. This supports the fact that asymptomatic carriers are likely to have more contacts in a community/nation when there is an outbreak of bacterial meningitis compared to ill individuals who will typically be bound to their beds during the acute phase of the disease, thereby lowering the probability of infecting susceptible individuals. Finally, the numerical simulation shows that the optimal (best) way of controlling the transmission of meningitis in Sub-Saharan Africa and the world at large is to encourage susceptible individuals to go for vaccination against meningitis at the health centers and also report any suspected symptoms of meningitis to health practitioners for early detection and immediate care (treatment).



---

*Source: 2657461-2018-03-27.xml*