Abstract
Tuberculosis (TB) and coronavirus (COVID-19) are both infectious diseases that globally continue affecting millions of people every year. They have similar symptoms such as cough, fever, and difficulty breathing but differ in incubation periods. This paper introduces a mathematical model for the transmission dynamics of TB and COVID-19 coinfection using a system of nonlinear ordinary differential equations. The well-posedness of the proposed coinfection model is then analytically studied by showing properties such as the existence, boundedness, and positivity of the solutions. The stability analysis of the equilibrium points of submodels is also discussed separately after computing the basic reproduction numbers. In each case, the disease-free equilibrium points of the submodels are proved to be both locally and globally stable if the reproduction numbers are less than unity. Besides, the coinfection disease-free equilibrium point is proved to be conditionally stable. The sensitivity and bifurcation analysis are also studied. Different simulation cases were performed to supplement the analytical results.
1. Introduction
The outbreak of coronavirus (COVID-19) was first reported by the World Health Organization (WHO) as pneumonia of unknown cause on December 31, 2019, observed in Wuhan, China. However, on the 10th of January 2020, the virus was determined to be the new family of novel coronavirus with the same category as the severe acute respiratory syndrome virus (SARS-CoV) and the Middle Eastern Respiratory Syndrome virus (MERS-CoV) that occurred, respectively, in Asia and Saudi Arabia [1]. COVID-19 is a respiratory virus that transmits mainly through droplets of saliva generated from an infected person through coughing or sneezing. The rapid expansion of the virus forced the WHO to declare it as a pandemic on March 11, 2020 [2]. The WHO report of COVID-19 indicated that on 10 August 2021, more than 202.14 million people were infected with the virus, and over 4.28 million have died. Different variants of COVID-19 occur in different countries, namely, , , , and . The variant is the highest contagion, which is more transmissible than the variant and with a higher risk of hospitalization [3].
Tuberculosis (TB) is a communicable disease which is caused by Mycobacterium tuberculosis. It is a major cause of death from a single infectious agent. The bacteria that cause TB are usually spread when an infected person coughs, sneezes, or any other forceful expiratory maneuver that shears respiratory secretions. It mostly affects the lungs but can also affect other organs such as bones, brain, kidneys, and glands [4]. About a quarter of the worldwide population is infected with Mycobacterium tuberculosis and thus at risk of developing TB disease. Globally, approximately 10 million people fell ill with TB, and 1.45 million people died from TB in 2018 [5].
Both TB and COVID-19 are infectious diseases that are transmitted mainly via close contact. The growing clinical evidence suggests that TB diseases are associated with COVID-19 outcomes, including an approximately two to three-fold increase in mortality and a 25% relative decrease in the possibility of recovery for COVID-19 coinfection with current TB disease [6]. Studies on COVID-19 and TB coinfection are limited. It was forecast that people infected with both TB and COVID-19 may have poorer treatment outcomes [7]. TB patients should take precautions as suggested by health professionals to take care of COVID-19 and continue their TB treatment as directed. People infected with COVID-19 and TB show similar symptoms like cough, fever, and difficulty breathing [7, 8]. The coinfection of the two diseases exists when individuals are infected with both the COVID-19 and TB at the same time.
Modeling the transmission dynamics of diseases is a way to formulate what is known about the transmissions and explore all possible features of the disease with mathematical techniques [9]. The dynamical model study is used to investigate the vital parameters, predict future trends, and evaluate control measures to provide decisive information for decision making [10]. Recently, many mathematical models have been developed to study the transmission dynamics of COVID-19, see for example [11–16] and references cited therein. Mathematical models for TB are also extensively studied in [17–21] and references cited therein. The studies in [22] described the effect of speeding the detection time of COVID-19 infected individuals and reducing the transmission rate using both deterministic and stochastic models. In [23], the authors developed a mathematical to describe the dynamics of COVID-19 by taking into account all modes of transmission and prosocial awareness. Moreover, an optimal control strategy was proposed in [24] that helps to reduce the contact between TB infectious and exposed individuals.
The coinfection of various diseases is widely studied in the mathematical literature, for instance, TB and HIV/AIDS coinfection [4, 25], Hepatitis E and HIV [26], and cholera-schistosomiasis coinfection in [27]. On the contrary, the study on COVID-19 and other coinfection are at its infant stage, and related mathematical literature is scant. For instance, COVID-19 and malaria coinfection are discussed in [28, 29], while TB and COVID-19 coinfection presented in [30, 31]. The studies in [30] attempted to estimate the number of TB patients co-infected with COVID-19 with and without public health interventions in New Delhi, India. Their study showed that the peak of COVID-19 and TB coinfected patients would occur on the day in the absence of public health interventions and on the day with interventions. Besides, the coinfection of the MERS-CoV and TB model is analyzed in [31]. In their model, the authors studied the bifurcation analysis of the model and performed a sensitivity analysis of the basic reproductive number for the model parameters. To the best of our knowledge, mathematical models considering the coinfection of TB and COVID-19 are rare and they prompt us for further investigation.
In this work, we propose a new mathematical model for TB and COVID-19 coinfection by classifying the total population into eight compartments. The paper is organized as follows. A detailed introduction to both diseases has been discussed in Section 1. The model formulation of COVID-19 and TB coinfection is investigated in Section 2. The analysis of the full model is given in Section 3, and the mathematical results of the submodels are also investigated. In Section 4, numerical results are obtained and discussed with detailed for each submodel as well as the coinfection model, while our work is concluded in Section 5.
2. Model Formulation
In this section, we formulate the governing mathematical model for the codynamics of TB and COVID-19 by dividing the total population into eight mutually exclusive compartments denoted as susceptible individuals (), latent TB infected individuals who have no symptoms of TB disease and are not infectious (), active TB-infected individuals (), COVID-19 infected individuals with no symptoms but are infectious (), COVID-19 infected individuals with clinical symptoms and are infectious (), both latent TB and COVID-19 coinfected individuals (), active TB and COVID-19 coinfected individuals (), and recovered population from both TB and COVID-19 (). With these in mind, the total population at time , denoted by , is given by
We assumed that the susceptible population is increased by the recruitment at a rate . All population in each compartment suffer from natural death rate . Susceptible individuals acquire TB through contact with active TB patients by a force of infection as expressed in equation (3). In this expression, denotes the effective contact rate for TB infection. The latent TB infected individuals are assumed to be asymptomatic that does not transmit the disease [32]. Similarly, susceptible individuals acquire COVID-19 infection following effective contact with people infected with COVID-19 at a force of infection defined as in equation (4). Here, represents the effective contact rate for COVID-19 disease. Moreover, we assumed that individuals in latent TB compartment () leave to compartment at a rate , and to latent TB and COVID-19 coinfected class at a force of infection and some portion recovered at a treatment rate . Besides, Individuals in TB-infected class () become recovered from the disease at a rate while the remaining portion transfers to the coinfection compartment at a rate or die due to TB-induced death at a rate .
The population in compartment () either progress to coinfected class at a rate or die with COVID-19-induced death rate . The remaining individuals are assumed to be transferred to either compartment at a constant multiple of as illustrated in Figure 1. That is, the population in class moves to class with a rate of , to compartment class rate of , and become recovered at a rate of . Furthermore, we assumed that individuals in compartment leaves to compartments , , or , respectively, at a rate of , , or while coinfection-induced death rate is denoted by . In addition, the COVID-19 exposed individual () has the chance to leave either to compartment , , or , respectively, at a rate of , , or . Similarly, the population in compartment () is recovered at a constant rate of or transferred to co-infection compartment at a rate of . The disease induced death rate in this compartment is denoted by . The flow diagram of the proposed model is illustrated in Figure 1.
Based on the flow diagram, it results in systems of the following nonlinear differential equations:
where and with nonnegative initial conditions , , , , , , , and . The model parameters are described in Table 1.
2.1. Positivity and Boundedness of Solutions
As we are dealing with human populations, all solutions must be positive and bonded in a feasible region. To assure these, we have the following theorem.
Theorem 1. The solutions of the system of equation (2) are positive, unique, and bounded in the region
Proof. All functions on the right-hand side of equation (2) are on . Thus, by the Picard–Lindelf theorem [33], equations (2) has a unique solution. Let the system of equations in (2) be written in the form , where , and is the right hand sides of the equation (1). Following the results of proposition A.1 in [34], the functions has the property of
whenever . As there exists a unique solution for the system of equations (2), it follows that for all whenever .
The change of total population is governed by
The solution for this linear first order ode is . Thus, for the initial data , we obtain
Hence, the solutions of the system of equations in (2) exists, unique, and bounded in a feasible region , which concludes the proof.
3. Model Analysis
To understand the dynamics of the proposed model, we find the equilibrium points of the system and investigate the dynamics of the equilibrium points. The analysis will be done by investigating the behavior of the submodels for TB and COVID-19 as well as the coinfection model.
3.1. COVID-19 Submodel
Without considering the infections of people with tuberculosis, the COVID-19 submodel is given by
3.1.1. Local Stability of Equilibrium Points
The equilibrium point of the COVID-19 submodel is a solution to the system of equations:
The disease-free equilibrium point is obtained by putting and in equation (9). Hence, the disease-free equilibrium point, denoted by , becomes .
The basic reproduction number of the submodel defined as the average number of secondary infections produced by a single COVID-19 infected individual in a completely susceptible population. It is computed using the next generation matrix as given in [35]. Let be the input rate of newly infected individuals in the compartment and and , respectively, be the input rates by other means and rate of transfer of individuals out of a compartment , where the ’s are our state variables and
Hence, from equation (8), it follows that
The Jacobian matrices of and , respectively, are
Hence, the next-generation matrix becomes
The eigenvalues of the next generation matrix are found to be , and . The largest spectral radius of the next generation matrix denotes the basic reproduction number for the submodel. Thus, we have the following result.
Theorem 2. The basic reproduction number of the COVID-19 sub-submodel is We use the Jacobean matrix to determine the local stability of the equilibrium points. For the submodel (8), the Jacobean matrix becomes The Jacobean matrix of the submodel at the disease-free equilibrium point is given as Here, the eigenvalues for are , and the remaining eigenvalues can be easily obtained from the submatrix: In order to determine the local stability of the disease-free equilibrium point, it is sufficient to show that the trace of is less than zero, and its determinant is greater than zero.
The trace of , which is less than zero if
. Moreover, the determinant of
This value is greater than zero if , that is, if and if . Hence, the disease-free equilibrium point of the COVID-19 submodel is stable if and unstable if .
The endemic equilibrium point of the COVID-19 submodel is the solution of the system of equation in (18).
Now, solving this system of equations in terms one can obtain the following expressions for the endemic equilibrium point where
Here, , , and . Using (19) in the expression for , it gives . After simplification, we have in compact form
This implies that the force of infection is positive at endemic equilibrium point only if . Thus, we have just proved the following theorem.
Theorem 3. If , then the submodel system (8) has a unique endemic equilibrium point.
3.1.2. Global Stability of the Endemic Equilibrium Point
The global stability analysis of the endemic equilibrium point is analyzed using Lyapunov function. For this, we define the following:
The function is greater than zero, and it is equal to zero at the endemic equilibrium point . Differentiating the function with respect to time, we obtain
This holds true for which corresponds the existence of . Thus, implies that the function is a strictly Lyapunov function which indicates that the endemic equilibrium point is globally asymptotically stable. Epidemically, this result implies that COVID-19 will continue to exist in the human population for a long duration.
3.1.3. Sensitivity Analysis of the Submodel
In this section, we perform the sensitivity analysis of the model parameters for the COVID-19 submodel given in (8). The sensitivity of a parameter reflects how the model behavior responds to a small change in a parameter value, and it is defined as [35]
In our context, the sensitivity analysis of each parameters for equation (8) becomes
Here, we can observe that the contact rate positively influences the spread of COVID-19. Further, if , the transfer rate from exposed compartment to infected class has a positive impact on the spread of the virus. That is, increasing the value of these parameters will increase the endemic. The remaining parameters, , , , and , have negative impact, which means increasing the value of these parameters will decrease the number of people infected with COVID-19. However, in the study of sensitivity analysis, it is unethically increasing human mortality rate to control disease epidemic and hence do not considered.
3.1.4. Bifurcation Analysis for COVID-19 Submodel
Next, we study the solution behavior of equations (8) by taking as the bifurcation parameter. Calculating the value of from , i.e., , we have . Following the result given in Theorem 4 [36], we can calculate the eigenvalues of the Jacobin matrix at the disease-free equilibrium point by substituting . Thus, substituting in equations (15), it gives zero eigenvalue. This means that the Jacobean matrix in equation (15) at has a left eigenvector (associated with the zero eigenvalue) which is calculated from . Here, , where
Similarly, the right eigenvector associated with the zero eigenvalue can be calculated from with where
Now, assume that represents the right-hand side of the equation in the COVID-19 submodel (8) and let denote the corresponding state variable for . Define
The local dynamics of equation (8) near the bifurcation point are then determined by the signs of two associated constants and with . Note that, in , the first zero corresponds to the DFE, , for the submodel (4). In other words, , for , if and only if the right-hand sides of equations (8) are equal to zero at . Moreover, from , we have when , which is the second zero component in .
For the submodel (8), the associated nonzero partial derivatives at the disease-free equilibrium are
Then, using the above expressions for , it follows that
For the sign of , it can be shown that the associated nonvanishing partial derivatives are
Clearly, it also follows from the above expressions that
From the computations of and , we observe that is always positive, and the sign of the bifurcation coefficient depends on the minimum value of reinfection that induces bistability (). Thus, using Theorem 4 of [36], the following result is proved.
Theorem 4. If the minimum value of the reinfection that induces bistability is less than unity, then there exists a forward bifurcation for the COVID-19 submodel (8).
3.2. TB Submodel
In the absence of COVID-19 disease, equation (2) automatically reduced to TB-submodel expressed as
In a similar argument, the disease-free equilibrium point of equation (32) becomes , and the corresponding basic reproduction number is given by
Theorem 5. The disease-free equilibrium point of the TB submodel (32) is locally asymptotically stable if .
Proof. The Jacobean matrix at the disease-free equilibrium point for the TB submodel is given by
Clearly, the two eigenvalues are equal to , while the remaining can be obtained from the reduced matrix
Here, the trace of the matrix is , and its determinant is . Thus, the determinant is greater than zero if . Therefore, the disease-free equilibrium point of the TB submodel is locally asymptotically stable if . This completes the proof.
The endemic equilibrium point of the TB-submodel is obtained by solving the system of equations:
where . Solving this equation, we have
The simplified expression for in terms of , , , and becomes
Note that the endemic equilibrium point exists if and only if since .
Theorem 6. The endemic equilibrium point of the TB submodel (32) is globally asymptotically stable if
Proof. To proof this theorem, we define the following function: Next, we show that this function is a Lyapunov function. Differentiating the function with respect to time , we obtain For , the endemic equilibrium point exists. Thus, implies that the function is strictly Lyapunov function which deduces that the endemic equilibrium point is globally asymptotically stable.
3.3. Analysis of the Coinfection Model
The equilibrium point of the coinfection model (2) is obtained by solving the system of equation (41): where and are as given in equations (3), and (4), respectively. The disease-free equilibrium point () of the system becomes
In the following, we calculate the basic reproduction number () of the coinfection model using the method of next generation matrix [37]. For this, we use the notation to denote the vector functions: representing the appearance of new infections and the transfer of individuals into and out of the infected compartments, respectively. The Jacobian matrices of and at the disease free equilibrium point are, respectively, where , , , , , and . The next-generation matrix is then given by where
The largest spectral radius of the next generation matrix denotes the basic reproduction number of the coinfection model. Clearly, four of the eigenvalues of the matrix are equal to zero. The remaining eigenvalues are obtained from the reduced matrix Thus, computing the eigenvalues of one can easily obtain that where . Hence, the basic reproduction number of the coinfection model (1) is given by
3.3.1. Local Stability Analysis of the Disease-Free Equilibrium Point
The local stability analysis of the disease-free equilibrium point for the coinfection model is analyzed using the method of linearization. The Jacobean matrix of the model equation (2) at disease-free equilibrium point represented as
The local asymptotic stability of is analyzed based on the sign of the eigenvalues. Expanding the characteristic equation with the first and eighth columns, and the sixth row, we obtain three eigenvalues and . Further, we obtain the remaining five eigenvalues from the characteristic equation of the reduced matrix:
Here, the characteristic polynomial is computed using elementary row operations method [38]. Thus, the row reduced form of the matrix becomes
Hence, the characteristic polynomial of our matrix is the product of the diagonal entries times , where is the number of row swaps. We perform three row swaps, so that the characteristic polynomial is
Simplifying it gives
That is, the eigenvalues are the solutions of the characteristic equation: where
Thus, using the Routh-Hurwitz stability criterion, the roots of the equation (54) have negative real parts if the following condition holds.
Hence, we have proved the following result.
Theorem 7. The disease-free equilibrium point of the coinfection model (2) is locally asymptotically stable if the condition given in equation (56) holds.
3.3.2. Global Stability of Disease-Free Equilibrium
Let us rewrite the coinfection model (2) as where and . Here, the components of denotes the number of uninfected individuals while represents the number of infected ones [39]. The disease-free equilibrium is denoted now as
The conditions and below must be met to guarantee global asymptotically stability:
(H1) For , the equilibrium point is globally stable
(H2) for , where is a Metzler matrix (the nondiagonal entities of are nonnegative), and is the biological feasible region of the model
From the coinfection model (2), we have
Hence,
Therefore,
Here, , and hence
Thus, . This shows that condition () is not satisfied. As a result, and then the disease-free equilibrium point cannot be globally asymptotically stable.
3.4. Impacts of TB on COVID-19
To analyze the impact of TB disease on COVID-19 (and vice versa), we began by expressing the basic reproduction numbers in terms of (and vice versa) [40]. Expressing the parameter in the equation (33) in terms of , we have
Solving for and simplifying, we obtain
Substituting this value of in (13), we obtain in terms of as where , and .
The partial derivative of with respect to is then given by with , and
Here, we observe that is positive. This implies that the expansion of TB infection exacerbates the pandemic of COVID-19.
Remark 8. If , the expansion of TB in the community has no significant impact on the spread of COVID-19. In the contrary, when , the increase in TB epidemic will negatively influence the spread of COVID-19.
Similarly, the impact of COVID-19 for TB cases can be addressed by expressing in terms of and identifying the sign of the partial derivative of with respect to .
4. Numerical Simulation
In the previous sections, we have discussed the analytical behaviors of the model equation (2) and the submodels. In the following, we perform numerical solutions of submodels as well as the coinfection model. The solutions of the system equations are integrated using the ode45 solver in MATLAB. The parameter values are also described in respective subsections.
4.1. Numerical Solutions of the COVID-19 Submodel
The analytical findings of the submodel (8) are illustrated in section 3.1. We now present the numerical solutions using initial conditions for state variables , , , , and the parameter values described in each figure.
Figure 2 presents the solution of the COVID-19 submodel using the parameter values described in the caption. The basic reproduction number is obtained as , which implies the disease free equilibrium point is stable, and no endemic equilibrium point exists. Using different initial conditions for each compartment, all solution curves converge to the disease-free equilibrium point as seen in Figure 3. In Figure 2, the susceptible populations approach to the value of , while the other compartments go to zero, which supports our analytical findings in the previous sections.
(a) Susceptible
(b) Infected
(c) Recovered
Figure 4 illustrates the time serious plot of the COVID-19 submodel, for which its basic reproduction number is greater than unity. In the analytical findings, we have shown that the disease-free equilibrium point is unstable, and the endemic equilibrium point is stable. Using different initial conditions for each compartment, all solution curves converge to the endemic equilibrium point as seen for the susceptible, infected, and recovered plots of Figure 5. As time goes, in Figures 5(a)–5(c), the solutions close to the endemic equilibrium point given in equation (19) which agrees with the analytical properties.
(a) Susceptible
(b) Infected
(c) Recovered
Figures 6 and 7 show the effects of the COVID-19 transmission rate and the recovery rate of the COVID-19 infected individuals. The other parameter values used in the simulation of the figures are , , , , , and in Figure 6, while in Figure 7, . When the transmission rates are increasing, the populations in each compartment (exposed, infected, and recovered) become increasing accordingly. If there are low contact rates, then the exposed and infected individuals are low and as a result, the number of recovered individuals is also low.
(a) Exposed
(b) Infected
(c) Recovered
(a) Exposed
(b) Infected
(c) Recovered
4.2. Numerical Solutions of the TB Submodel
The numerical solutions of the TB submodel with the basic reproduction number of are given in Figure 8. The solutions are executed using the initial conditions of the state variables , , , and . We observe that all solution curves approaches to their endemic equilibrium point. This supports the analytical results that corresponds to , where the endemic equilibrium point of the TB submodel is globally asymptotically stable. Figure 9 shows the impact of the contact rate () of TB endemic. As it increases, the number of individuals in the compartments, latent TB, active TB, and recovered ones is also increasing. If there is a low contact rate, then the transmission of TB disease will decrease as shown in Figures 9(a)–9(c).
(a) Latent TB
(b) Infected TB
(c) Recovered
4.3. Numerical Solutions of the Coinfection Model
In this subsection, we present the numerical solutions of the coinfection model (2). In doing so, we use initial values of the state variables for the coinfection model given in Table 2, and all parameter values are illustrated in Table 3.
The time series graphs for the numerical solutions of the coinfection model (2) are plotted in Figure 10 with parameter values given in Table 3. The impact of treatment rate for the coinfected individuals with both diseases is illustrated in Figure 11. It was shown that an increase in the treatment rate for the disease will decrease the number of individuals in that class.
Figure 12 displays the impact of transferring rates to coinfected with TB and COVID-19 from each actively infected individual of the two diseases. In Figure 13, we have shown that the impact of contact rates of the disease with infected compartments of system (2). Thus, we conclude that the states , , , and increase as the transmission coefficients increase.
(a) Coinfected
(b) Infected COVID-19
(a) Coinfected
(b) Infected TB
(c) Infected COVID-19
5. Conclusion
In this study, a nonlinear deterministic mathematical modeling for the coinfection of TB and COVID-19 is proposed to examine the possible spread of these diseases. The biological meaningfulness of the model is proved by showing the existence, uniqueness, positivity, and boundedness of solutions in a given region. Then, the equilibria are computed for the coinfection model and each submodels independently. The stability analysis of these equilibria is also presented with the help of the calculated respective basic reproduction numbers. The analytical results in both submodels reveal that the disease-free equilibrium points are stable if the basic reproduction numbers are less than one otherwise unstable. Furthermore, there exists a stable endemic equilibrium point for the basic reproduction numbers greater than one. Besides, the stability analysis of coinfection disease-free equilibrium point is investigated via corresponding basic reproduction number. Different simulation cases were performed to supplement the analytical results, and it is observed to be in good agreement. Furthermore, the simulation result reveals that an increase in contact rate exacerbated the coinfection of both diseases, while a decrease in the effective contact rates and increase in treatments would minimize the spread of co-infection diseases.
Data Availability
All the necessary data was included in the main text.
Conflicts of Interest
The authors would declare that they have no known competing interests that could have appeared to influence this work.
Acknowledgments
The first author acknowledges Adama Science and Technology University for the financial support under grant number ASTU/SP-R/120/21.