Modeling the impact of COVID-19 on future tuberculosis burden.

Background: The ongoing COVID-19 pandemic has greatly disrupted our everyday life, forcing the adoption of non-pharmaceutical interventions in many countries and putting public health services and healthcare systems worldwide under stress. These circumstances are leading to unintended effects such as the increase in the burden of other diseases.
Methods: Here, using a data-driven epidemiological model for tuberculosis (TB) spreading, we describe the expected rise in TB incidence and mortality if COVID-associated changes in TB notification are sustained and attributable entirely to disrupted diagnosis and treatment adherence.
Results: Our calculations show that the reduction in diagnosis of new TB cases due to the COVID-19 pandemic could result in 228k (CI 187-276) excess deaths in India, 111k (CI 93-134) in Indonesia, 27k (CI 21-33) in Pakistan, and 12k (CI 9-18) in Kenya. Conclusions: We show that it is possible to reverse these excess deaths by increasing the pre-covid diagnosis capabilities from 15 to 50% for 2 to 4 years. This would prevent almost all TB-related excess mortality that could be caused by the COVID-19 pandemic if no additional preventative measures are introduced. Our work therefore provides guidelines for mitigating the impact of COVID-19 on tuberculosis epidemic in the years to come.

Introduction

Tuberculosis (TB) is an infectious disease caused by the bacterium Mycobacterium tuberculosis (M.tb.) that usually affects the lungs. It is a preventable but complex disease with a high global burden that requires early detection and long treatments. Despite the global effort to eradicate TB and recent decreases in its burden due to the implementation of strategies aimed at optimizing diagnosis and treatment[1,2], it remains one of the greatest threats to public health worldwide, being the deadliest single-agent persistent infectious disease nowadays. According to the 2021 Global TB Report by the World Health Organization (WHO)[3], ten million people developed TB and nearly 1.5 million people died because of TB infection in 2020, and for the first time in a decade, there is an increase in TB-caused deaths. In the last decades, the WHO has deployed a series of global strategies that have since been the backbone of the global fight against TB. In 1995, the Directed Observed Treatment Strategy was introduced, which significantly strengthened the capacity of national programs to diagnose and treat TB cases. Later, the Stop TB Strategy, announced in 2006, was the first of such plans at setting a TB elimination horizon, defined as a reduction of incidence levels under one case per million and year by 2050. A redefinition of the eradication goal took place in 2014 when the previous objective was moved forward to 2035 within the End TB Strategy.

If the elimination target set by the End TB strategy was already an ambitious goal[4], the emergence of the COVID-19 pandemic caused by the new coronavirus SARS-CoV-2 sheds significant concerns on whether these goals are still reachable. During the acute stages of the COVID-19 pandemic, economic and human resources were redirected to control and mitigate the emergency caused by the pandemic, which led to a great reduction in the diagnosis of new cases of other diseases, as already documented for cancer, or malaria[5,6]. Interventions such as long lockdowns and mobility restrictions have exacerbated shortages in resources otherwise destined for the care of patients suffering these, and other pathologies. Moreover, COVID-19 has greatly affected healthcare workers[7-9], thus creating additional pressure to healthcare systems.

TB diagnosis and patient care are no exceptions, as reported in previous literature[3,10,11]. As a primary and immediate effect of COVID-19 spreading onto TB transmission dynamics, a reduction in the case notification ratio has been observed during and immediately after lockdowns and periods of high COVID-19 incidence and saturation of healthcare facilities[3]. We hypothesize that this disruption alone will lead to a surge in TB burden in the next years, even before the more complex, and less predictable effects of the COVID-19 pandemic on TB management and transmission dynamics can be properly characterized. For example, drastic drops in laboratory capacity needed to support TB diagnosis are expected along with interruptions in the supply of drugs, which could result in shortages of medications and could delay the start of treatments until the supply chain is reestablished[12-14]. Moreover, as suggested by Cilloni et al.[15], even temporary stoppages might cause long-term increases in TB incidence and mortality, and a peak in TB burden is to be observed as a consequence of the healthcare system disruption.

In this work, we assess the impact of COVID-19 on the expected TB burden until the year 2035, which marks the target horizon of the End TB Strategy. Specifically, we incorporate the observed drop in TB diagnosis and treatment compliance rates caused by COVID-19 into a mathematical model that produces long-term forecasts of TB burden[16]. This allows us to: (1) quantify the effect of the COVID-19 stoppage in relation to a baseline scenario in which no pandemic happened, and (2) compute the effect that a rapid response to the uprising TB burden in the following years, in the form of a compensatory intervention aiming at boosting TB diagnosis rates as soon as the COVID-19 pandemic ends, has over long-term TB goals. Our results show that an effort focused on increasing TB diagnosis capabilities once the pandemic is over could revert the effect of the pandemic in the long term.

Methods

Model calibration and diagnosis rate

In this study, we have capitalized on the detailed M.tb. transmission model developed by Arregui et al.[16,17] (see Supplementary Methods and Supplementary Fig. 1). Conceptually, this model is an age-stratified compartmental model that describes TB dynamics within a whole, closed population, stratified into 15 age groups during periods of the order of several decades. The model is detailed enough to include demographic evolution and aging, along with heterogeneous contact patterns among age groups that have been adapted from empirical survey studies (see Supplementary Fig. 2).

Here, the model is calibrated to reproduce TB incidence and mortality rates in each country under study for the period 2000–2019, using the burden estimates provided by the WHO. The calibration process gives the diagnosis rate d(t) and the scaled infectiousness β(t), which are modeled as half-sigmoid-like curves, and, among other parameters, are country-specific. This allows the model to reproduce different epidemiological scenarios. Specifically, the diagnosis rate is defined as:

Therefore, the diagnosis rate is parameterized by two quantities (d0, d1), where d0 is the value at the beginning of the calibrating window (i.e., year 2000 in this study), and d1 defines its evolution, either increasing or decreasing with time depending on d1’s sign. In the case of a decreasing evolution, the diagnosis rate is bounded to be greater than zero, while in the case of increasing evolution the upper bound is (year−1)[16]. This latter upper bound corresponds to a minimum diagnosis period of 1 month, assuming that, with a conservative lower boundary, the main symptom of TB is a continuous cough lasting for 3 weeks, followed by a time to diagnose estimated to last at least 10 days[18]. For further details regarding the specific values of epidemiological parameters, calibration processes, and uncertainty estimates, the reader is referred to the original source[16] and the Supplementary Methods.

Once the model is calibrated, we use it to produce forecasts until 2035 under two different scenarios: the baseline scenario, namely, a scenario in which there is no COVID-19 pandemic and thus, no disruption in healthcare systems is introduced, and another one in which a disruption is introduced at the start of 2020 up to the end of 2021, which is the pandemic scenario. During the duration of the pandemic, the diagnosis rate drops according to the reduction observed in the notifications of TB cases that are reported by WHO online, in the last global TB report, and also by the Nikshay program in India. Therefore, the drops in diagnosis rate are country-specific. These drops in TB notifications are fitted to a bump-like asymmetric function, as described through dred(t) in Eq. (2). This function reproduces the real data and is then applied to the model-calibrated diagnosis rate to produce the diagnosis function under the pandemic scenario. The fitting procedure is a Levenberg-Marquardt Nonlinear Least-Squares using minpack.lm R’s package[19], where Eq. (2) is applied to the data normalized by the 2019 mean for each country:

The bump-like function described in Eq. (2) serves as a multiplier to the model-calibrated diagnosis rate, thus, being the diagnosis rate under the pandemic scenario D(t) = d(t) dred(t) with d(t) ≠ 1 only during the COVID-19 pandemic, and dred(t) = 1 otherwise. In Table 1 we report the fitted values of each parameter involved in the bump-like description of the TB notification drops.

Table 1

Fitted parameters for diagnosis reduction in selected countries.

Bump parameters
Country	Θ = {h, t1, t2, k1, k2}
Indonesia	Θ = {0.494,   0.391,   17.31,   6.226,   117.4}
Pakistan	Θ = {0.398,   0.257,   18.04,   3095.3,   826.8}
Kenya	Θ = {0.248,   0.859,   5.218,   0.905,   46.51}
India	Θ1 = {0.393,   0.272,   1.25,   134.7,   0.645}
	Θ2 = {0.594,   0.139,   1.069,   3.111,   114.12}

The fitting procedure is a Levenberg-Marquardt Nonlinear Least-Squares using minpack.lm R’s package[19], where Eq. (2) is applied to the WHO data normalized by the 2019 mean. For Indonesia, Pakistan and Kenya, one bump is enough for reproducing the data, whereas in India two separate bumps need to be concatenated, and are denoted here as Θ1 and Θ2 respectively. h, k1 and k2 are dimensionless quantities, whereas t1 and t2 have units of year−1.

Finally, during the period of recovery, interventions are aimed at compensating for the drop in diagnosis rates during the pandemic years. We modeled this by multiplying the diagnosis by a scale parameter, as already discussed. Once the recovery period is over, we assume that the diagnosis rate goes back to its original value as given by d(t) up to the end of the simulation.

Modeling decisions about the disruption

A conceptually deep limitation of this study that needs to be stressed is that we only describe the effects of COVID-19-induced reductions in TB diagnosis rates and treatment adherence as the main drivers of the interaction between both processes. Admittedly, the effects of the COVID-19 crisis on TB dynamics are more complex than what is described here, and will most likely include alterations in transmission dynamics, effects mediated by economic impact, and long-term damages to health care quality standards beyond diagnosis rates; all these being aspects that lie out the scope of our study, mainly because the relevant data needed to describe the effects of them on TB dynamics are yet to be produced.

Although some of the non-pharmaceutical interventions adopted worldwide have proven their efficacy in reducing COVID-19 spreading[20], their effectiveness highly depends upon general public adherence and proper knowledge about the pandemic risks. Whereas some studies[21-24] show that the knowledge, attitude, and practice toward COVID-19 basic preventive strategies and conducts are in general positive, there is a great variation between communities and, for example, in India, between socioeconomic levels. Specifically, rural populations, as well as individuals with lower education, and unskilled occupations, are associated with lower scores of knowledge, attitude, and practice toward the basic preventive strategies against COVID-19, which would in turn be expected to contribute to halting TB transmission too[23]. This lack of adherence in the lower socioeconomic levels[25] suggests that it might be misleading to assume that the implementation of countermeasures induces a reduction in the TB force of infection. On the other hand, the changes in mobility due to lockdowns and other restrictions indicate that most of the interactions happen in residential areas (e.g., households) while these interventions are in place (see Supplementary Notes 1 and Supplementary Fig. 3). Admittedly, this could be at the root of some recent observations that report that the number of children diagnosed with TB has increased and that non-pharmaceutical public health interventions likely reduced influenza transmission, but have a lesser effect on M.tb. transmission during 2020[26-28].

To contextualize our findings in broader scenarios where changes in TB transmission—either toward enhanced or reduced spreading—are considered, we show the results of the basic burden outcomes, incidence, and mortality, in each country, for scenarios in which the transmission is either reduced or enhanced. We also considered an alternative scenario in which the treatment availability is higher than the one adopted in the main text, thus, exploring the effect of an overestimate of the disruption over the treatment (see Supplementary Notes 2 and Supplementary Figs. 4 and 5).

First-line treatment reduction

According to previous reports[3,15], first-line TB treatment completion has dropped effectively as a consequence of the COVID-19 pandemic, with interruptions in the supply of drugs that delay the start of the treatment in those cases in which the remaining medical capabilities have been enough to diagnose the disease. This inconvenience could not only worsen the expected treatment outcome for the patient but also drive secondary infections even in diagnosed patients if they are not able to quarantine until the treatment could be carried out. We modeled this situation in terms of the epidemiological model by including a fraction of under-treatment pulmonary TB individuals (T) in the expression of the force of infection (λ(t)). On the baseline scenario and without disruptions, those T individuals are not able to contribute to λ(t) as we assume that they are under control by the healthcare system, thus, being controlled and either under quarantine or, later on, medicated with TB drugs that greatly reduce their infectiousness. This means that, under normal circumstances, diagnosed individuals are expected not to be a risk for the rest of the population. However, when disruptions in the supply chain appear, a drop is observed in the first-line and second-line treatments completion[15], and then diagnosed individuals who are not able to either start the treatment or quarantine could become a risk. For this reason, we obtain an estimate of the fraction of T individuals that contribute to λ(t), , from Cilloni et al.[15] as:where η = 0.788. This value attempts to capture this kind of impact in countries like India and Kenya. It is based on expert opinion in the Stop TB Partnership and USAID about the side effects of the COVID-19 pandemic on TB treatment completion. We assume it to be a good proxy for the real value for the other countries included in this study.

Table 2

Cumulative number of averted deaths (in thousands) in 2035 with a post-COVID-19 intervention initiated in 2022 in each of the countries studied.

Number of averted deaths (in thousands) by 2035
Country	Diagnosis effort	T = 2 years		T = 3 years		T = 4 years
Indonesia	1.15	36	(31–43)	52	(44–61)	66	(57–78)
	1.30	67	(57–79)	93	(80–110)	118	(102–139)
	1.45	93	(79–109)	128	(110–150)	160	(139–187)
Kenya	1.15	9	(7–13)	13	(10–19)	17	(12–24)
	1.30	18	(13–25)	25	(18–35)	32	(23–45)
	1.45	25	(18–35)	35	(26–50)	45	(32–63)
India	1.15	121	(98–156)	173	(141–224)	223	(182–288)
	1.30	223	(178–294)	315	(254–418)	402	(326–532)
	1.45	309	(246–419)	432	(347–586)	547	(444–735)
Pakistan	1.15	14	(12–16)	20	(17–23)	27	(23–30)
	1.30	25	(22–29)	36	(31–41)	47	(41–53)
	1.45	35	(30–39)	49	(43–56)	64	(56–71)

The values in the table are computed by calculating the difference between the model forecast for mortality with the pandemic scenario and with non-pharmaceutical interventions of different intensities of diagnosis effort and duration of the recovery period. Values are the median of the outcome and figures in parentheses are the 95% CI of the model projections.