Mathematical Modelling of the Epidemiology of Tuberculosis with Silicosis Coinfection

Mathematical Modelling of the Epidemiology of Tuberculosis with Silicosis Coinfection

Innocent Cleopas Eli*

Department of Mathematics and Statistics, Federal University Otuoke, Bayelsa State, Nigeria

*Corresponding Author

DOI: https://doi.org/10.51584/IJRIAS.2025.10020051

Received: 10 February 2025; Accepted: 14 February 2025; Published: 20 March 2025

ABSTRACT

The study presents an innovative mathematical model analysing the epidemiology of Tuberculosis with silicosis coinfection. It effectively integrates epidemiological factors and historical theoretical research with well-structured model formulation and numerical verification through MATLAB. The use of partial differential equation, Jacobian matrix, deterministic techniques as well as Routh Hurwitz algebraic criteria plays significant role in the stability of disease-free equilibrium point and stability of the endemic equilibrium point analytically which indicates locally stable system asymptotically which indicates locally stable system asymptotically as equally demonstrated by the reproduction number. The solutions of the model equations are integrated using the Range Kutta Fourth order method in MATLAB and observed the impact of 

 which proves that the endemic equilibrium point increased for the recorded population meanwhile, decreased for the coinfected population as increases. Since the Ro < 1, it shows that the disease-free equilibrium point is stable beyond 2500 days and no endemic equilibrium point exists. It is equally observed that the solution trajectories of the silicosis only sub-model converge to a single point believe to be disease free equilibrium point also known as silica free movement

Keywords: Mathematical Modelling, Epidemiology, Tuberculosis, Silicosis and Coinfection

INTRODUCTION

Tuberculosis is an infectious disease caused by bacteria that mostly affects the lungs, it remains a global health issues due to its high mortality and it is the leading cause of death in the majority of Sub-Sahara Africa countries [1] and South-East Asia countries.

Tuberculosis cases increase in the year 2013 to 2015. Nigeria came third behind India and China in the new tuberculosis census [2]. The effect of tuberculosis is still high and devastating till date [3]. Tuberculosis is caused by the bacterium mycobacterium tuberculosis and primarily affects the lungs. It is transmitted through the inhalation of airborne droplets containing the bacteria, making it highly contagious. In recent years, an emerging concern in the epidemiology of tuberculosis is its co-occurrence with silicosis.

Silicosis is a debilitating occupational lung disease caused by inhalation of crytallin silica dust [4]. Workers in various industries such as mining, construction and manufacturing are at risk of developing silicosis due to their exposure to silica dust. The link between silicosis and TB is well-established, as silicosis weakens the immune system, making individuals more susceptible to TB infection and increasing the risk of progression from latent TB infection to active TB disease. Furthermore, the co-existence of TB and silicosis present unique challenges in diagnosis treatment and prevention. The co-effective poses a significant challenge for healthcare systems as the treatment and management of individuals with TB and silicosis require a more comprehensive and nuanced approach [5]. [6] studied the dynamical behaviour of epidemiological models with non-linear incidence rates. Their studies show that models with non-linear incidence have a much wider range of dynamical behaviours than those with bilinear incidences rates, far these models, there is a possibility of multiple attractive basins in phase space and because of that, the disease survival depends not only upon the parameters but also upon the initial conditions. [7] studied the SEIR model with non-linear incidence rates in epidemiology. The purpose of their paper is to prove the global stability of the non-trivial equilibrium for their SEIR method. [8] also studied the global dynamics of a SEIR model with varying population size using the homogeneity of the vector field of the model to analyse the derived system of the fractions in determining the behaviour of the population sizes and the total population. [9] investigate a SEI TB model with immigration that includes infected (both Latent and Infectious) individuals. The model assumed constant recruitment with fixed fractions entering each class, however, they were able to prove that under certain restrictions on the parameters (including the treatment rates, disease transmission rate and TB induced death rate) the disease will approach a unique endemic level. [10] investigate the impact of immigration on the transmission dynamics of tuberculosis and showed that disease never dies out but becomes endemic in host areas, also, disease will persist in the population if there is a prevalence of TB immigrants and there will be no disease free-equilibrium [11] presented a SEIR tuberculosis model which incorporated treatment of infectious individuals and chemoprophylaxis to show that treatment of infectious individuals is more effective in the first years of implementation as it cleared active TB immediately. [12] investigated the importance of Heterogeneity to the epidemiology of Tuberculosis using mathematical modelling to simulate the spread of TB in a population focusing on accounting for heterogeneity in TB dynamics. Their results highlight that TB transmission is not uniform across the population and understanding these variations is critical for effective intervention strategies. A study by [13] focuses on the healthcare seeking behaviour of TB patients and the time for them to access care, it sheds light on the challenges and delays in the diagnostic and treatment process in a specific region. [14] presented a projected effects of tobacco smoking on worldwide tuberculosis control, their work predicted that smoking would produce an excess of 18 million tuberculosis cases and 40 million death from TB between 2010 and 2050 if smoking trends continued along current trajectories. [15] presented a mathematical model in 2019, to assess vaccination and effective contract rate impact in the spread of tuberculosis via deterministic epidemic model (SV ELI – (Susceptible, Vaccinated, Early Latent, Late Latent, Infectious). Using Lyapunor-Lasalle method to analyse the epidemic system (SVELI) around the equilibrium (disease free and endemic), they discovered that the global asymptotic stability of the unique endemic equilibrium whenever is proved, and when it is less than 1, TB can be eradicated. Also a TB date found in a literature related to Cameroon, shows that vaccination coverage is not sufficient to control TB, effective contact rate has a high impact in the spread of TB. [16] presented mathematical models of the population dynamics of TB for the effects of efficient treatment Versus Incomplete Treatment. His studies analyzed the spread, asymptotic behaviour and possible eradication of the disease versus persistence in TB. [17] study and presented a mathematical model to evaluate the impact of the response of TB cells and macrophages in the control of MTB (multi-drug resistant tuberculosis). Their analysis reveals the existence of two equilibrium states, infection-free equilibrium and endemically infected equilibrium which can represent a state of latent or active infection, depending on the amount of bacteria. [18] in their study, modelled the qualitative behaviour of a system of ordinary differential equations and a system of differential integral equations for the dynamics of disease transmission for tuberculosis and discusses. The possibility of a person infected with TB may develop active TB as a result of endogenous infection. [19] studies an optimal control TB mathematical model this include the presence of exogenous reinfection in the dynamics of the disease, modifying the model by [18] above. [20] presented different mathematical models and biological scales in understanding the immune response in tuberculosis. In their study, they use four different mathematical tools to explore both the global immune response as well as the more local one (granuloma formation) and compare and contrast results obtained using these methods. [21] presented a system of ordinary differential equations modelling the population of infective. [22] Investigate a first time comprehensive review of work on within host TB model that describe the immune response of the host to infection, including the format of lung granulomas. The survey application of this models to TB therapy and prevention suggest future directions to impact this global disease TB. [23] focus on the study of an age-structure model for TB transmission dynamics in populations subjected to a vaccination program. They use the theoretical results to vaccination policies to determine the optional age or ages at which an individual should be vaccinated, proving the existence of an endemic steady state when commonly used method does not apply and showing how to compute the optimal vaccination strategies in such situations. [24] presented a mathematical model and simulation to control the spread of Multidrug-resistant Tuberculosis, analyzing a suitable strategy in controlling the development of susceptible individuals to active TB and even multidrug-resistant TB. But in this paper, we will consider the Mathematical modelling of the epidemiology of Tuberculosis with silicosis co-infection.

Mathematical Formulation

Model Description 

Our proposed model is an extension of the models proposed by [25] and [26]. 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 and silicosis through contact with active TB and silicosis patients (including dust inhalation) by contact rates 

 respectively. Individuals exposed to TB and silicosis are transferred into the co-exposed class at rates  respectively. Exposed TB individuals could either recover from TB disease at rate  or become actively infected with TB ate rate  Active TB infected individuals could either recover from the disease at rate  or acquire silicosis infection at rate  Since Silicosis infection does not have a cure the co-infected population recover from TB only at a rate of  We assume that silica is found throughout the earth’s crust and it is harmless until disturbed in a way that creates dust by a constant dust production rate M. Silica dust is lost a rate of  Silicosis exposed individuals are transferred to the infected class at a rate  while infected silicosis individuals become co-infected at a rate  Death due to TB, silicosis and the co-infection are denoted by 

Flow Diagram of the Proposed Model

The flow diagram of the proposed model is given in figure 2.1 below. The parameters used in the model description are also given below.

Where;

 Effective contact rate for the Silicosis submodel

 Effective contact rate for the TB

 Recruitment rate into Susceptible class

 Fraction of TB exposed individuals to be TB infected

 Fraction of Silicosis exposed individuals to be Silicosis infected

 Silicosis infection rate from TB infected individuals

 TB infection rate from Silicosis infected individuals

 Rate at which TB exposed individuals become exposed with Silicosis

 Rate at which Silicosis exposed individuals become exposed with TB

 

 Rate at which individuals leave the co-exposed class

 Recovery rate of TB infected individuals

 Recovery rate of the coinfected class

 Natural death rate

 Death rate due to TB disease

 Death rate due to Silicosis infection

 Death rate due to the coinfection

 Silica dust deposition rate

 Silica dust production rate

METHOD OF ANALYSIS

The following system of nonlinear differential equations is derived from the flow diagram.

         (3.1)

         (3.2)

        (3.3)

       (3.4)

        (3.5)

        (3.6)

        (3.7)

        (3.8)

         (3.9)

    (3.10)

With 

 and 

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 behaviour of the sub-models for TB and Silicosis as well as the coinfection model.

TB Sub-Model

Without considering the infections of people with Silicosis, the TB sub-model is given as

 

         (3.11)

         (3.12)

         (3.13)

          (3.14)

        (3.15)

Disease-Free Equilibrium Point (DFEP)

The disease free equilibrium point are the points

To get S*, we solve the system (3.11 

 3.14) when

We get,

          (3.16)

         (3.17)

         (3.18)

          (3.19)

From equation 3.16, we get

Then the DFEP becomes

        (3.20)

Stability of the Disease-Free Equilibrium Point

Theorem 3.3.1: The DFEP is locally asymptotically stable if the basic reproduction number is lesser than one and unstable otherwise.

Proof: 

Let 

 represent the basic reproduction number of the TB sub-model at DFEP.

We obtain 

 by the next generation matrix proposed by [24].

Let 

 be the terms which contains only secondary infections and  be the other terms which do not contain secondary infections. Then  is the spectral radius of 

Where 

 and 

 , and  

At DFEP

 and    

Then,

 

           

 

 

The DFEP equilibrium is asymptotically stable if

Endemic Equilibrium Point

The endemic equilibrium point of the TB sub-model is the solution of the system of equation (3.16 

 3.19).

From equation 3.18,

          (3.21)

Substitute equation 3.21 into equation 3.17

 then 

And

            (3.22)

Substitute equation 3.22 into equation 3.16, we get

Then,

       (3.23)

Substitute equation 3.23 into equation 3.21, we get,

       (3.24)

Substitute equation 3.23 into equation 3.19, we get,

      (3.25)

Stability of Endemic Equilibrium Point

Theorem 3.5.1: The endemic equilibrium point is asymptotically stable if all eigenvalues of its characteristic polynomial are negative.

Proof:

Let,

Then, the Jacobian matrix is given as,

        (3.26)

At EEP,

            (3.27)

On simplifying equation 3.26 we get,

        (3.28)

Where,

By Routh Hurwitz stability criterion, equation 3.27 has all eigenvalues 

 since all coefficients satisfy  and 

Silicosis Sub-Model

Without considering the infections of people with TB, the Silicosis sub-model is given as,

          (3.29)

         (3.30)

         (3.31)

         (3.32)

Disease-Free Equilibrium Point

From the dynamics of Silicosis, if there is no silica production in the community, then it means that there is no Silicosis patients in the community and clearly M = 0. Then 

          (3.33)

         (3.34)        (3.35)                                                                                                                               (3.36)

From equation 3.32

The DFEP is given as 

Stability of the Disease Free Equilibrium Point

Theorem 3.7.1: The disease free equilibrium point is locally asymptotically stable if all eigenvalues of the characteristic polynomial of its Jacobian matrix are all negative.

Proof:

Let,

By equation 3.26,

Since all eigenvalues are negative, the DFEP is locally asymptotically stable.

Theorem 3.7.2: The DFEP is locally asymptotically stable if the basic reproduction number is lesser than 1.

Proof: 

Let 

 represent the basic reproduction number of the Silicosis sub-model at DFEP. We obtain  by the next generation matrix method.

  and 

 and 

The DFEP will be locally asymptotically stable if 

Determination of Endemic Equilibrium Point

The endemic equilibrium point is the solution to the system of equations (3.33 

 3.36)

From equation 3.36

           (3.37)

Substitute equation 3.37 into equation 3.33

On simplifying, we get a quadratic equation of the form

       (3.38)

We solve equation 3.38 and take only the positive root, we get,

 

     (3.39)

Substitute equation 3.39 into equation 3.37, we get

     (3.40)

Substitute equation 3.39 and equation 3.40 into equation 3.34, we get

          (3.41)

Substitute equation 3.41 into equation 3.31, we get,

         (3.42)

Stability of the Endemic Equilibrium Point

Theorem 3.9.1: The endemic equilibrium point is locally asymptotically stable whenever all eigenvalues of the characteristic equation of its Jacobian matrix are negative.

Proof: 

On solving and simplifying this determinant equation, we get,

        (3.43)

Where,

 

 

By Routh Hurwitz stability criterion, equation 3.27 has all eigenvalues 

 if all coefficients satisfy  and 

Analysis of the Coinfection Model

Determination of the Disease Free Equilibrium Point

The disease free equilibrium will occur when there is no cause of infection, and by consequence

We let, 

Equations (3.1 

 3.9) becomes

          (3.44)

         (3.45)

         (3.46)

        (3.47)

        (3.48)

         (3.49)

        (3.50)

        (3.51)

          (3.52)

At DFEP,

DFEP = 

Stability of the Disease Free Equilibrium

Theorem 3.12.1: The disease free equilibrium point is locally asymptotically stable if the basic reproduction number 

Proof: we obtain 

 by the next generation matrix method.

Then 

GRAPHICAL RESULTS

We now present graphical results using Runge Kutta fourth order method in MATLAB and the parameter values are described in the respective subsections.

Figure 4.1: Solutions of the TB sub-model with parameter values 

 with basic reproduction number 

Figure 4.2: The solution curves of the susceptible, infected, and recovered populations for different initial values of the TB only model at DFEP

Figure 4.3: Solutions of the TB sub-model with parameter values 

with  basic reproduction number 

Figure 4.4: The stable solution curves of the susceptible, infected, and recovered populations for different initial values of the TB only model using parameter values of the Figure 4.3 in which its basic reproduction number is greater than unity ( 

 ).

Silicosis Sub-Model

We present the numerical solutions to the Silicosis sub-model using initial conditions for state variables 

 and the parameter values described in each figure.

Figure 4.5: Numerical simulations showing the solution trajectories of the Silicosis sub-model converging to a single point believed to be the disease-free equilibrium (or silica-free environment).  (a) M = 0 and 

 (b) M = 0 for infected population at DFEP (c)

Effect of Parameter M on the fraction of infected population at endemic equilibrium. The parameter values are, 

TB-Silicosis Coinfection Model

In this subsection, we present the numerical solutions of the coinfection model. In doing so, we used initial values of the state variables for the coinfection model given as 

and  The parameter values used are  

Figure 4.6: Solutions of the TB-Silicosis coinfection model at DFEP with M = 0 and basic reproduction number 

Figure 4.7: Solutions of the TB-Silicosis coinfection model at EEP with M = 0.0001 and basic reproduction number 

Sensitivity Analysis of the Coinfection Model

In this section, we perform the sensitivity analysis of the model parameter 

 for the TB-Silicosis coinfection model. The sensitivity of a parameter reflects how the model behaviour responds to a small change in a parameter value, and it is defined as in [24].

Figure 4.8: Effect of change in 

 on the recovered and coinfected populations of the combined model.

DISCUSSION

Figure 4.1 presents the solution of the TB sub-model using the parameter values described in the caption. The basic reproduction number is obtained as 

 which implies the disease free equilibrium point is stable beyond 2500 days, 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 4.2. In Figure 4.1, the susceptible populations approach to the value of Λ/μ, while the other compartments go to zero, which supports our analytical findings in the previous sections.

Figure 4.3 illustrates the time serious plot of the TB sub-model, for which its basic reproduction number is greater than one. 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 4.4. As time goes, in figures 4.4(a)–4.4(c), the solutions get close to the endemic equilibrium point, which agrees with the analytical properties.

In figures 4.5(a) and 4.5(b), it is observed by testing different initial conditions, that the solution trajectories of the Silicosis only sub-model converge to a single point believed to be the disease free equilibrium point (or silica-free environment) while in figure 4.5(b)  compared different silica dust production rate on the silica and also observed that all curves converge towards the endemic equilibrium point.

Figure 4.6 shows the solution curves of the coinfection model where the disease-free equilibrium is stable, while figure 4.7 shows the solution curves where the disease-free equilibrium is unstable and the endemic equilibrium is stable.

Figures 4.5 and 4.6 show the effects of the TB transmission rate 

 on the recovered population and TB-Silicosis coinfected individuals.

CONCLUSION

This work formulates mathematical model to study the transmission pattern of tuberculosis and silicosis coinfection. Through analytical studies, it is determined that the disease free equilibrium point and endemic equilibrium point through which establishment of stability of the TB-only sub-model, Silicosis-only sub-model and the coinfection model were achieved. It is also determined that the reproduction number of each sub-model and the coinfection model, which was also key in establishing stability. The investigation went further to numerically simulate the results of the proposed model and obtained the time series plots for each sub-model and the coinfection model. It is observed that the numerical results agreed with  analytical findings. Lastly, sensitivity analysis on the model parameter 

 to determine the model response to it was carried out and observed that the endemic equilibrium point increased for the recovered population but decreased for the coinfected population as   increases.

REFERENCES

1. Ojo, M. M; Benson, T. O; Peter, O. J. & Goufo, E.F.D (2022), Non-linear optimal control strategies for a mathematical model of COVID-19 and influenza co-infection. Phys A state Mech. Its Appl. 607,128173.
2. World Health Organization (WHO) (2021). Global Tuberculosis Report, WHO Report.
3. Parsons, L.M., Somoskovi, Aji Gutierrez, C.; Lee, E. (2011). Laboratory diagnosis of tuberculosis in resource-poor countries: challenges and opportunities. Clin. Microbiol. Rev. 24(2): 314.
4. Rees D, & Murray J. Silica, (2007). Silicosis and tuberculosis. Int. J. Tuber Lung Dis. 11(5): 474484.
5. Yang Q; Lin, M.; He, Z.; Liu, X.; Xu, Y.; Wu, J.; Sun, F.; Jiang, T.; Gao, Y.; Huang. X.; Zhang, W.; Ruan, Q.; & Shao, L. (2022). Mycobacterium tuberculosis Infection among 1,659 silicosis patients in Zhejiang Province, China. Microbiology Spectrum, 10(6).
6. Liu, W.M., Hetchote, H.W. & Levin, S.A (2022) Influence of non-linear incidence rates upon epidemiological models. J. Matt. Biol., 23: 187-204.
7. Li, M.Y. & J. Muldowney, (2019). Global stability for the SEIR model in epidemiology Math. Biosci, 12, 155-164.
8. Li, M.Y; Graef JR.; Wang L.C; & Karsai J. (2020). Global dynamics of a SEIR model with a varying total population size. Math. Biosci 160: 191-213.
9. McChuokey Connell, C. & Van Den Driessche, P. (2023). Global Analysis of Two Tuberculosis Models. Journal of Dynamics and Differential Equations, 16(1), 0100-0139.
10. Jia ZW, Jia XW, Liu YX, Dye C, Chen F, & Chen F, & Chen CS (2016) Spatial analysis of tuberculosis cases in migrants and permanent residents. Emerg. Infect Dis; 14:1413-9.
11. Bhunu, C.P.; Garira, W.; Mukandavire, Z.; & Zimba, M. (2008). Tuberculosis Transmission Model ith chemoprophylaxis and Treatment Bulletin of Mathematical Biology, 70, 1163-1191.
12. Traver, J.M; Dodd, P.J. & Gomes, M.G.M. (2019). The importance of Heterogeneity to the Epidemiology of Tuberculosis. Clinical infectious Diseases 69(1), 159-166.
13. Yimer, S.A., Bjune, G.A.; & Holm-Hansen, C. (2019). “Time to first consultation, diagnosis and treatment of TB among patients attending a referral hospital in northwest Ethiopia”. BMC Infectious Diseases, 19(1), 814.
14. Basu, S.; Stuckler, D.; Bitton A; Glantz, S.A. (2011). Projected effect of tobacco smoking on worldwide tuberculosis control: Mathematical modelling analysis. Bmj, 343.
15. Nkamba, L.N., Manga, T.T; Agouanet, F., & Mann Manyombe, M.L. (2019). Mathematical model to assess vaccination and effective contact rate impact in the spread of tuberculosis. Journal of Biological Dynamics 13(1), 26-42.
16. Ayodeji O. Adebiyi (2016). Mathematical modelling of the population dynamics of tuberculosis. University of the Western Cape, 1-92.
17. Lourdes Esteva, Eduardo Ibarguen-Mondragon, Leslie Chavez-Galan (2011). Mathematical biosciences and engineering & (4), 973-986.
18. Zhilan Feng, W. Hung & Carlos Castillo-Charez (2001) on the role of variable latent periods in mathematical model for Tuberculosis. Journal of Dynamic and Differential Equations 13(2), 425-452.
19. Rachik, M.; Hattaf. K.; Saddis, S. Tabit Y. & Yousfi N. (2009), optimal control of tuberculosis with exogenous reinfection. Appl. Math. Sci. 3(5): 231-240.
20. David Gammack, Suman Garrigue, Simeon Mario, Jose Segovia-Juarez, & Denise E. Kirschner (2005). Understanding the immune response in Tuberculosis using Different Mathematical Models and Biological Scales. Society for Industrial and Applied Mathematics, 3(2), 312-345. DOI:101128/CMR. 00059-10.
21. Okuonghae, C. & Aihie, V.U. (2010). Optional control measures for Tuberculosis. Mathematical models including immigration and isolation of infective. Journal of Biological systems 18(1), 17-54.
22. Elseje Pienaar, Denise Kirschner, Simeon Marino & Jennifer J. Linderma (2017). A Review of Computational and Mathematical modeling contributions to our understanding of Mycobacterium tuberculosis within-host infection and treatments. PMC Pub.Med, 3, 170-185.
23. Carlos Castillo-Chavez, & Z. Feng. (1998). Global stability of an age-structure model for TB and its applications to optimal vaccination strategies. International Journal of Mathematical Biosciences 151, 135-154.
24. Van den Driessche, P.; & Watmough, J. (2002). Reproduction number and sub-threshold endemic equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180(1-2), 29-48.
25. Tenkam, H. M., Goufo, E. D., Tsanou, B., Hassan, A. S., Hussaini, N., & Terefe, Y. A. (2020). Classical and fractional analysis of the effects of Silicosis in a Mining Community. Alexandria Engineering Journal, 59(4), 2683-2694.
26. Abdul Halim, N. (2013). Tuberculosis Model: A Mathematical Analysis (Doctoral dissertation, University of Malaya).