Abstract
The main goal of this article is to devise the spatial-temporal spread of TB, in multiple neighboring domains, taking into account the epidemiological diversity of their populations. However, since both the environment and any population are spatially heterogeneous, it is obviously desirable to include spatial structure into an epidemic model. Individuals with tuberculosis can spread the disease by moving from one area to another. In addition, people travel by air between cities, so diseases can be spread quickly between very distant places (as was the case with the COVID-19). In our model, each region’s studied population is divided into five compartments S, L1, I, L2, and R. Further, we introduce in our discrete systems three control variables which represent the effectiveness rates of vaccination, travel-blocking operation, and treatment. We focus in our study to control the outbreaks of an epidemic that affects a hypothetical population belonging to a specific region. Firstly, we analyze the epidemic model when the control strategy is based on the vaccination control only, and secondly, when the travel-blocking control is added, we finish with the introduction of the treatment control. The optimal control theory, based on Pontryagin’s maximum principle, is applied thrice in this paper, for the characterizations of the vaccination, travel-blocking, and treatment controls. The numerical results associated with the multipoint boundary value problems are obtained based on the forward-backward sweep method combined with progressive-regressive Runge–Kutta fourth-order schemes.
1. Introduction
Conventional mathematical models have focused only on the temporal spread of infection but ignored or neglected spatial dynamics, although the spatial spread of the epidemic was actually observed and many infectious diseases were transmitted from one region to another. Since an infectious case is first found at one location and then the disease spreads to other areas, the spatial variation of the population plays an important role in the control and eradication of this disease. However, due to the large mobility of people within a country or even worldwide, models that ignore the effect of space are not sufficient to give a realistic picture of the spread of the diseases. The epidemic could spread over large and remote areas. It can reach continents, and among the examples in human history there is the black death (plague) which spread in Europe in the 13th century, then measles and smallpox between the 16th and 17th centuries which spread in the New World [1, 2], HIV/AIDS in the year 1981, the West Nile virus that emerged in North America in the late 1990s [3], SARS that appeared in Asia in 2003 [4], and recently the new corona virus (COVID-19) that appeared at the end of 2019 [5]. It continues to spread still in the world and jumps from one continent to another, especially due to the explosive development of means of air, sea, and land transport, which contributed to a wide spread. The overlap of many different components makes modeling the spatial spread of infectious diseases a complex task. And spatial heterogeneity requires formulating a small model for each region (city, country, and so on). One possible approach is to consider people traveling between distinct geographic regions (subfields of the field of study), assuming that contagion does not occur while traveling. With regard to diseases which affect animals and insects (vectors) and which can also be transmitted to humans [6], it is very important to define a precise and clear pattern of not only the contact/mixing between humans and different species but also the movement patterns of each species involved, including humans. The importance of modeling the spatial spread of infectious diseases can be highlighted by the outbreaks of corona disease in many countries around the world. Before this disease was diagnosed at the end of 2019, infected people were moving from city to city, and after diagnosing some infected cases and identifying the problem, it took some time for the political decision to prevent all movement, especially from the affected areas, and the regions were actually divided according to the number of infected cases it contains (from one affected continent to another, or from one affected country to another, or from one affected city to another …) [7]. Thus, in order to control spatial spread of the disease, it was necessary to follow all precautions and take all necessary measures to reduce displacement to and from the affected areas.
The mathematical modeling approach is one of the useful tools and can be used to study the dynamical behavior of a disease in a better way and to develop some useful strategies for possible eradication and control. In the last few years, several transmission models based on different approaches have been suggested to analyze the dynamics of TB disease in various regions. Most of the researchers who have studied mathematical models of tuberculosis have focused on the continuous-time models described by ordinary differential equations (see [8–15]). We observe that most of those researchers focused on the continuous-time models described by the differential equations. It is noted that, in recent years, more and more attention has been given to discrete-time models (see [16–19] and the references cited therein). The reasons for adopting discrete modeling are as follows. Firstly, the statistical data are collected at discrete moments (day, week, month, or year). So, it is more direct and more accurate and timely to describe the disease using discrete-time models than continuous-time models. Secondly, the use of discrete-time models can avoid some mathematical complexities such as choosing a function space and regularity of the solution. Thirdly, the numerical simulations of continuous-time models are obtained by the way of discretization. The deterministic, spatial disease models have the advantage that they allow for the investigation, either analytically or numerically, of various spatiotemporal control strategies. The oldest control method for the spatial spread of infectious diseases is quarantine through travel restrictions.
A study of the impact on disease dynamics of this type of quarantine was undertaken by Sattenspiel and Herring [20], as applied to the model of Sattenspiel and Dietz [21]. Numerically, this intervention method was shown to have little impact on the ultimate course of the epidemic.
There is increasing interest in studying the spatial spread of infectious diseases [22–25]. Most of the models studied have been partly continuous because of their mathematical tractability. In Allen et al. [26], a discrete-time, age-independent SIR-type epidemic model with n subgroups is formulated and analyzed and then applied to a measles epidemic model on a university campus. In similar pattern, Zakary et al. [7] devised a discrete-time SIR model that describes the propagation of a disease in a population of individuals who travel between p regions (domains). In the same way, in our article, we will study the propagation of tuberculosis in discrete-time case between several regions. Tuberculosis is a contagious disease, secondary to infection with “bacillus of Koch” (Mycobacterium tuberculosis). This bacterial agent is transmitted by air via the droplets contaminated by the bacterium, which are suspended in the exhaled air by patients, especially during coughing. Inhaling a small number of contaminated droplets is enough to infect an individual. The displacement of populations (travelers and refugees) has largely contributed recently in the spread of the disease in the world. People who are more likely to acquire TB infection are the following:(1)People recently exposed to someone who has symptomatic TB disease.(2)People who live in congregate settings with high-risk persons.(3)People who live or have lived in countries where TB is common.(4)People who are healthcare workers who are in contact with TB patients when proper infection control procedures are not followed.
Many people who acquire TB infection do not have symptoms and may never develop TB disease. These people have latent TB infections (LTBIs) [27]. After exposure to the tuberculosis bacillus, some people developed a primary infection, which is controlled by the immune system in 90% of cases: tuberculosis is then said to be “latent.” The bacillus remains in the body, but the immune system prevents it from multiplying. In 10% of those infected, the bacillus is not sufficiently controlled by the immune system, and these people develop a form of so-called “active” tuberculosis, which will cause disease and complications.
The Moroccan Ministry of Health registers nearly 30,000 cases of tuberculosis each year, including new cases and cases of relapse. The incidence rate is around 87 cases per 100,000 inhabitants, of which pulmonary tuberculosis represents half. They also noted that this disease affects young people aged 15 to 45, adding that among the causes of the spread of this disease are socio-economic factors, especially housing conditions, poverty, and malnutrition [28, 29].
In our model, the population total is divided into five categories:(i)(S): susceptible, who have never encountered the Mycobacterium.(ii)(): early latent, that is, individuals recently infected (less than two years) but not infectious.(iii)(): infected, that is, individuals who have active tuberculosis and are infectious.(iv)(): persistent latent, that is, individuals who were infected and remain latent.(v)(): recovered, that is, individuals who were previously infected and treated.
We are interested in studying the dynamics of Koch bacillus spread at several region in a discrete-time model. In addition, in order to find the best strategy to reduce the number of susceptible, infected who have active MTB or recently and persistent infected latent, and to rise the number of recovered, we will use three control, namely, vaccination, travel-blocking, and treatment programs. The paper is organized as follows: in Section 2, the model. S, L1, I, L2, and R is described for a multiregion discrete model. In Sections 3 and 4, we present the optimal control problem for the proposed model where we give some results concerning the existence of the optimal controls and we characterize these optimal controls using Pontryagin’s maximum principle in discrete time; illustrative numerical simulations are given at the end. Finally, we conclude the paper in Section 5.
2. Formulation of the Multiregion Discrete Mathematical Model
We assume that there are p geographical regions (domains) of the domain studied (in Figure 1, is the region of oriental with p = 8).
Let , and let be the population of domain at time i, i.e., the number of individuals who are physically present in , both residents and travelers.
According to the disease transmission mechanism, the host population of is grouped into five epidemiological compartments, and let and be the number of individuals in the susceptible, early latent, infective, persistent latent, and removed compartments of at time i, respectively. In addition to the death and recruitment, there are population movements among those five epidemiological compartments from time unit i to time i + 1. We assume that the recruited individuals (by birth and immigration) are constant and enter the susceptible compartment. To make the model a little more realistic, but in order to work with a constant overall population, we suppose that birth and death occur with the same rate. In addition, we suppose that individuals who are out of their domain do not give birth, and so birth occurs in the home domain at a per capita rate . And death takes place anywhere with a per capita rate . After one unit time, the susceptible individuals may stay in the susceptible compartment, or get infected and move to the infectious compartment, or die. Disease transmission is assumed to occur between individuals present in a given domain . The individuals in the infectious compartment may continue to be infected, or recover and be transferred to the recovered compartment, or die. Moreover, here we have assumed that there is no mortality due to infection. For instance, since the average life span of Moroccan people is 76.45 years, we have [30]. Individuals in the early latent compartment can progress either to active disease with rate or to a persistent latent infection with rate , following the approach in [31]. Parameter reflects that only 5% of infected individuals will ever develop active TB [15]. We choose d such that the progression rate from early infections to active disease is , which roughly approximates the data by Styblo [32], describing the proportions of disease development after conversion. For the rates of reactivation we adopt for untreated latent infections [33] and for those who have undergone a therapeutic intervention. As in Gomes et al. [34], the partial susceptibility factor affecting the rate of exogenous reinfection of untreated individuals, s, is fixed at 0.25, in accordance with the highest estimates of protection conferred by BCG vaccination (see [35]). In treated patients, this factor becomes , for which several exploratory values are adopted. Treatment of different infection stages is implemented at specific rates: applies to active TB and represents the rate of recovery (typically as a result of treatment, though here it also accounts for the infrequent natural recovery); and apply, respectively, to the latent classes and as the rates at which chemotherapy or a postexposure vaccine is applied. The rate is fixed at , corresponding to an average duration of infectiousness of 6 months, while and are considered at different exploratory values (see [15]). The values of the model parameters presented in control system (2) are given in Table 1. The disease transmission in a given domain at time i is modeled using standard incidence, given bywhere the disease transmission coefficient is the proportion of adequate contacts in domain between a susceptible from and an infective from another domain . The multiregion discrete-time TB model associated with is written as follows:where is the population size corresponding to the domain at time i. It is clear that the population size remains constant for all ; in fact
3. The Model with Vaccination and Travel Blocking
3.1. Presentation of the Model
Let and let us denote by the set of indexes of domain at high risk; we introduce a control variable that characterizes the travel-blocking operation, in order to restrict movements from domains to , where
We also introduce a control variable that characterizes the effectiveness of vaccination in the aforementioned model (2). Then, for a given region the model is given as follows:
3.1.1. An Optimal Control Approach
We are interested in controlling the population of region . The problem is to minimize the objective functional given bywhere Bk > 0, Cj > 0, > 0, and > 0 are the weight constants of controls, the infected, and the removed in region Ωj, respectively. where and . Our goal is to minimize the infected and the cost of applying controls and increase the number of recovered in .
In other words, we are seeking optimal control and such thatwhere and are the control sets defined by . .
By using Pontryagin’s maximum principle [37], we derive necessary conditions for our optimal controls. For this purpose, we defined the Hamiltonian as
Theorem 1. Given an optimal control and solutions and there exist . The adjoint variables satisfy the following equations:where ; and are the transversality conditions.
In addition, we have
Proof. Using Pontryagin’s maximum principle [37, 38], we obtain the following adjoint equations:with .
To obtain the optimality conditions, we take the variation with respect to controls and and set it equal to zero:and then we obtain the optimal control pairBy the bounds in and of the controls, it is easy to obtain and in the following form:
3.1.2. Numerical Results
We now present numerical simulations associated with the aforementioned optimal control problem. We write a code in and simulate our results using different data. The optimality systems are solved based on an iterative discrete scheme that converges following an appropriate test similar to the one related to the forward-backward sweep method (FBSM). The state system with an initial guess is solved forward in time and then the adjoint system is solved backward in time because of the transversality conditions. Afterwards, we update the optimal control values using the values of state and costate variables obtained at the previous steps. Finally, we execute the previous steps till a tolerance criterion is reached.
The multiregion TB epidemic model we suggest here is applicable for any number of regions p. For example, in order to show the importance of our work, we choose , i.e., we consider three regions (domains) , and with different parameters based on [12] and given in Table 1. Other epidemiological and numerical parameters are presented in Tables 1 and 2, respectively. The initial value of each category, ; and , is given in Table 2, which is based on [39] (Figures 2–13).
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(e)
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(e)
(a)
(b)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
Since the BCG vaccine in Morocco is mandatory, we considered the BCG vaccination control in the three regions given by (11). We are also interested in controlling and , by assuming that is at high risk of infection (see Figure 2); then, . We introduce the second control given by (10), in order to reduce entry of infected people from region .
In Figure 2, we can observe that in the absence of control and in the presence of an epidemic which spreads in three regions characterized by different parameters, the number of susceptible individuals varies from , and as initial conditions to 5088, 12094, and 19230, respectively, in the three regions; on the other hand, when the BCG vaccination controls , given by (11), are introduced in the three regions and we blocked the movement of infected individuals to the Ω2 region by introducing the controls given by (10), we observe that the number of susceptible subjects decreases from the initial state in the three regions to 703, 545, and 3556, respectively (see Figure 4).
To understand well, we compared each state variable in the three regions without control (see Figure 3), and with the introduction of vaccination and blocking the movement (see Figure 5), we can see clearly in Figure 5 that the curves of the susceptible decrease in all regions especially in second region after 2 years. On the other hand, the recovered curves increase markedly after the implementation of the first strategy which proves the impact of these controls.
In Figure 4, we can clearly see that after two years, there is effectiveness of vaccination on all the regions and especially on the second region controlled by the blockage of displacement, which is seen through the decrease in the curves related to the number of susceptible. At the same time, we see that we have a rapid increase in the number of recovered, especially in the second region.
In Figure 6(a), it can be seen that the higher the value of the vaccination control (), the more the number of susceptible is reduced, and this is clearly seen in Figure 5. On the other hand, in Figure 6(b), we find that when the travel-blocking control takes small values, i.e., when the contact between susceptible and infected individuals is too low, we obtain good results. This can be seen by comparing regions one and three for the number of infected (see Figure 5).
4. The Model with Travel Blocking, Vaccination, and Treatment
4.1. Presentation of the Model
We are now going to introduce a third control variable that characterizes the treatment of the infected in order to reduce them.
Then, for a given region , the model is given as follows:
4.1.1. An Optimal Control Approach
We are interested in controlling the population of region . The problem is to minimize the objective functional given by where Bk > 0, Cj > 0, Dj > 0, , and are the weight constants of controls, the infected, and the removed in region , , and .
Our goal is to minimize the infected and the cost of applying controls and increase the number of recovered in .
In other words, we are seeking optimal control and and such thatwhere , and are the control sets defined by
By using Pontryagin’s maximum principle [37], we derive necessary conditions for our optimal controls. For this purpose, we defined the Hamiltonian as
Theorem 2. Given optimal controls and solutions and , there exist . The adjoint variables satisfy the following equations:where ; and are the transversality conditions.
In addition, we have
Proof. Using Pontryagin’s maximum principle [37], we obtain the following adjoint equations:with .
To obtain the optimality conditions, we take the variation with respect to controls , and and set it equal to zero:and then we obtain the optimal controlsBy the bounds in , and of the controls, it is easy to obtain , and in the following form:
4.1.2. Numerical Results
In this section, we will also introduce the third control (a very important control); it is the treatment control, and of course we force all the infected to take anti-tuberculosis drugs; it is observed, from these subfigures, that using this control strategy, there is a significant decrease in the number of infected individuals treated compared to those untreated, especially in the second region (Figure 7).
In Figure 7, we can clearly see that after first year of disease spread, there is effectiveness of treatment in all the regions and especially in the second region controlled by the blockage of displacement. We observe that the number of infected individuals increases from , and as initial conditions to 347, 691, and 1036 after two years, respectively, in the three regions, in the absence of controls; on the other hand, when we introduced the treatment control in the three regions, the number of infected individuals decreases to 124, 114, and 269 after 2 years, respectively. This is a direct result of the introduction of the treatment in the three regions, and for that we notice in the second region that the red curve reaches its maximum before two years.
In Figure 8, we notice that the number of susceptible decreases from , and as initial conditions to 215, 296, and 833 in the three regions, respectively; it is because of the control effect especially in the second high-risk region. Similarly, the number of individuals cured in the second region increases from 166 people without control to 15267 with the effect of vaccination and blockage of movement and finally treatment. The impact of the introduction of the treatment has given good results with the decrease in the number of infected individuals.
In Figure 9, we present the populations of susceptible in each region in three cases without control (blue curves) and with control in two strategies: with vaccination and travel blocking (red curves) and with vaccination, travel blocking, and treatment (green curves); we see that the number of susceptible decreases quickly after implementing the two strategies.
In Figure 10, we present the populations of early latent in each region in three cases without control (blue curves) and with control in two strategies: with vaccination and travel blocking (red curves) and with vaccination, travel blocking, and treatment (green curves); we see that the number of early latent decreases, and of course, the best strategy is the second.
In Figure 11, when we compare the infected in each region without control and with control in two strategies, we clearly see that the green curves decrease rapidly to zero. This is due to the effect of the treatment on the infected individuals.
In Figure 12, the comparison of persistent latent shows us the efficiency of the introduction of controls; it is therefore seen that the infected with long controlled duration does not increase like those uncontrolled.
In Figure 13, by comparing the curves corresponding to the recovered for each region, we can see that there is a significant difference between the individuals recovered with and without control; this is due to the impact of the control strategies used, and for the second region, we have a large increase due to the blocking of movement towards this last region and the majority of those infected are treated.
In conclusion, when we compare the curves obtained in the numerical simulations, we see that the strategy which includes vaccination, treatment, and blocking of trips is more effective than the first strategy (vaccination and blocking of trips).
5. Conclusion
Tuberculosis represents a considerable burden on health services all over the world because the process is long and very expensive, so it is impossible to immunize the entire population, especially for developing countries. In this case, it is necessary to target the regions most likely to catch tuberculosis, i.e., “regions of fragility”; in this way, we will reduce the cost and the effort deployed by public authorities. This is the reason why we chose to work on multiregion models. This paper analyzes two different optimal control issues for a TB multiregion discrete-time model. In the first part, we introduce two controls that characterize the effectiveness of vaccination and blocking trips, where we vaccinated some of the susceptible people in all regions; except that in the high-risk region, we vaccinated as many susceptible people as possible, in order to reduce the number of susceptible and increase the number of recovered; at the same time, the displacement to high-risk areas was prohibited by controlling contact between susceptible people and infected individuals at high risk of infection. The numerical results associated with this first part assured us of the effectiveness of the vaccination and the impact of the travel blocking on reducing the spread of the disease from one region to another. We extended the first part to a second part by adding a processing control for the infected individuals, and above all, we concentrated the efforts on the high-risk region and we treated most of the individuals in this region. Our numerical results show that this strategy (the strategy that combines vaccination, treatment, and blocking trips) is more effective than the first strategy (see Figure 11). In conclusion, we have shown the interest of each case studied by a comparison between the numerical simulations provided. We have demonstrated the importance and effectiveness of our approach in controlling the spread of infection. Control programs based on these strategies can effectively reduce the number of infected cases and increase the number of individuals removed from a specific region. As a perspective of this work and in order to make our model more realistic, we will introduce the notion of delay to model the latency period.
Data Availability
The data used to support the findings of this study can be obtained from the Moroccan Ministry of Health.
Conflicts of Interest
The authors declare that they have no conflicts of interest.