The authors have declared that no competing interests exist.
There is an emerging consensus that achieving global tuberculosis control targets will require more proactive case finding approaches than are currently used in high-incidence settings. Household contact tracing (HHCT), for which households of newly diagnosed cases are actively screened for additional infected individuals is a potentially efficient approach to finding new cases of tuberculosis, however randomized trials assessing the population-level effects of such interventions in settings with sustained community transmission have shown mixed results. One potential explanation for this is that household transmission is responsible for a variable proportion of population-level tuberculosis burden between settings. For example, transmission is more likely to occur in households in settings with a lower tuberculosis burden and where individuals mix preferentially in local areas, compared with settings with higher disease burden and more dispersed mixing. To better understand the relationship between endemic incidence levels, social mixing, and the impact of HHCT, we developed a spatially explicit model of coupled household and community transmission. We found that the impact of HHCT was robust across settings of varied incidence and community contact patterns. In contrast, we found that the effects of community contact tracing interventions were sensitive to community contact patterns. Our results suggest that the protective benefits of HHCT are robust and the benefits of this intervention are likely to be maintained across epidemiological settings.
Screening household members of newly detected tuberculosis cases is an efficient method for finding previously undiagnosed cases in high-burden settings. Despite the intuitive appeal of this approach, randomized trials examining the population-level effects of these interventions in settings with sustained community transmission have shown mixed results. One explanation for these inconclusive findings is that household transmission is responsible for a varying proportion of overall tuberculosis burden between locations, with the impact of household transmission being a function of both the overall incidence and the relative intensity of disease-transmitting contacts in the community and the household. In this manuscript, we use an individual-based network model to explore how local incidence levels and patterns of community contact impact the effectiveness of household-based approaches for interrupting tuberculosis transmission. Our analyses suggest that protective benefits of household-based interventions are maintained across a wide range of epidemiological settings. Our findings provide evidence for the robustness of household-based interventions and suggest that variable results from trials may be primarily due to implementation challenges rather than inherent limitations of these interventions.
Despite recent progress in the development of new tuberculosis (TB) diagnostics [1], drugs [2], and vaccines [3], the decline in TB incidence remains far too slow to meet global targets for TB control. We need more effective case detection strategies to extract the maximum benefit available from existing tools. There is an emerging consensus that passive TB case-finding (i.e. waiting for individuals with symptoms consistent with TB to seek medical care), is insufficient and must be augmented with more active approaches that allow cases to be detected and treated as early as possible [4].
Household contact tracing (HHCT) has been advocated as an efficient approach to TB treatment and prevention [5, 6], both because household contacts of known TB cases are at high risk of transmission, and because they are relatively easy to identify. These interventions typically begin when a household index case is diagnosed with TB after presenting for care, i.e. ascertained by passive case-finding. Household contacts may then be screened and treated for active TB and latent TB infection (LTBI). HHCT is routinely applied in high-income, low-incidence settings. However, evidence regarding its efficacy in reducing population-level TB risk in higher-incidence settings is mixed. While HHCT has been demonstrated to improve the yield of TB case finding over passive detection [7], conclusions about the impact of HHCT on population-level TB incidence from cluster-randomized trials are varied, with some demonstrating improved TB control at the community level [8] while others failing to show robust population-level effects [9].
One potential explanation for these heterogeneous outcomes is the relative concentration of TB transmission in households and the community. Intuitively, where the fraction of transmission concentrated within households contributes to a greater share of overall TB burden (e.g., in low-incidence settings), we might expect HHCT to be more effective than in settings where there is a greater community burden of TB and household transmission represents a smaller fraction of the total (e.g., high-incidence settings). The average number of contacts and physical distance between contacts (i.e., social mixing patterns) can also affect the distribution of TB transmission in households and the community. For example, longer distances between community contacts will lead to fewer shared contacts between household members. This will likely result in a smaller fraction of co-prevalent household cases and HHCT therefore being less effective. In order to better understand the complex relationship between endemic incidence and social mixing patterns on the impact of HHCT, we developed a spatially-explicit network model of coupled household and community TB transmission. We then used this model to characterize the relationship between endemic incidence levels and community contact patterns on the impact of HHCT.
In this section, we outline the components of our individually-based network transmission model, which extends a previously published individual-based TB model of coupled household and community transmission [10] by introducing spatially structured community contact networks adapted from [11]. Specifically, we extended the Gaussian community contact network from [11] to include household transmission. Further, we extended the natural history model from [10] to include TB transmission across the network and added in separate intervention (e.g., treatment) compartments.
Utilizing a network representation of community contact allows us to assign discrete contacts at the individual level, which will then be used for the contact tracing interventions described below. Contacts represent individuals who interact with each other and are therefore potential tuberculosis transmission links. They may be close contacts (e.g., within households) or casual contacts (e.g., within the community). Our model represents a population consisting of 100,000 individuals divided evenly into 20,000 households (5 individuals per household). Each household is placed uniformly on a 2-dimensional grid. The spatial extent of the grid is implicitly set by the network density, which we fixed to equal 1. Individuals within households are all connected to each other, while community contacts are formed according to a Gaussian connectivity kernel in which the probability of connection between individuals varies as a function of physical distance between their households. On the network level, this kernel controls both the average number of contacts each individual has (average degree) as well as the average physical distance between community contacts (average connection radius). For instance, in networks with a higher average degree, individuals have more potential transmission links from the community. Furthermore, networks with a lower average connection radius have contacts (potential transmission links) that are closer together and may lead to transmission that is more clustered and less widely dispersed throughout the community. For simplicity, we assumed that household and community contacts are static for the duration of the simulation. See S1 Text for more details on the construction of community contact networks.
To explore the robustness of HHCT to different intensities and configurations of community contact (i.e., to examine the impact across different community contact settings), we conducted simulations across a wide range of connectivity kernel parameters (see Table 1). To account for random variation in contact network structure we generated 10 network realizations for each input parameter set. We also calculated the global clustering coefficient, C, (S1 Text) for each network [12] among community contacts, which measures the degree to which individuals have overlapping contacts. Specifically, for a given generated network, the clustering coefficient was calculated by finding all triads (groups of 3 households) with ≥2 connections (i.e., households were deemed ‘connected’ if there was at least one community contact between them). Then dividing all triads with 3 ties by the total number triads [12]. See Fig 1 for a schematic of household and community network structure and example networks with different clustering coefficients and S1 Fig for features of generated networks.


Network structure.
(A) Schematic of Network Structure (top): All individuals are fully connected within their households. Individuals form community contacts based on a Gaussian (normally distributed) connectivity kernel [11]. Networks consist of 100,000 individuals divided evenly into 20,000 households. (B & C): Example networks consisting of 1000 nodes (i.e., households) and similar average degrees with minimal (B; bottom left) and high levels (C; bottom right) of community clustering.

| Parameter | Description | Range | Source |
|---|---|---|---|
| n | Average degree i.e., the total number of community contacts | 25 to 200 | A large range is set to explore variation in community contact. Actual range of generated networks was slightly different. |
| σ | Average connection radius | 0.5 to 5 | Range set to obtain predetermined average degree distribution. Actual range of generated networks was slightly different. |
| ρ | Network density: i.e., density of households in grid | Fixed at 1 | Assumption |
The natural history of TB is characterized by a multi-stage latency period, in which recently infected individuals are at highest risk of progression to active disease, but may also enter a state of long-term latency from which they may progress to disease many years after infection, with most infected individuals (∼80%) never developing active TB. See S2 Fig for the distribution of latent progression times and the fraction progressing to active TB across different parameterizations of the model. Because early detection is key to the success of HHCT, our model must provide an adequate representation of these transitions. To accomplish this, we adapted a previously published model (from [10]) which includes 5 key disease states. Individuals are born as uninfected and susceptible (S). Upon infection, individuals enter an early latent period (EL) in which the annual rate of progression to infectious active TB (I) decreases over five years before entering the late or long-term latent (LL) state from which a small subset of individuals may progress to active TB. Infectious individuals may then enter the recovered state (R) spontaneously or as a function of treatment. Finally, individuals may exit the model by death as a function of TB mortality or the background mortality rate. When estimating the rate at which individuals become infected i.e., the force of infection (FOI) on individual i in household j at time t (λij(t)), we account for variable intensities of household and community contact as follows:



TB transmission model.
Schematic of TB Transmission Model. Where S is susceptible; EL is early latent (divided in sub-states to represent the decreased rate of progression to I with time since infection); LL is late latent; I is infectious active TB; R is recovered; T is treatment; IPT represents individuals who are currently taking preventive therapy. Individuals transitioning to the treated states (T, and IPT) are represented by dot-dashed lines. The rate at which individuals transition to the IPT states depends on the number of contacts. Births and deaths are represented by dotted lines.
Where possible, we obtained point estimates and uncertainty intervals for natural history parameters from published sources, including reviews of historical studies conducted prior to the advent of TB chemotherapy, which provide information on rates of death and spontaneous recovery from untreated active TB [15]. We also consulted systematic reviews e.g., for the rate of progression from LL to active TB [16], and modeling analyses e.g., for the rate of progression from EL to active TB [17]. Household and community transmission rates were calibrated to reproduce a range of incidence levels (20-400 cases per 100,000 person-years). We chose this range to focus on epidemiological settings in which there was sustained endemic community transmission. Furthermore, in the most countries, TB incidence levels are within this range. For instance, Bulgaria has an estimated annual incidence of 22 per 100,000 person-years, Peru has an estimated annual incidence of 123 per 100,000 person-years, and Eswatini has an estimated annual incidence of 329 per 100,000 person-years [18]. Following results indicating a higher per-contact rate of transmission from household vs. community contacts [19], in all simulations, we constrained the per-contact community transmission rate to be less than or at most equal to the household transmission rate, i.e. βC ≤ βHH. However, since individuals are likely to have more community than household contacts, this does not exclude the possibility that an infectious individual will cause substantially more infections in the community than their household. Because we were primarily interested in examining how the relative proportions of TB transmission occurring in households and the community alters the impact of interventions, we allowed for a large range of household and community transmission rate values and constrained other key drivers of incidence e.g., the rate at which individuals with active TB seek care. See Table 2 for a full list of transmission model parameters and their sources.

| Parameter | Description (units) | Range | Source/Explanation |
|---|---|---|---|
| Natural History | |||
| βHH | Household transmission parameter | 0 to 1.25 | Upper bound is set empirically based on which βHH value when equal to βuC results in the target incidence levels (20 to 400 cases per 100,000 person-years) in a sufficient number of simulations. |
| βuC | Unscaled community transmission parameter | 0 to βHH (i.e., βuC ≤ βHH) | Derived by multiplying βHH by a scaling factor between 0 and 1. βuC is then divided by the average degree of community contacts (n − 4; see Table 1) to obtain βC which is used in the FOI. This parameter also has strong impact on overall TB incidence. (See Eq 1) |
| ω | Amount of immunity conferred by current state (%) | For S = 0%; EL, LL, and R = 80%; I, T, IPT = 100% | [13, 14] |
| θ | Life expectancy (years) | 72.38 | Global average life expectancy. [20]. The mortality rate used in the model is the inverse of the life expectancy |
| ϵ | Early latency progression rates (/yr) | A value between 0.0817 to 0.0905 is sampled and then multiplied by (1, 0.41, 0.13, 0.086, 0.028) to derive the rate of progression for each EL sub-state | [17] |
| τ | Late latency progression (/yr) | 0.0005 | [16] |
| γ | Recovery rate (/yr) | 0.09 to 0.15 | [10, 15] |
| κ | Active TB mortality rate (/yr) | 0.05 to 0.4 | [10, 15] |
| TB Screening and Treatment | |||
| cdr | The rate at which individuals with active TB self-present for care (/yr) | 1 | Individuals are detected an average of 1 year after progressing to active TB. We assumed: and used global numbers from [21] to obtain a disease duration of ∼ 1 year. This is consistent with other analyses e.g., [22, 23]. |
| txd | Treatment duration (months) | 6 | [24] |
| iptd | Preventive therapy duration (months) | 6 | [25] |
In all simulations, household index cases were ascertained via passive case finding [24, 26], representing current practice in most TB endemic settings [27]. We assumed that on average, it takes one year for an individual with active TB to be detected by this mechanism. Individuals with active TB who are found through passive case finding are given treatment, assumed to no longer be infectious, and eventually recover. In our intervention simulations, passive detection will then trigger the active case finding (ACF) interventions outlined below.
For each case discovered through passive case finding, all 4 household contacts of the index case were screened (during the same one-month time step) for active TB and LTBI. This scenario was adapted from [10, 28]. Contacts found to have active TB were placed on treatment, while those with LTBI (EL or LL) were given preventive therapy. While receiving preventive therapy, we assumed that individuals cannot become reinfected or progress to active TB; at the end of the treatment period, individuals enter the R state.
To determine whether any effects of HHCT could be attributed specifically to focusing ACF at the household level, we simulated two additional interventions in which the same number of individuals were screened per detected case: (1) community contact tracing (community CT) in which 4 community contacts of the index case were screened and (2) community-wide ACF, in which 4 randomly selected individuals were screened, following [10, 28]. See S1 Text for additional details of these alternative ACF interventions.
To account for parameter uncertainty and explore parameter values of interest (e.g. βC), we ran the model with 9,000 parameter sets obtained using Latin Hypercube Sampling [29] from predefined ranges based either on published estimates or ranges set based on our interpretation of the literature. At the beginning of each simulation run, we selected a network at random and ran the transmission model 5 times using different random number seeds to allow stochastic variation across realizations. We ran the model with passive case finding only until it reached endemic equilibrium. Next, we implemented a simulated ACF trial for 5 years representing a plausible time horizon to evaluate the impact of screening interventions. We ran all interventions and a passive case finding only scenario with each parameter set and random number seed combination. For more details see S1 Text. We generated networks for the model using R version 3.4.1 [30]. Our individual-based network model was run on Matlab version R2019a [31]. We ran our simulations in parallel on a high performance computing cluster and each individual simulation run took between 2 and 10 minutes and used ∼4.5 gigabytes. Model code can be accessed at https://github.com/jhavumaki/network_tb_ibm.
To assess the impact of HHCT, we calculated the 5-year cumulative incidence rate at the end of the ACF trial period. Incidence was defined as the number of new active TB cases per 100,000 person-years over the previous 60 one-month time steps. We then calculated rate ratios (RRs) to quantify the protective benefit of HHCT for each parameter set and random number seed combination. Specifically, we divided the incidence rate (at the end of the ACF trial period after HHCT) by the passive case finding only scenario incidence rate at the end of the simulation (see S3 Fig for overview of the simulation workflow). Both comparator interventions were also compared to the passive case finding only scenario.
RRs were compared across all parameter sets and also, within strata of network parameters (average degree and average connection radius), and initial incidence-levels (i.e., immediately before the active screening interventions were implemented) to measure the protective benefit conferred by HHCT in different transmission and community contact (i.e., network) settings.
To ensure that our model behaved as expected and to evaluate the impact of variation in community contact patterns on transmission dynamics, we examined all model runs from the passive case finding only scenario (without restricting incidence to our target range). The proportion of individuals with LTBI was approximately as expected for different TB prevalence levels (S4 Fig). Additionally, a greater proportion of infections were caused by community transmission in higher incidence settings compared with lower incidence settings (S5 Fig). Furthermore on clustered networks, community transmission contributed more to overall TB burden than on less clustered networks. Despite the fact that higher TB incidence levels were driven by more community transmission and that community transmission contributed more to overall TB burden among clustered networks, we observed that as the network clustering coefficient increased, incidence levels were generally lower (S6 Fig). This is explained by the fact that more overlapping contacts (on clustered networks) likely created local contact saturation and depletion of the pool of susceptible individuals [32]. The remainder of the analysis was conducted among simulations within our target range of incidence levels i.e., between 20 to 400 cases per 100,000 person-years.
Overall, the protective benefits conferred by HHCT were robust across all settings (i.e., varying incidence level, average degree, average connection radius, and clustering coefficient) with median RRs equaling ∼0.7. See Fig 3 and S13 Fig and S1–S4 Tables for more details.


Effectiveness of HHCT across different settings.
Fitted splines representing relationship between all RRs and (A) the incidence rate immediately before HHCT (per 100,000 person-years) (top left), (B) the average degree (top right), (C) the community clustering coefficient (bottom left), and (D) the average connection radius (bottom right). Lines are splines calculated using the LOESS method in R [33]. Among model runs with incidence rates between 20 and 400 cases per 100,000 person-years. Shaded regions represent 95% confidence intervals.
To ensure that the impact of HHCT was not sensitive to input parameter combinations and endemic incidence levels (e.g., to account for effect modification), we plotted the density of RR values within all pairwise distributions of incidence and network parameter strata (e.g., incidence level by average degree, or average degree by average connection radius). This enabled us to determine whether e.g., average degree affects the impact of HHCT for a given incidence level (S8–S12 Figs). Overall, there did not appear to be any additional emergent trends. See for details.
Contrary to the robustness of HHCT across settings, community CT was sensitive to average degree, average connection radius, and clustering coefficient with the effect of incidence level being similar to HHCT. For example, as average degree increased, community CT became less effective (Fig 4). Notably although community CT was sensitive to network parameters, standard deviations associated with the RRs were greater than the differences between strata. See S14 Fig and S5–S8 Table for all results.


Effectiveness of community CT varying average degree.
Fitted splines representing relationship between all community CT RRs and the average degree. Lines are splines calculated using the LOESS method in R [33]. Among model runs with incidence rates between 20 and 400 cases per 100,000 person-years. Shaded regions represent 95% confidence intervals.
With respect to community-wide ACF, its effects were close to null and it was most sensitive to incidence level. See S15 Fig and S9–S12 Tables for results.
We made additional comparisons between interventions to further assess the performance of HHCT. First, we examined the infectious period duration for each intervention and found that the median infectious period was ∼1 month shorter for HHCT compared with community CT, community-wide ACF and passive surveillance only scenarios (S16 Fig). Next, we estimated the number of secondary cases averted among household contacts due to preventive therapy. As expected, the modeled HHCT intervention prevented more cases than the modeled community CT which in turn, prevented more cases than the modeled community-wide ACF (S17 Fig). We then examined the relative number of preventive therapy administrations compared with treatment administrations for each screening scenario. HHCT led to substantially more preventive therapy administrations compared with the other screening scenarios indicating that it may be a more effective approach to identify and treat individuals with recent infection at highest risk of progression (S18 Fig). Finally, we examined the prevalence of LTBI and active TB among household and community contacts. We found that HHCT resulted in lower prevalence of active TB among both household and community contacts, but the prevalence of LTBI was similar between interventions (S19 Fig). The prevalence of LTBI was similar because it represents a cumulative lifetime exposure status and therefore changes substantially more slowly than the prevalence of active TB. Therefore, because we only ran the intervention for 5 years, we did not see major changes in the LTBI prevalence.
We conducted additional sensitivity analyses to investigate the extent to which different sources of variability affected the estimated screening scenario impacts. To do this we examined the relative contributions of stochasticity in the disease transmission model, using results from multiple network realizations, and variability in parameter inputs. We constructed a hierarchical regression model to examine the relative contributions of these inputs to variation in estimated RRs across simulations (S13 Table). We found that stochasticity and network realization had only modest effects on RR values, contributing to ∼5% and ∼3.5% of variability in RRs, respectively. Finally, the overall model R2 was ∼85.6% indicating that the parameterization of the model contributed substantially to overall changes in RRs.
Next, we incorporated imported TB cases into the force of infection to understand if a low background level of risk for the entire population might affect the projected impact of interventions. Results from this analysis revealed that although the impacts of all screening interventions were reduced, the main conclusions from our analysis did not change (S22–S24 Figs).
Finally, to explore whether variation in the implementation of HHCT might lead to less robust results, we conducted an additional sensitivity analysis varying the coverage of HHCT. We found that in higher incidence settings (>100 cases per 100,000 person-years), lower HHCT coverage levels resulted in much smaller population-level benefits (S21 Fig).
Our analyses suggest that the protective benefits of HHCT are likely to be robust across diverse settings characterized by variation in incidence level and community contact patterns (Fig 3). As compared to the consistent effectiveness of HHCT, community CT appeared to be more sensitive to different network parameters, suggesting that its utility is more limited to specific scenarios. For instance, as average degree increased, the protective benefits of community CT decreased (S14 Fig). Overall, this suggests that HHCT is robust to different relative proportions of TB transmission occurring in households and the community while community CT is not.
Despite the relative effectiveness of HHCT in reducing community incidence, community transmission was the dominant mode of infection in all but the lowest-incidence settings (see Table 2 and S5 and S7 Figs). This can be partially explained by contact saturation within households [32]. Additionally, transmission causing more infections in higher-incidence settings is consistent with molecular epidemiology studies which have revealed that, in these settings, the majority of co-prevalent TB cases within households are genetically discordant [34, 35]. Thus, both the transmission rates and the number of contacts are important drivers of TB transmission. One might expect HHCT to appear less effective in settings where most TB transmission occurs in the community. However, our results challenge this intuition, as we found that the benefits of HHCT were maintained despite differences in the proportion of transmission events occurring in the community (S20 Fig). This finding is consistent with a recent model-based analysis indicating that the effects of HHCT are robust to observed variation in local community-level exposure [28].
Our finding of the robustness of the projected benefits of HHCT across settings suggests that differences in community and household forces of infection probably do not explain mixed results of randomized control trials of HHCT to reduce TB incidence [8, 9]. Others have suggested that HHCT may have limited effects in high transmission settings where evidence suggests that a minority of transmission events occur within the home [34–36]. Our results reveal that HHCT has a substantial impact even when within-home transmission does not account for the majority of infection events because household contacts represent a higher-risk population given the likelihood of clustering of exposure risk. This suggests that other explanations may be needed for the mixed results of trials of HHCT, including the coverage and quality of case finding activities. Our sensitivity analysis examining the impacts of varying HHCT coverage confirms that this might be a potential explanation (S21 Fig)
Our approach has some limitations that should be considered in interpreting these findings. For the sake of parsimony, we did not include additional factors in our model that could have altered the risk of TB at the individual level. For example, TB can cluster due to shared risk factors (aside from transmission) like alcohol use and malnutrition [37]. Future models exploring the impact of individual-level heterogeneity on the robustness of HHCT will be helpful in clarifying the real-world utility of this approach. Next, our model does not include an interaction between the intensity of symptoms and the infectiousness of cases, as has been done in earlier models [38]. Including these dynamics would have likely affected the relative performance of all interventions in the same way so it would have not changed our overall conclusions regarding the efficacy of HHCT. At the same time, models including the time-dynamics of infectiousness will be extremely helpful for guiding the timing of the proposed interventions.
It is important to note that the implementation of contact patterns were simplified to reduce the computational demands of our model. First, we assumed static network connections. Dynamic network connections would have allowed infection to escape more easily from local clusters [39], dampening differences across more and less clustered networks. To ensure that we accounted for this, we conducted our analyses across different networks (i.e., with different clustering coefficients) and this likely accounts a wide range of distributions of TB. We also assumed that number of individuals per household is fixed (i.e., at 5 people). Although this is unrealistic, keeping the household size fixed allowed us examine how variation in community contact alone might impact the effectiveness of different screening interventions. Examining the impact of variation in two dimensions (i.e., on both the household and community levels), may lead to different effects. Finally, our births reshuffling scheme based on [39] is not realistic, however, it allows for the population size to be constant and prevents our model from placing susceptible individuals into areas of high endemic transmission.
These results suggest important future directions this research can take to understand the real-world applicability of HHCT in high-burden settings. Our model focused on the effects of heterogeneity across networks with varying average degrees and average connection radii but relatively low within-network variability in connectivity. Future work should consider the effects of heterogeneity within a single network in more detail to capture the impact of variable household sizes, variable average degrees, super-spreading individuals, spatial hotspots, and other causes of right-skewed community contact distributions. Adding within network heterogeneity will alter the local risk of TB. Although the effects of HHCT have been found to be robust to observed variation in local community-level exposure in one setting [28], it is unclear what a systematic exploration of different within network heterogeneities would reveal. Additionally, mechanisms that impact susceptibility and the risk of progression to active TB (e.g., HIV, malnutrition, diabetes) may also lead to more community-level heterogeneity.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39