Stability and Bifurcation Analysis of Measles Model
- O. A. Adeboye
- O. A Odebiyi
- J. K. Oladejo
- E.O. Elijah
- A. A. Taiwo
- O. A. Olajide
- 507-522
- Oct 21, 2024
- Mathematics
Stability and Bifurcation Analysis of Measles Model
O. A. Adeboye 1*, O. A Odebiyi2, J. K. Oladejo3, E.O. Elijah4, A. A. Taiwo5, O. A. Olajide6
1,2,3,5Department of Pure and Applied Mathematics, Ladoke Akintola University of Technology, PMB 4000, Ogbomoso, Nigeria.
4,6Department of Mathematics, Federal University of Technology, P.M.B. 65, Minna Niger state, Nigeria.
DOI: https://doi.org/10.51584/IJRIAS.2024.909045
Received: 08 September 2024; Accepted: 18 September 2024; Published: 21 October 2024
ABSTRACT
A six compartmental deterministic mathematical model, governed by a system of ordinary differential equations for measles was formulated in other to study and analyze the transmission dynamics of measles in human population. The model was shown to be mathematically and epidemiologically meaningful. The basic reproduction number of the model was obtained and global stability of the disease-free equilibrium and endemic equilibrium were obtained and shown to be asymptotically stable, whenever the basic reproduction and unstable if otherwise. More so, if then the endemic equilibrium of the model equation is globally asymptotically stable The effect of some parameters of the model relatives to the basic reproduction number was calculated using the normalized forward sensitivity indices, and it was shown that increase in the parameters with negative indices will reduce the value of the basic Reproduction number, while increase in those with positive indices will increase the value of basic reproduction number. The bifurcation analysis was also carried out and the model was shown to exhibit backward bifurcation which indicates that is no longer sufficient for effective disease control. The numerical result shows that isolation of infective plays a major role in reducing the transmission of the disease in the population.
Keywords: Bifurcation Analysis, Basic reproduction number, sensitivity analysis, stability analysis
INTRODUCTION
Measles is a highly contagious viral illness caused by the measles virus. It’s a serious disease that can lead to serious complications, especially in children. It is highly contagious, serious airborne disease caused by a virus that can lead to severe complications and death. It is an infectious disease and highly contagious respiratory disease through person to person transmission mode, with over 93% transmission rates among susceptible persons. It is the worst eruptive fever during childhood. It also shows characteristics of reddish rash, fever, and leads to serious and fatal complications including, diarrheal, pneumonia, and encephalitis [1]. It can affect anyone but most common among children. The mortality rate and the incidence rate of different infectious disease (e.g. measles, cholera, tuberculosis), play a major public health concern in the developing countries. Most infectious diseases are caused by micro-organisms such as virus, parasites, bacteria, fungi. There modes of transmission is through direct or in direct contact from one-person to another. Transmission spread of these diseases to human population depends on various factors, such as number of susceptible, infective, exposed, modes of transmission (social, ecological, geographical) conditions. These diseases are spread by direct contact between infective and susceptible from droplet of an infected individual by talking, sneezing, coughing, drinking, kissing, or body contact. Diseases such as measles are mostly spread by respiration, while others are spread by vectors or bacteria [2].
Over ten million two hundred thousand (10.2 million) deaths annually attributed to infectious diseases and most of these diseases occur in developing countries [3].
Many infected children suffer blindness, impaired vision or deafness. Measles confers a lifelong immunity from further attacks. Measles vaccination have been very effective, and it’s been prevented by MMR (Measles Mumps Rubella) vaccine. Before the vaccination program, an estimated value of five million to six million (5million-6 million) people are infected annually with 6000-7000 confirmed deaths, also with 52,000 hospitalized who develops a chronic disability from measles encephalitis. However, the global vaccination has helped in reducing the global incidence, but measles remains a public health problem in developing countries. Measles stands as one of the leading vaccine-preventable killer of many children in undeveloped countries like Democratic Republic of Congo (DRC), Madagascar and Nigeria. As at January 2022, 254 cases of measles were confirmed in Nigeria. At the end of 2021, there were over 10,000 such cases in the country [4].There have been 11 outbreaks (defined as more related cases) reported in 2024, and 67% of cases 101 of 151 are outbreak- associated. For comparison, 4 outbreaks were reported during 2023 and 48% 0f cases (28 of 58) were outbreak-associated [5]. In 2018, about 140,000 people died from measles worldwide [4]. The overall case fatality rate for children below 5 years was 12.6%, for unvaccinated children below 5, 16.2% and among children below 9 months, 24%.
Measles can be transmitted from an infected person to another through a contagious respiratory disease or through body contact with an infected person. Measles can be contacted through kissing, hugging, exchange of sweat from an infected person or close respiratory contact with an infected individual. Measles cannot be contacted through handshaking, dishes, door knobs or drinking glasses. An infected person can spread the virus at any stage of the infection. Some people develop the measles symptoms shortly after been infected, but these symptoms quickly show in the children between the ages of 2-8 years. Early detention of the virus can help to reduce complications by using medications and vaccinations. Infact, the World Health Organization African Region established a 2020 measles elimination goal [18].
Some researchers have worked on the model of measles and few among them are:
[1] developed a mathematical analysis of effect of measles. The paper presents a robust compartmental mathematical model of (SVEIR). The model has a disease-free equilibrium which is globally asymptotically stable (GAS). There also exist a unique endemic equilibrium point which is locally stable whenever the association threshold quantity R0 exceed unity. Runge – Kutta of order (4) was used to solve the model numerically. [6] developed a mathematical model for control of measles epidemiology. They used SEIR model to determine the impact of exposed individuals at latent period through the stability analysis and numerical simulation. [7] worked on the dynamical analysis of a model for measles infection. His model determined the required vaccination coverage and dosage that will guarantee eradication of measles within a population. [8] proposed a mathematical model of measles dynamics with vaccination by considering the total number of recovered individuals either from natural recovery or recovery due to vaccination. Numerical simulation of the model shows that vaccination is capable of reducing the number of exposed and infectious population. In the research of [9], a deterministic SIR model was employed to simulate the spread of measles under different vaccination scenarios in a population with a specific size and age distribution. The model accurately forecasted a measles outbreak in 1997, which played a crucial role in informing the decision to launch a comprehensive MMR (Measles, Mumps and Rubella) vaccination campaign in New Zealand that year.
However, this work presents the global stability and bifurcation analysis of measles reoccurrence in vaccinated population the paper is organized as follows. In section 2, we formulate and explain the model positivity. In section 3, we explore existence of disease free equilibrium point, the endemic equilibrium point and the global stabilities of their equilibrium were analyzed using Lyapunov functions, the computation of sensitivity analysis and bifurcation analysis were also investigated. In section 4, the paper ends with some numerical simulations to support and compliment the theoretical finding.
Model Formulation
A mathematical model for measles was formulated in order to study and analyze the transmission dynamics of measles in a human population. The mathematical model was governed by a system of ordinary differential equations which was subdivided into six mutually exclusive classes, namely: Susceptible human \( S_h(t) \) , Exposed human \( E_h(t) \) , Isolated human \( J_h(t) \), Infectious human \( I_h(t) \) , Vaccinated human \( V_h(t) \) , Recovered human \( R_h(t) \)
In this research, a population size \( N(t) \) was partitioned into six subclasses of individuals with sizes denoted by \( S \), \( V \), \( E \), \( I \), \( J \), and \( R(t) \), respectively, such that:
Figure 1: Schematic Diagram of the Model
The following system of ordinary differential equations of the proposed model is therefore considered:
\[
\frac{dS}{dt} = (1 – p)\pi – \frac{\lambda S}{A} – \mu S + \omega V + \sigma R
\]
\[
\frac{dE}{dt} = \frac{\lambda S}{A} – (k + \mu)E
\]
\[
\frac{dI}{dt} = kE – (\epsilon + \tau_1 + \mu + \delta)I \tag{1}
\]
\[
\frac{dJ}{dt} = \epsilon I – (\tau_2 + \mu + \delta)J
\]
\[
\frac{dV}{dt} = P\pi – (\omega + \mu)V
\]
\[
\frac{dR}{dt} = \tau_1 I + \tau_2 J – (\sigma + \mu)R
\]
where:
\[
\lambda = \beta \eta_d I
\]
The model parameters are defined as follows
Table 1: Variables and definitions as used
Variable | Definition |
S (t) | Number of Susceptible at times (t) |
E (t) | Number of Exposed at times (t) |
I (t) | Number of infected at times (t) |
J (t) | Number of Isolated Individuals at times (t) |
V (t) | Number of Vaccinated at times (t) |
R (t) | Number of Recovered at times (t) |
Table 2: Model parameters, definitions, value and source used.
Parameter | Definition | Value | Source |
---|---|---|---|
P | Vaccine Rate | 0.05263 | [1] |
\(\pi\) | Birth Rate | 0.004 | [1] |
\(\beta\) | Contact Rate | 0.2 | [1] |
A | Area per square meter | 1.0 | [1] |
\(\mu\) | Natural death rate | 0.02 | [1] |
\(\omega\) | Vaccine waning rate | 0.1 | [1] |
\(\sigma\) | Loss of immunity | 2.00 | [1] |
\(\tau_1, \tau_2\) | Treatment for infected and isolated individuals | 0.8, 0.6 | [1] |
\(\epsilon\) | Rate of Isolation | 0.6 | [1] |
k | Progression rate to infectious class | 0.3 | [1] |
\(\delta\) | Disease induced death rate | Modification parameter based on Area | 0.09 |
\(\eta_d\) | Disease induced death rate | 0.9 | [1] |
The Invariant Region
Lemma 1: The feasible region of the measles model given by:
\[
D = \left\{ (S,V,E,I,J,R) \in \mathbb{R}_+^6 : S + V + E + I + J + R \leq \frac{\pi}{\mu} \right\} \tag{2}
\]
is positively invariant and attracting.
Proof: Let the total human population of the model be denoted by \( N(t) \). Adding all the parameters of the model (1) together, then, the rate of change of the total human population gives:
\[
\frac{dN}{dt} = \pi – \mu N – \alpha(I + J) \tag{3}
\]
So that,
\[
\frac{dN}{dt} \leq \pi – \mu N \tag{4}
\]
Then, using the method of integrating factor:
\[
\frac{dN}{dt} – \mu N \leq \pi
\]
Solving this gives:
\[
N(t) \leq N(0)e^{-\mu t} + \frac{\pi}{\mu} – \frac{\pi}{\mu} e^{-\mu t}
\]
\[
N(t) \leq N(0)e^{-\mu t} + \frac{\pi}{\mu} \left( 1 – e^{-\mu t} \right) \tag{5}
\]
Since \( N(t) \leq \frac{\pi}{\mu} \) wherever \( N(0) \leq \frac{\pi}{\mu} \), then the region \( D \) is positively invariant. Further, if \( N(0) > \frac{\pi}{\mu} \), either the solution enters \( D \) in finite time or the total population tends to the limit \( \frac{\pi}{\mu} \) and the infected classes tend to zero. Therefore, the feasible region \( D \) is attracting, which implies that all solutions initiated in \( \mathbb{R}_+^6 \) eventually enter \( D \). Therefore, it is sufficient to study the dynamics of measles in the feasible region \( D \), where the model is considered to be mathematically and epidemiologically well-posed.
Positivity and Boundedness of Solution of the Model
Lemma 2: The solution set \( \{ S(t), V(t), E(t), I(t), J(t), R(t) \} \) of the measles model (1) with positive initial data in \( D \), remains positive in \( D \) for all time \( t > 0 \).
Proof: From the first compartment of the model:
\[
\frac{dS}{dt} = (1 – p)\pi – \frac{\lambda S}{A} – \mu S + \omega V + \sigma R
\]
\[
\frac{dS}{dt} + \left( \frac{\lambda}{A} + \mu \right) S(t) \geq 0
\]
Using the integrating factor method gives:
\[
I.F. = e^{\int_0^t \left( \frac{\lambda}{A} + \mu \right) dt}
\]
Simplifying further yields:
\[
S(t) = S_0 \geq e^{\int_0^t \left( \frac{\lambda}{A} + \mu \right) dt} \geq 0
\]
Thence,
\[
S(t) \geq 0 \quad \text{for all} \quad t \geq 0
\]
The remaining variables can be solved following the same procedure to be positive. Therefore, all the state variables are non-negative.
Mathematical Analysis of the Model
Disease-Free Equilibrium
At DFE, it is assumed that there is no infection, i.e., \( E = I = J = 0 \), but the DFE of the model is denoted by \( \varepsilon^0 \) such that at critical points, \( \frac{dS}{dt} = \frac{dV}{dt} = \frac{dE}{dt} = \frac{dI}{dt} = \frac{dJ}{dt} = \frac{dR}{dt} = 0 \). Then, the model equation (1) becomes:
Disease-free Equilibrium at disease-free, \( E = I = J = R = 0 \):
\[
\therefore \varepsilon^0 = (S^0, V^0, E^0, I^0, J^0, R^0) = \left( \frac{1}{\mu} \left[ (1 – p)\pi + \frac{\omega p \pi}{\omega + \mu} \right], \frac{p \pi}{\omega + \mu}, 0, 0, 0, 0 \right) \tag{6}
\]
Endemic Equilibrium Point
At endemic equilibrium, disease is present, but the endemic equilibrium point is denoted by \( \varepsilon^{**} \). At steady states:
\[
\varepsilon^{**} = (S^{**}, V^{**}, E^{**}, I^{**}, J^{**}, R^{**})
\]
Hence, endemic equilibrium for model (1) is obtained as:
\[
V^{**} = \frac{p \pi}{\omega + \mu}, \quad I^{**} = \frac{k E^{**}}{\epsilon + \tau_1 + \mu + \delta}, \quad J^{**} = \frac{k \epsilon E^{**}}{(\epsilon + \tau_1 + \mu + \delta)(\tau_2 + \mu + \delta)}
\]
\[
R^{**} = \frac{E^{**}}{\sigma + \mu} \left( \frac{\tau_1 k}{\epsilon + \tau_1 + \mu + \delta} + \frac{\tau_2 k \epsilon}{(\epsilon + \tau_1 + \mu + \delta)(\tau_2 + \mu + \delta)} \right)
\]
\[
S^{**} = (1 – p)\pi + \frac{\omega p \pi}{\omega + \mu} + \frac{\sigma + \omega}{\sigma + \omega} \left( \frac{\tau_1 k}{\epsilon + \tau_1 + \mu + \delta} + \frac{\tau_2 k \epsilon}{(\epsilon + \tau_1 + \mu + \delta)(\tau_2 + \mu + \delta)} \right)
\]
\[
E^{**} = \frac{\mu \sigma A (\epsilon + \tau_1 + \mu + \delta)^2 (\tau_2 + \mu + \delta)(\sigma + \omega)(k + \mu)(R_0 – 1)}{\beta \eta_d k \sigma \left[ (\sigma + \omega)(k + \omega)(\epsilon + \tau_1 + \mu + \delta)(\tau_2 + \mu + \delta) – \sigma k \left( \tau_1 (\tau_2 + \mu + \delta) + \tau_2 \epsilon \right) \right]} \tag{7}
\]
The endemic exists only if \( R_0 > 1 \). If \( R_0 < 1 \), no endemic equilibrium exists, and it also coincides with the disease-free equilibrium.
Basic Reproduction Number
The basic reproduction number \( R_0 \) is a fundamental concept in epidemiology that represents the average number of secondary cases generated by a single infectious person in a completely susceptible population during the early stages of an outbreak [10].
Considering the infection-related compartment:
\[
\frac{dE}{dt} = \frac{\lambda S}{A} – (K + \mu)E
\]
\[
\frac{dI}{dt} = KE – (\epsilon + \tau_1 + \mu + \delta)I
\]
\[
\frac{dJ}{dt} = \epsilon I – (\tau_2 + \mu + \delta)J
\]
\[
\frac{dR}{dt} = \tau_1 I + \tau_2 J – (\sigma + \mu)R
\]
Where:
\( K_1 = K + \mu \)
\( K_2 = \epsilon + \tau_1 + \mu + \delta \)
\( K_3 = \tau_2 + \mu + \delta \)
\( K_4 = \sigma + \mu \)
\[
F = \begin{pmatrix} \frac{\beta \eta_d I S_0}{A} \\ 0 \\ 0 \end{pmatrix}, \quad V = \begin{pmatrix} K_1 E \\ K_2 I – KE \\ K_3 J – \epsilon I \end{pmatrix} \tag{8}
\]
\[
F = \begin{pmatrix} 0 & \frac{\beta \eta_d S_0}{A} & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix}, \quad V = \begin{pmatrix} K_1 & 0 & 0 \\ -K_1 & K_2 & 0 \\ 0 & -\epsilon & K_3 \end{pmatrix} \tag{9}
\]
\[
FV^{-1} = \begin{pmatrix} 0 & \frac{\beta \eta_d S_0}{A} & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix} \begin{pmatrix} \frac{1}{K_1} & 0 & 0 \\ \frac{K}{K_1 K_2} & \frac{1}{K_2} & 0 \\ \frac{k \epsilon}{K_1 K_2 K_3} & \frac{\epsilon}{K_2 K_3} & \frac{1}{K_3} \end{pmatrix}
\]
\[
= \begin{pmatrix} \frac{\beta \eta_d K S_0}{A K_1 K_2} & \frac{\beta \eta_d S_0}{A K_2} & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix}
\]
\[
|FV^{-1} – \lambda I| = 0, \quad \text{where} \quad I \text{ is a 3×3 identity matrix}
\]
\[
\left| \begin{pmatrix} \frac{\beta \eta_d K S_0}{A K_1 K_2} – \lambda & \frac{\beta \eta_d S_0}{A K_2} & 0 \\ 0 & -\lambda & 0 \\ 0 & 0 & -\lambda \end{pmatrix} \right| = 0 \tag{10}
\]
Thus,
\[
\left( \frac{\beta \eta_d K S_0}{A K_1 K_2} – \lambda \right) \lambda^2 = 0
\]
\[
\lambda = 0 \quad \text{(twice)} \quad \text{or} \quad \lambda = \frac{\beta \eta_d K S_0}{A K_1 K_2} = R_0
\]
The basic reproduction number \( R_0 = \rho(FV^{-1}) \), is the spectral radius of the dominant eigenvalue of (10).
Hence,
\[
R_0 = \frac{\beta \eta_d K}{\mu A (K + \mu)(\epsilon + \tau_1 + \mu + \delta)} \left[ (1 – p)\pi + \frac{\omega p \pi}{\omega + \mu} \right] \tag{11}
\]
Stability Analysis
Global Stability of the Disease-Free Equilibrium
Theorem 2: If \( R_0 \leq 1 \), then the disease-free equilibrium is globally asymptotically stable, and unstable otherwise.
Proof: Let \( F \) be a candidate Lyapunov function for the model, such that:
\[
F = (S – S^* – S^* \ln \frac{S}{S^*}) + \frac{K E}{(K + \mu)(\epsilon + \tau_1 + \mu + \delta)} + \frac{I}{\epsilon + \tau_1 + \mu + \delta} \tag{12}
\]
Where \( S^* = \frac{1}{\mu} \left[ (1 – P)\pi + \frac{\omega P \pi}{\omega + \mu} \right] \) is the value of \( S(t) \) at DFE.
Obviously, the second and third terms on the RHS of (12) are positive. For the first term, \( S^* \leq S \) (since \( S^* \) is an equilibrium point of \( S \)). Then \( S – S^* – S^* \ln \frac{S}{S} \) is also positive.
Therefore, \( F = F(S,E,I) \) is positive definite. Now for the time derivative of \( F \) along the solutions of the model equation (1), we have:
\[
\frac{dF}{dt} = \left( 1 – \frac{S^*}{S} \right) \frac{dS}{dt} + \frac{K}{(K + \mu)(\epsilon + \tau_1 + \mu + \delta)} \frac{dE}{dt} + \frac{1}{\epsilon + \tau_1 + \mu + \delta} \frac{dI}{dt} \tag{13}
\]
Substituting \( \frac{dS}{dt}, \frac{dE}{dt}, \) and \( \frac{dI}{dt} \) from (1) gives:
\[
\frac{dF}{dt} = (1 – \frac{S^*}{S}) \left[ (1 – P)\pi – \frac{\beta \eta_d IS}{A} – \mu S + \omega V + \sigma R \right]
+ \frac{K}{K + \mu}(\epsilon + \tau_1 + \mu + \delta) \left[ \frac{\beta \eta_d IS}{A} – (K + \mu)E \right]
+ \frac{1}{\epsilon + \tau_1 + \mu + \delta} \left[ KE – (\epsilon + \tau_1 + \mu + \delta)I \right]
\]
\[
\frac{dF}{dt} = (1 – \frac{S^*}{S}) \left[ \frac{\beta \eta_d I^* S^*}{A} – \frac{\beta \eta_d IS}{A} + \mu (S^* – S) + \omega (V – V^*) + \sigma (R – R^*) \right]
+ \frac{K}{K + \mu} (\epsilon + \tau_1 + \mu + \delta) \left[ \frac{\beta \eta_d IS}{A} – (K + \mu) \cdot \frac{\beta \eta_d I^* S^*}{A (K + \mu)} \right] + \frac{1}{\epsilon + \tau_1 + \mu + \delta} \left[ \epsilon + \tau_1 + \mu + \delta (I^* I) \right]
\]
At DFE, \( I^* = R^* = 0 \), so \( \frac{dF}{dt} \) becomes:
\[
\frac{dF}{dt} = -\frac{\beta \eta_d IS}{A} \left( 1 – \frac{S^*}{S} \right) – \mu (S – S^*) – \omega (V – V^*) + \sigma R
+ \left[ \frac{K \beta \eta_d S}{A (K + \mu)(\epsilon + \tau_1 + \mu + \delta)} – 1 \right] I
\]
At disease-free equilibrium, \( \varepsilon^0 = (S^0, V^0, E^{0*}, I^0, J^0, R^0) = \left( \frac{1}{\mu} \left[ (1 – P)\pi + \frac{\omega P \pi}{\omega + \mu} \right], \frac{P \pi}{\omega + \mu}, 0, 0, 0, 0 \right) \), (14).
Therefore,
\[
\frac{dF}{dt} = -\frac{\beta \eta_d IS}{A} \left( \frac{S – S^*}{S} \right) – \mu (S – S^*) – \omega (V – V^*) + \left\{ \frac{\beta \eta_d k}{\mu A (K + \mu)(\epsilon + \tau_1 + \mu + \delta)} \left[ (1 – P)\pi + \frac{\omega P \pi}{\omega + \mu} \right] – 1 \right\} I
\]
\[
\frac{dF}{dt} = -\frac{\beta \eta_d IS}{A} \left( \frac{S – S^*}{S} \right) – \mu (S – S^*) – \omega (V – V^*) + (R_0 – 1) I \tag{15}
\]
Obviously, from (15), \( \frac{dF}{dt} < 0 \) if \( R_0 \leq 1 \).
\( \frac{dF}{dt} = 0 \) if and only if \( S = S^*, V = V^*, I = 0 \).
Thus, \( (S, V, E, I, J, R) \to \left( \frac{1}{\mu} \left[ (1 – P)\pi + \frac{\omega P \pi}{\omega + \mu} \right], \frac{P \pi}{\omega + \mu}, 0, 0, 0, 0 \right) \) as \( t \to \infty \), and the target compact invariant set is the singleton \( \{ \varepsilon^* \} \). So, by LaSalle’s invariance principle [16], every solution of the model system (1) with initial conditions in \( D \) approaches \( \varepsilon^* \) as \( t \to \infty \) whenever \( R_0 \leq 1 \). Hence, the disease-free equilibrium is globally asymptotically stable whenever \( R_0 \leq 1 \), and unstable otherwise.
Global Stability of Endemic Equilibrium
Theorem 4: If \( R_0 > 1 \), then the endemic equilibrium point of the model equation is globally asymptotically stable in \( \Omega \), provided \( S \geq S^*, V \geq V^*, E \geq E^*, I \geq I^*, J \geq J^*, R \geq R^* \).
Proof: To establish the global stability of the endemic equilibrium \( E^* \), we analyze by constructing the following quadratic Lyapunov function \( L \), following the approach of [14], such that:
\[
L = \frac{1}{2} \left[ (S – S^*) + (V – V^*) + (E – E^*) + (I – I^*) + (J – J^*) + (R – R^*) \right]^2
\]
By direct calculation of the time derivatives \( L(t) \) along the solutions of the system (1), we obtain:
\[
\frac{dL}{dt} = \left[ (S – S^*) + (V – V^*) + (E – E^*) + (I – I^*) + (J – J^*) + (R – R^*) \right] \left( \frac{dS}{dt} + \frac{dV}{dt} + \frac{dE}{dt} + \frac{dI}{dt} + \frac{dJ}{dt} + \frac{dR}{dt} \right)
\]
Substituting the appropriate solutions of the system (1) into the derivative of \( L(t) \) gives:
\[
\frac{dL}{dt} \leq \left[ (S – S^*) + (V – V^*) + (E – E^*) + (I – I^*) + (J – J^*) + (R – R^*) \right] \frac{dN}{dt}
\]
\[
\frac{dL}{dt} \leq \left[ (S – S^*) + (V – V^*) + (E – E^*) + (I – I^*) + (J – J^*) + (R – R^*) \right] (\Lambda – \mu N)
\]
\[
\leq (N – \frac{\pi}{\mu})(\Lambda – \mu N) \tag{17}
\]
We obtain the result by rearranging and simplifying (17):
\[
\leq -\frac{1}{\mu} (\pi – \mu N)^2
\]
Let \( \xi = \pi – \mu N \). Then,
\[
\frac{dL}{dt} \leq -\frac{1}{\mu} \xi^2
\]
Hence, \( \frac{dL}{dt} \leq 0 \) and \( \frac{dL}{dt} = 0 \) if and only if \( S = S^*, V = V^*, E = E^*, I = I^*, J = J^*, R = R^* \).
Therefore, if \( X < Y \), then \( \frac{dL}{dt} \) will be negative definite, implying that \( \frac{dL}{dt} = 0 \) if and only if \( S = S^*, V = V^*, E = E^*, I = I^*, J = J^*, R = R^* \). Therefore, the largest positive invariant set in \( \{ (S^*, V^*, E^*, I^*, J^*, R^*) \in \Omega : \frac{dL}{dt} = 0 \} \) is a singleton \( \{ E_1 \} \), where \( E_1 \) is globally asymptotically stable in the set \( \Omega \) in accordance with LaSalle’s invariance principle [16]. It then implies that \( E_1 \) is globally asymptotically stable in \( \Omega \).
Sensitivity Analysis
Using the approach of [15], the normalized forward sensitivity index of a variable “P” that depends differentiably on a parameter “q” is defined as:
\[
X_v^{R_0} = \frac{\partial R_0}{\partial v} \cdot \frac{v}{R_0} \tag{18}
\]
As we have an explicit formula for \( R_0 \) in equation (18), we derive an analytical expression and the associated numeric values for the sensitivity of \( R_0 \) as \( X_v^{R_0} = \frac{\partial R_0}{\partial v} \cdot \frac{v}{R_0} \) with respect to each of the parameters involved in \( R_0 \) in Table 2, which is depicted by the bar chart in Figure 2.
Table 3: Parameter Values, Sensitivity Expression, and Sensitivity Indices of \( R_0 \)
Parameter | Sensitivity Expression | Sensitivity Index Value |
---|---|---|
\( \beta \) | 1 | 1 |
\( A \) | -1 | -1 |
\( \pi \) | 1 | 1 |
\( \delta \) | \( \frac{\delta}{k_2} \) | -0.07438 |
\( \tau \) | \( \frac{-\tau}{k_2} \) | -0.661157 |
\( \eta_d \) | 1 | 1 |
\( \mu \) | \( \frac{-\mu(k_1 + k_2 + k_1 k_2 (1 + \omega \rho \pi))}{\mu k_1 k_2 (\omega + \mu)^2} \) | -0.92288 |
\( P \) | \( \frac{-p \mu}{\mu(1 – p) + \omega} \) | 0.05263 |
\( K \) | \( \frac{\mu}{k + \mu} \) | 0.06250 |
\( \omega \) | \( \frac{\omega \rho \mu}{(\omega + \mu)(\mu(1 – p) + \omega)} \) | 0.04386 |
\( \epsilon \) | \( \frac{-\epsilon}{k_2} \) | -0.24793 |
Figure 2: Sensitivity indices value chart
Interpretation of Sensitivity Indices
Table 3, Figure 2, and Figure 8 represent the sensitivity index for the baseline parameter values. They show that when the parameters with positive sensitivity index values increase while the other parameters remain constant, the value of \( R_0 \) will also increase, implying that they increase the endemicity of the disease as they have positive indices and should be targeted by intervention strategies.
Moreover, when the parameters with negative sensitivity index values increase while keeping other parameters constant, the value of \( R_0 \) will decrease, implying that they decrease the endemicity of the disease as they have negative indices as depicted by plots in Figures 10, 11, 12, and 13 and should equally be targeted by intervention strategies. For instance, \( X_\beta^{R_0} = +1.0 \) means that increasing or decreasing \( \beta_1 \) by 10% increases (or decreases) \( R_0 \) by 10%, while \( X_{\tau_1}^{R_0} = – 0.661157 \) means that increasing or decreasing \( \tau_1 \) by 10% decreases (or increases) \( R_0 \) by 6.61157%. We can easily calculate for every other parameter following a similar procedure.
Bifurcation Analysis
Existence of Backward Bifurcation
Bifurcation refers to a change in the stability or qualitative behavior of a system. The backward bifurcation refers to a sudden change in the system’s behavior, where the disease can persist and spread even if the basic reproduction number is less than 1, which typically indicates disease extinction.
The existence of backward bifurcation is explored using the center manifold made popular by [11]. This has been widely used in the study of some epidemiological models [12-13]. For understanding and convenience, we consider the following change of variables:
Let \( S = x_1 \), \( E = x_2 \), \( I = x_3 \), \( J = x_4 \), \( V = x_5 \), and \( R = x_6 \). Also, let \( x = (x_1, x_2, x_3, x_4, x_5, x_6)^T \), so the system (1) can be rewritten in the form \( \frac{dx}{dt} = F(x) \), where \( F = (f_1, f_2, f_3, f_4, f_5, f_6)^T \), as follows:
\[
\frac{dx_1}{dt} = (1 – p)\pi – \frac{\lambda x_1}{A} – \mu x_1 + \omega x_5 + \sigma x_6 = 0
\]
\[
\frac{dx_2}{dt} = \frac{\lambda x_1}{A} – (k + \mu) x_2 = 0
\]
\[
\frac{dx_3}{dt} = kx_2 – (\epsilon + \tau_1 + \mu + \delta) x_3 = 0 \tag{19}
\]
\[
\frac{dx_4}{dt} = \epsilon I – (\tau_2 + \mu + \delta) x_4
\]
\[
\frac{dx_5}{dt} = P\pi – \omega x_5 – (\omega + \mu) x_5 = 0
\]
\[
\frac{dx_6}{dt} = \tau_1 x_3 + \tau_2 x_4 – (\sigma + \mu) x_6 = 0
\]
Taking \( \beta_1 = \beta_1^* \), where \( \beta_1^* \) is the chosen bifurcation parameter and the case when \( R_0 = 1 \), then:
\[
\beta_1^* = \frac{\mu A (K + \mu)(\epsilon + \tau_1 + \mu + \delta)(\omega + \mu)}{\eta_d K \pi (\omega + \mu(1 – \rho))}
\]
The Jacobian of system (1), evaluated at the disease-free equilibrium \( E_0 \), and \( \beta_1^* \) is given by:
\[
J(E_0, \beta_1^*) = \begin{pmatrix} -\mu & 0 & -\frac{\beta \eta_d S_0}{A} & 0 & \omega & \sigma \\ 0 & -(K + \mu) & \frac{\beta \eta_d S_0}{A} & 0 & 0 & 0 \\ 0 & K & -(\epsilon + \tau_1 + \mu + \delta) & 0 & 0 & 0 \\ 0 & 0 & \epsilon & -(\tau_2 + \mu + \delta) & 0 & 0 \\ 0 & 0 & 0 & 0 & -(\omega + \mu) & 0 \\ 0 & 0 & \tau_1 & \tau_2 & 0 & -(\sigma + \mu) \end{pmatrix} \tag{20}
\]
The associated right eigenvectors \( w \) corresponding to the zero eigenvalue, where \( w = (w_1, w_2, w_3, w_4, w_5, w_6)^T \), can be obtained from \( J(E_0, \beta_1^*) w = 0 \). Hence, we have:
\[
\begin{aligned}
&w_1 = \frac{w_3 \left( \beta^* \eta_d S_0 K_3 K_4 + \sigma (K_3 \tau_1 + \epsilon \tau_2) \right)}{A K_3 K_4} \\
&w_2 = \frac{w_3 K_2}{K} \\
&w_3 = w_3 > 0 \\
&w_4 = \frac{\epsilon w_3}{K_3} \\
&w_5 = 0 \\
&w_6 = \frac{w_3 (K_3 \tau_1 + \epsilon \tau_2)}{K_4 K_3}
\end{aligned} \tag{21}
\]
Similarly, the Jacobian, \( J(E_0, \beta_1^*) \), has a left eigenvector \( v = (v_1, v_2, v_3, v_4, v_5) \) given by:
\[
\begin{aligned}
&v_1 = 0 \\
&v_2 = v_2 > 0 \\
&v_3 = \frac{\beta^* \eta_d S_0 v_2}{A K_3} \\
&v_4 = 0 \\
&v_5 = 0 \\
&v_6 = 0
\end{aligned} \tag{22}
\]
Satisfying \( v \cdot w = 1 \), the product and substitution give:
\[
\frac{v_2 w_3 K_2}{K} = 1 \quad \Rightarrow \quad v_2 w_3 = \frac{K}{K_2}
\]
Thus, the preceding equality is satisfied if we choose the values of:
\[
v_2 = \frac{1}{K_2}, \quad w_2 = K \tag{23}
\]
Computation of \( a \) and \( b \)
The coefficients \( a \) and \( b \), as defined by Castillo-Chavez and Song [17], are given by:
\[
a = \sum_{i,j,k=1}^{\infty} v_k w_i w_j \frac{\partial^2 f_k (E_0, \beta^*)}{\partial x_i \partial x_j}, \quad b = \sum_{k,i,j=1}^{\infty} v_k w_i \frac{\partial^2 f_k (E_0, \beta^*)}{\partial x_i \partial \beta^*}
\]
They are algebraically computed as follows, considering only the nonzero components \( f_i \) for \( i, j, k = 1, 2, 3, 4, 5, 6 \). With these definitions, the following were obtained:
\[
a = \frac{2 \beta^* \eta_d K^2 \left( \beta^* \eta_d K_3 K_4 S_0 + A \sigma (K_3 \tau_3 + \epsilon \tau_2) \beta^* \right) (\omega + \mu (1 – p))}{A \mu K_2 K_3 K_4 (\omega + \mu)} \tag{24}
\]
\[
b = \frac{\eta_d K^2 (\omega + \mu (1 – p))}{A \mu K_2 (\omega + \mu)} \tag{25}
\]
Clearly, from (24) and (25), \( a > 0 \) and \( b > 0 \), then the model exhibits a backward bifurcation around \( R_0 = 1 \), that is, there exists bi-stability of equilibrium points (one stable disease-free and stable endemic equilibrium for \( R_0 < 1 \)). This implies that the disease-free state is not globally stable. From an epidemiological point of view, the existence of backward bifurcation implies that the classical necessary requirement of \( R_0 < 1 \) is insufficient to control the spread of measles in the population.
Numerical Simulations and Discussion
This section presents the numerical simulation results of modeling measles reoccurrence in vaccinated infants. The numerical simulations of the models are done using parameter values obtained from the literature as shown in Table (2), and this is done via MAPLE 18 software with initial conditions \( S(0) = 2885 \), \( V(0) = 4192 \), \( E(0) = 15784 \), \( I(0) = 1645 \), \( J(0) = 1050 \), and \( R(0) = 14215 \).
Figure 3. Behavior of total population at Disease Free Equilibrium
Figure 4. Behavior of total population at Endemic Equilibrium
The numerical results of total population at disease free equilibrium point for measles virus model is showed in Figure 2. Similarly, result of total population at Endemic Equilibrium is also showed in Figure 4. Additionally, Figure 5 shows that the disease-free equilibrium exist and it is globally stable as shown by the plot of I(t) at different initial values.
Figure 5. Global Stability of Disease-Free Equilibrium
The portrait of global stability of the disease free equilibrium with various initial conditions as illustrated in Figure 5. It simply suggest that if the basic reproduction number , the measles virus can be eradicated from the population, regardless of how many people are initially infected and the disease will eventually die out, regardless of the initial number of infective individuals whenever.
Figure 6. Backward Bifurcation Diagram of Basic Reproduction number
Backward bifurcation is a mathematical modeling phenomenon where a system’s dynamics shift, resulting in the coexistence of multiple stable states as indicated in Figure 6. In the context of measles, this means that a critical threshold (Ro) exists, often below 1, beyond which the disease exhibits complex behavior. Below this threshold, isolation and vaccination efforts are likely to be effective, but beyond it, controlling the disease becomes increasingly difficult.
Figure 7: The Vaccinated Human against time
Figure 7 shows the graph of Vaccinated Human (I) against time (t) at varying vaccine waning rate, It shows the effect of waning of vaccine on the population of the Infectious class. It was noted that at when =0.6, there is an outright decrease in Vaccinated Human population, the figure shows that there is the possibility of high degree of infections as individuals lose the protection that the vaccine offers recorded within the first five to seven years. The figure shows that there is gradual decrease in the number of Vaccinated human population rate at =0.1 and =0.3 respectively. This indicated that the rate of Vaccinated human vanishes at the same point at time (t)= 18,10 years respectively.
Figure 8(a): Graph of \( \beta \) versus Area per Square Meter \( A \) ,Figure 8 (b): Graph of \( \beta \) versus \( \tau_1 \)
It is observed that as the Area per square meter (A) increases, the number of transmission is reduced. This therefore play an important role in the spread and transmission of measles virus.
CONCLUSION
This study was formulated and analyzed as a mathematical model for measles in order to gain more insights and understanding of the epidemiological features of some parameters on the transmission dynamics of measles in human population. The mathematical model was shown to be mathematically and epidemiologically well-posed through the theory of positivity and boundedness of solutions in D. The disease-free and endemic equilibrium points of the model were obtained. The basic reproduction number of the model was calculated using the next generation matrix method and the stability of the disease-free equilibrium was investigated and shown to be locally asymptotically stable whenever the basic reproduction number is less than unity and unstable if otherwise. The global asymptotic stability of the model was shown to be globally asymptotically stable whenever the associated threshold parameter is less than unity and unstable if otherwise.
The effect of some parameters of the model relative to the basic reproduction number was calculated using the normalized forward sensitivity indices and it was established that increase in the parameters with negative indices will reduce the value of the basic reproduction number while increase in those with positive indices will increase the value of the basic reproduction number. The model was shown to exhibits backward bifurcation which implies that the existence of the classical necessary requirement of is insufficient to control the spread of measles in the population. Therefore, epidemiological features such as the vaccination rate, effective contact rate, progression rate and treatment rate should be given great attention in order to effectively control and prevent the dynamical spread of measles infections in human population.
The results of the bifurcation analysis suggest that the model exhibit a backward bifurcation, which has significant implications for measles control. This occurs when a system exhibits a sudden change in behavior, resulting in a stable equilibrium coexisting with an unstable equilibrium, allowing it to maintain a foothold in the population even with low values.
Conflict of interest: The authors declares no conflict of interest
Availability of Data: All data used are declared in the work
REFERENCE
- Adewale, S.O. Mohammed I.T. and Olopade. I.A. (2014). Mathematical analysis of the dynamical spread of measles, IOSR J. Engineering. 4(3): 43–57.
- Shukla V, Maan HS, Dhole TN. Identification of different lineages of measles virus strains circulating in Uttar Pradesh, North India. Virology Journal.2012; 9: 237.
- Walker, C.L.F., Munos, M.K., Black R.E. (2013) Quantifying the indirect effects of key child survival interventions for pneumonia, diarrhoea, and measles. Epidemiology and Infection 141 (1), 115-131.
- W.H.O (2024). Measles cases in Nigeria 2018-2022, by status.
- Centers for disease control and prevention (2024). Measles cases and outbreaks. https://www.cdc.gov>data-res….
- Momoh, AA, Ibrahim MO, Uwanta I.J., and Manga S.B., (2013): Mathematical model for control of measles epidemiology, Intentional Journal Pure Appl. Math. 87(5) (2013), pp. 707-717.
- Bolarin, G. (2014) On the dynamical analysis of a new model for measles infection, Intentional Journal Mathematics Trends Technology. 7(2):144–154
- Robert M.G. and Tobias M.I. (1999), Predicting and preventing measles epidemics in New Zealand: application of a mathematical model. Ag Research, WallaceŠille Animal Research Centre, Upper Hutt, New Zealand.
- Roberts, M.G. (2000) – The elimination of childhood diseases by mass vaccination. Journal of Applied Mathematics and Decision Sciences, 2000 – dml.mathdoc.fr
- O. Diekmann, J.A. Heesterbeek, and J. A.J. Metz, (1990). On the definition and the computation of the basic reproduction ratio in models for infectious diseases in heterogeneous population. J. maths. Bio.28:265-382
- Castillo-Chavez and B. Song. (2024) Dynamical Models of Tuberculosis and their Applications, Math. Biosci. Eng. 1(2) (2004), pp.361-404.
- Olaniyi, S., K.O. Okosun, S.O. Adesanya and R.S. Lebelo. (2020) Modelling malaria dynamics with partial immunity and protected travelers: Optimal control and cost effectiveness analysis. Journal of Biological Dynamics, vol. 14, No. 1, 90-115. http://doi.org/10.1080/17513758.2020.1722265
- Adeyemi, M.O., and T.O. Oluyo; (2023). Mathematical Modeling for the Control of Fly Borne Mastitis Disease in Cattle. Frontiers in applied mathematics and statistics. June, 2023 DOI: 10.3389/fams.2023.1171157
- Odebiyi, OA; Oladejo, JK; Yahaya, AA and Elijah, EO (2024). Analysis of HIV/AIDS Model with nonlinear Incidence Function,. International Journal of Research and Innovation in Applied Science. 9 (4): 317-341. DOI: 10.511584/IJRIAS/
- Chitnis N, Cushing J, M and Hyman M (2008). Determining Important parameters in the spread of Malaria: Through the sensitivity analysis of Mathematical model. Bulleting of Mathematical Biology, 70:1272-1296
- LaSalle, J.P. the stability of dynamical systems, regional conference series in Applied Mathematics. SIAM, Philadelphia,1976.
- Castillo-Chavez and B. Song. (2024) Dynamical Models of Tuberculosis and their Applications, Math. Biosci. Eng. 1(2) (2004), pp.361-404.
- Center for Disease Control (CDC) (2023). Morbidity and Mortality weekly Report (MMWR). Progress Toward Measles Elimination- African Region, 2017-2021. Weekly/September 8, 2023/72(36):985-991