Monocyte Transcriptional Responses to Mycobacterium tuberculosis Associate with Resistance to Tuberculin Skin Test and Interferon Gamma Release Assay Conversion

Abstract

ABSTRACT

Heavy exposure to Mycobacterium tuberculosis, the etiologic agent of tuberculosis (TB) and among the top infectious killers worldwide, results in infection that is cleared, contained, or progresses to disease. Some heavily exposed tuberculosis contacts show no evidence of infection using the tuberculin skin test (TST) and interferon gamma release assay (IGRA); yet the mechanisms underlying this “resister” (RSTR) phenotype are unclear. To identify transcriptional responses that distinguish RSTR monocytes, we performed transcriptome sequencing (RNA-seq) on monocytes isolated from heavily exposed household contacts in Uganda and gold miners in South Africa after ex vivo M. tuberculosis infection. Gene set enrichment analysis (GSEA) revealed several gene pathways that were consistently enriched in response to M. tuberculosis among RSTR subjects compared to controls with positive TST/IGRA testing (latent TB infection [LTBI]) across Uganda and South Africa. The most significantly enriched gene set in which expression was increased in RSTR relative to LTBI M. tuberculosis-infected monocytes was the tumor necrosis factor alpha (TNF-α) signaling pathway whose core enrichment (leading edge) substantially overlapped across RSTR populations. These leading-edge genes included candidate resistance genes (ABCA1 and DUSP2) with significantly increased expression among Uganda RSTRs (false-discovery rate [FDR], <0.1). The distinct monocyte transcriptional response to M. tuberculosis among RSTR subjects, including increased expression of the TNF signaling pathway, highlights genes and inflammatory pathways that may mediate resistance to TST/IGRA conversion and provides therapeutic targets to enhance host restriction of M. tuberculosis intracellular infection.

IMPORTANCE After heavy M. tuberculosis exposure, the events that determine why some individuals resist TST/IGRA conversion are poorly defined. Enrichment of the TNF signaling gene set among RSTR monocytes from multiple distinct cohorts suggests an important role for the monocyte TNF response in determining this alternative immune outcome. These TNF responses to M. tuberculosis among RSTRs may contribute to antimicrobial programs that result in early clearance or the priming of alternative (gamma interferon-independent) cellular responses.

INTRODUCTION

Aerosol exposure to Mycobacterium tuberculosis results in three basic outcomes in humans: (i) failure of the bacterium to establish infection or its early clearance, (ii) infection within alveolar macrophages and other myeloid cells that is contained but not sterilized by an adequate host cellular response, or (iii) progressive pulmonary or disseminated infection that results in human disease and potentially ongoing transmission. The critical host-pathogen events that result in sterilizing immunity are largely unknown.

To better understand the immunologic mechanism of natural resistance to M. tuberculosis infection, we (1) and others (2–4) have profiled immune cell signatures from individuals who fail to convert the tuberculin skin test (TST) and interferon gamma (IFN-γ) release assay (IGRA) despite heavy exposure to M. tuberculosis, which we label the “resister” (RSTR) phenotype (5–7). We examined household contacts of pulmonary tuberculosis (TB) cases in Uganda who completed serial TST and IGRA testing over 8 to 10 years of follow-up and male gold miners in South Africa who have extreme occupational M. tuberculosis exposure (8–10). Resistance to TST/IGRA conversion in RSTRs may be explained by multiple mechanisms, including early clearance through innate immune responses (2, 11), IFN-γ-independent T cell responses (1, 12), or humoral responses (1, 3, 5, 6). However, we do not currently know whether RSTR subjects have distinct monocyte/macrophage responses upon M. tuberculosis infection and whether the cellular inflammatory response is heightened or diminished compared to latent TB infection (LTBI) subjects.

We recently compared RSTR and LTBI transcriptomes in whole blood and in monocytes to identify any cellular pathways that may distinguish these clinical phenotypes or provide clues into their unique immune responses (13). Across heterogeneous cohorts, RSTR individuals had strong enrichment of gene sets in monocyte metabolic pathways. However, these studies were performed in the absence of any stimulation of M. tuberculosis antigens or live infection. In the current study, we compare the monocyte transcriptional responses between RSTR and LTBI subjects following ex vivo M. tuberculosis infection. We discovered multiple pathways that distinguish these clinical phenotypes across study sites. Surprisingly, we detected increased expression of genes involved with tumor necrosis factor (TNF) signaling among RSTR monocytes, which may contribute to the unique immunologic outcomes to M. tuberculosis infection in these highly exposed cohorts.

RESULTS

Inflammatory and metabolic gene set enrichment distinguishes RSTR transcriptomes across populations.

To evaluate whether resistance to TST/IGRA conversion correlates with distinct monocyte transcriptional responses to M. tuberculosis infection, we examined two cohorts of HIV-negative subjects with heavy M. tuberculosis exposure in Uganda and South Africa and compared transcriptomes from those with persistently negative TST/IGRA (RSTR) versus LTBI controls with concordant positive tests. We enriched CD14⁺ monocytes from cryopreserved peripheral blood mononuclear cells (PBMCs) from each cohort (n = 49 RSTR and 52 LTBI subjects in Uganda; n = 20 RSTR and 29 LTBI subjects in South Africa) and infected them with M. tuberculosis (H37Rv) or medium for 6 h before transcriptome (RNA-seq) analysis (Fig. 1). No demographic or epidemiologic features were associated with RSTR or LTBI phenotypes in Uganda (see Table S1 in the supplemental material), whereas ancestry was associated with the RSTR phenotype in South Africa (Table S2). We stratified subsequent South Africa analyses according to the dominant ethnic group.

FIG 1

Uganda household contact and South Africa gold miner recruitment and monocyte RNA-seq analyses. Peripheral blood mononuclear cells from each cohort were thawed and rested in M-CSF overnight before CD14⁺ magnetic bead column isolation, plated to allow for adherence overnight, and then infected with M. tuberculosis (H37Rv) for 6 h (MOI, 1). RNA was isolated for library preparation and RNA sequencing. The indicated subjects were excluded due to failed RNA yield or quality and contamination in Uganda. In South Africa, additional samples were excluded due to TST/IGRA conversion or reversion (3 RSTR and 1 LTBI subjects) and based on ancestry related to presumed lower M. tuberculosis exposure (6 RSTR subjects). The experiment was performed twice (Exp A and Exp B) among partially overlapping subjects from Uganda, resulting in a final 101 unique subjects (49 RSTR and 52 LTBI subjects) and once from 49 South Africa gold miner subjects (20 RSTR and 29 LTBI subjects). QC/MDS, quality control/multidimensional scaling.

Using gene set enrichment analysis (GSEA) (14, 15), we identified multiple “hallmark” gene sets (16) that differentiated RSTR and LTBI monocytes across Uganda and South Africa sites. These enriched RSTR transcriptional signatures included metabolic and inflammatory gene sets among unstimulated monocytes (Fig. 2A), a finding that we reported previously (13). Upon M. tuberculosis stimulation, several gene expression patterns emerged (groups i to iv in Fig. 2B). First (group i), genes within the “TNFα signaling via NF-κB” (TNFα hallmark) and the “inflammatory response” hallmark gene sets had lower expression at baseline and higher expression following M. tuberculosis among RSTR than LTBI monocytes. Second (group ii), collective expression of the “oxidative phosphorylation” and “adipogenesis” gene sets remained higher among RSTR subjects regardless of stimulation. M. tuberculosis infection resulted in decreased expression of these metabolic gene sets in both RSTR and LTBI (Fig. S1), consistent with a shift away from oxidative phosphorylation and toward aerobic glycolysis (Warburg effect) induced by mycobacterial stimuli (17, 18), but the relative RSTR versus LTBI expression patterns suggest that RSTR monocytes are less susceptible to this shift. Third (group iii), collective expression of genes in the “IFNγ response” gene set was lower in RSTR (versus LTBI) regardless of M. tuberculosis stimulation. At least for the M. tuberculosis-infected condition, this finding is expected, considering our phenotypes are defined by an M. tuberculosis-specific IFN-γ response (i.e., IGRA). The “allograft rejection” and “IFNα response” gene sets also followed this pattern. Overlap of the IFN-α and IFN-γ leading-edge (LE) genes (Fig. S2), which is a subset of genes that contributes to the core enrichment of a gene set (15), suggests that the IFN-α enrichment is largely driven by IFN-γ signaling.

FIG 2

Gene set enrichment analysis identifies pathways that distinguish RSTR and LTBI phenotypes under unstimulated and M. tuberculosis-infected conditions. RSTR and LTBI global transcriptional differences under unstimulated (A) or M. tuberculosis-infected (B) conditions were probed using FGSEA and “hallmark” curated gene sets (MSigDB). The effects of M. tuberculosis on global expression differences were similarly compared (Fig. S1). Gene sets that were significantly enriched (FDR, <0.1) both by phenotype (A and B) and by stimulation condition (Fig. S1A and B) across Uganda and South Africa are shown. Global expression patterns for each gene set can be grouped by direction of enrichment across analyses (i to iv), most of which are highly concordant across clinical sites (i to iii). For each analysis, a normalized enrichment score (NES) is plotted according to the direction of enrichment for the Uganda (square) and South Africa (triangle) data sets, with color indicating significant enrichment (FDR, <0.1). NS, nonsignificant (FDR, ≥0.1).

Since RSTR monocytes have greater upregulation of the TNFα hallmark and inflammatory response pathways in response to M. tuberculosis infection, we hypothesized that phenotype-specific modulation of these pathways is involved in resistance to TST/IGRA conversion. We first looked for candidate resistance genes that were differentially expressed among RSTRs at each study site. Using an interaction model, which includes an interaction term between the M. tuberculosis infection (“M. tuberculosis group” versus “medium”) and clinical phenotype (RSTR versus LTBI) main effects, we identified 260 differentially expressed genes (DEGs) in Uganda and five DEGs in South Africa (false-discovery rate [FDR], <0.2) (Fig. 3; Tables S3 and S4). Principal-component analysis showed clustering by infection status and not phenotype (Fig. S3), and post hoc pairwise contrasts of these DEGs indicated expression differences between RSTR and LTBI phenotypes were largely restricted to the M. tuberculosis-stimulated condition (Fig. S4). However, none of these candidate resistance genes overlapped across sites. To identify candidate RSTR genes shared across populations, we next explored gene-level features of the TNFα hallmark pathway limiting the analysis to genes that either were DEGs or were included in the leading-edge subset that defines the core enrichment from each GSEA (Fig. 4A and B). The leading-edge subsets showed substantial overlap across Uganda and South Africa with consistent enrichment under the medium-stimulated condition (RSTR expression decreased), M. tuberculosis-stimulated condition (RSTR expression increased), or both (Table 1). Interestingly, TNF was neither differentially expressed nor a member of the overlapping leading-edge subset across populations, suggesting the enrichment of TNFα hallmark genes among RSTRs is not strictly a function of TNF expression. Six Uganda DEGs were leading-edge members under both the medium- and M. tuberculosis-stimulated conditions, including ABCA1, DUSP2, NR4A2, and TNFAIP3 (Fig. 4C). Overall, our results identify RSTR genes related to TNF-α signaling, yet distinct from TNF, that we hypothesize play important roles in mediating resistance to TST/IGRA conversion.

FIG 3

Differentially expressed genes distinguish RSTR and LTBI phenotypes, but are not shared across Uganda and South Africa. Differentially expressed genes were identified using an interaction model that incorporates an interaction term (“Mtb:RSTR”) in addition to the main effects phenotype (RSTR versus LTBI) and stimulation (medium versus M. tuberculosis). Using this model, 260 DEGs were identified in Uganda and 5 were identified in South Africa (FDR, <0.2). Volcano plots for the Uganda (top) and South Africa (bottom) analysis indicate changes in gene expression in response to M. tuberculosis stimulation (M. tuberculosis versus medium, log₂ fold change) that contrast RSTR and LTBI.

FIG 4

Core Uganda enrichment of TNF-α signaling gene set contains multiple differentially expressed genes. Mean relative expression (log₂ fold change) between RSTR and LTBI is plotted for each gene in the “TNFα response via NF-κB” hallmark gene set using pairwise contrasts for the medium-stimulated (A) and M. tuberculosis-stimulated (B) conditions for Uganda (top) and South Africa (bottom) monocyte analyses. Leading-edge genes are indicated by darker dots. Gene labels reflect DEGs identified in Uganda from the interaction model that incorporates main effects and the stimulation:phenotype term. Gene label shading indicates a DEG that is also a leading-edge member (dark label) versus DEGs that did not contribute this core enrichment (light label). (C) Box plots for select Uganda DEGs ABCA1, DUSP2, NR4A2, and TNFAIP3, each of which contributed to the GSEA core enrichment are shown: the median (line), interquartile range (box), and 1.5× interquartile range (whiskers) of normalized expression (log₂) for each subject (dot) were plotted.

TABLE 1

Overlap of TNFα gene hallmark leading-edge subsets

Gene in subseta
Medium (RSTR down)		M. tuberculosis (RSTR up)
Uganda (33/98 LE genes)	South Africa (35/92 LE genes)	Uganda (37/96 LE genes)	South Africa (12/47 LE genes)
CCL20	CCL20	CCL20	CCL20
IL7R	IL7R	IL7R	IL7R
INHBA	INHBA	INHBA	INHBA
PTX3	PTX3	PTX3	PTX3
TNFAIP6	TNFAIP6	TNFAIP6	TNFAIP6
TNIP1	TNIP1	TNIP1	TNIP1
ABCA1	ABCA1	ABCA1
AREG	AREG	AREG
CCL4	CCL4	CCL4
CXCL1	CXCL1	CXCL1
CXCL2	CXCL2	CXCL2
CXCL3	CXCL3	CXCL3
DUSP2	DUSP2	DUSP2
DUSP4	DUSP4	DUSP4
EHD1	EHD1	EHD1
GPR183	GPR183	GPR183
HES1	HES1	HES1
IFIT2	IFIT2	IFIT2
	IL1B	IL1B	IL1B
IL23A	IL23A	IL23A
	MARCKS	MARCKS	MARCKS
NFKB2	NFKB2	NFKB2
NFKBIA	NFKBIA	NFKBIA
NR4A1	NR4A1	NR4A1
NR4A2	NR4A2	NR4A2
OLR1		OLR1	OLR1
PDE4B	PDE4B	PDE4B
PFKFB3	PFKFB3	PFKFB3
PNRC1	PNRC1	PNRC1
	PTGS2	PTGS2	PTGS2
	SLC2A6	SLC2A6	SLC2A6
SOCS3	SOCS3	SOCS3
SOD2	SOD2	SOD2
TNC		TNC	TNC
TNFAIP3	TNFAIP3	TNFAIP3
ZC3H12A	ZC3H12A	ZC3H12A
ZFP36	ZFP36	ZFP36

a

Genes are listed if they were members of the leading-edge (LE) subset in at least 3 of 4 gene set enrichment analyses (columns) using the hallmark “TNFα signaling via NF-κB” gene set (MsigDB). Expression of these LE genes was lower (medium condition) or higher (M. tuberculosis condition) in RSTR than in LTBI monocytes. Of the 200 genes in this gene set, the LE subset consisted of 98 genes (Uganda medium), 92 genes (South Africa medium), 96 genes (Uganda M. tuberculosis), and 47 genes (South Africa M. tuberculosis) in the respective analyses. Gene order is alphabetical, and genes are grouped by overlap across all 4 analyses (above line) or 3 of 4 analyses (below line). Boldface genes are included among the 260 DEGs in Uganda (interaction model).

Polymorphisms in Uganda DEGs associate with the RSTR phenotype.

Genetic variations at multiple chromosomal loci have previously been associated with resistance to TST conversion (19–23). We next used a candidate gene association study to examine whether single nucleotide polymorphisms (SNPs) of ≤5 kb cis to the 260 DEGs were associated with the RSTR clinical phenotype among Ugandan participants (n = 74 RSTR and 189 LTBI subjects). Among these 5,248 examined SNPs, we identified 11 SNPs in 10 genes that associated with the RSTR or LTBI phenotype (P < 0.005) (Table S4). These 10 genes included CIITA (rs6498130), which was also a leading DEG (log₂ fold change, −0.79; FDR, 1.5E−04), for which expression increased in LTBI monocytes but decreased in RSTR monocytes following M. tuberculosis stimulation (Fig. S5). Among the other associated SNPs, those for EPB41L3 (rs1719945), ZNF184 (rs1883216), AKT3 (rs12144559), and KLHL29 (rs1530045) are known expression quantitative trait loci (cis-eQTL) in various tissues represented in a public database (GTExPortal; https://www.gtexportal.org). This supports a functional role for the proteins encoded by these genes. None of the 5,248 SNPs queried was significant after correction for false discovery (FDR, <0.2), and while this method may be overly stringent considering linkage disequilibrium reduces the actual number of total comparisons, these SNP associations require confirmation.

Networks of Uganda DEGs highlight TNF and inflammatory signaling pathways.

Several DEGs in Uganda with the largest fold change expression values are related to proinflammatory signaling, including CIITA, CXCL9, IFNG, and IL-2 (Fig. 3; Fig. S5). To explore biologic pathways that relate to these and the other 260 Uganda DEGs, we used STRING network analysis to assign connections among the DEGs that are included in the STRING database (https://string-db.org). We found a large cluster (n = 77 genes) of highly interconnected DEGs, five smaller networks, and 155 DEGs that had two or fewer connections (Fig. 5A). To assign biologic function to each network, we used topGO enrichment analysis, where DEGs were labeled according to terminal branches of the hierarchical Gene Ontology (GO) network. The core of the large DEG network related to regulation of cytokine production, TNF signaling, and IFN-γ signaling. Furthermore, this core included CIITA and AKT3, each of which has SNPs that independently associated with clinical phenotype. Additional subclusters within this network showed enrichment for leukocyte degranulation/homeostasis, cell-cell junction organization and other GO terms. When the STRING networks are labeled by relative expression (RSTR versus LTBI), the IFN-γ subcluster clearly shows higher expression among LTBI monocytes compared to RSTR monocytes following M. tuberculosis stimulation (e.g., IL2, IFNG, and CXCL9), whereas DEGs in the remaining subclusters, including TNF signaling and leukocyte degranulation, have higher expression among RSTR monocytes (e.g., TNFAIP3, TNFRSF1B, FOXO3, and LIMS1) (Fig. 5B). The smaller networks were enriched for GO terms related to RNA and rRNA processing, glycoprotein synthesis, and protein polyubiquitination. One network included ST3GAL1, which also had an SNP associated with clinical phenotype, but generally the genes in these small networks had modest expression differences between RSTR and LTBI phenotypes. Overall, these DEG networks reinforce conclusions from our global gene set enrichment analysis demonstrating transcriptional differences between RSTR and LTBI monocytes that include elevated expression of IFN-γ signaling among LTBI donors whereas other inflammatory pathways, including TNF, have higher expression among RSTR monocytes following M. tuberculosis infection.

FIG 5

Network analysis of Uganda DEGs identify multiple biologic pathways, including TNF and IFN-γ signaling, that correlate with RSTR and LTBI relative expression. Six gene networks among the 260 Uganda DEGs were identified by STRING. To identify biologic functions of these gene networks, we used topGO enrichment analysis to separately analyze the large network (n = 77 DEGs) and each smaller network of 3 or more genes individually. Significantly enriched gene sets (FDR, <0.05) for each cluster were identified using Fisher’s exact test and plotted in topGO, and Gene Ontology (GO) terms were selected from terminal or near-terminal branches of the GO hierarchical networks to capture the most specific GO terms and reduce redundancy. Each DEG was then colored according to its most specific (lowest on tree) enriched GO terms (A) or according to log₂ fold change (RSTR versus LTBI) expression values in M. tuberculosis-stimulated monocytes (B).

DISCUSSION

The immunologic mechanisms underlying resistance to TST/IGRA conversion after heavy exposure to M. tuberculosis remain poorly defined. By leveraging rare patient cohorts with rigorous epidemiologic characterization, our study explores pathways enriched across RSTR subjects from two distinct populations. We identify differentially expressed genes and pathways that distinguish transcriptional responses to M. tuberculosis in RSTR and LTBI monocytes. Our results highlight differences related to inflammatory signaling pathways, including the TNF response, which surprisingly was enriched with higher expression in M. tuberculosis-infected RSTR monocytes than in LTBI controls.

TNF activates antimicrobial programs in M. tuberculosis-infected macrophages. TNF is essential for protection from M. tuberculosis disease in mice (24), and the use of multiple anti-TNF biologic agents is associated with an increased risk for TB disease in humans (25–27). Whether this TNF-mediated protection from disease extends to a protective effect at the time of exposure, as might be measured in RSTR subjects, is not currently known. In a previous study of Indonesian household contacts followed for 14 weeks, those who remained IGRA negative (termed “early clearers”) had increased TNF and proinflammatory cytokine responses after ex vivo whole-blood stimulation compared to IGRA converters, suggesting that heightened proinflammatory responses at early times postexposure may distinguish RSTR subjects from those who become infected (2). Genetic studies have also associated the TNF response with resistance to TST conversion. Genome-wide linkage studies identified a locus on chromosome 11p14 (TST1) associated with TST negativity in regions of high (19) and low (20) endemicity for tuberculosis. This locus was later found to be a quantitative trait locus (QTL) for TNF secretion in response to Mycobacterium bovis BCG (28) and supports the hypothesis that RSTR leukocytes have a distinct TNF response to mycobacterial exposure. TNF may also arise from antigen-stimulated T cells and contribute to IFN-γ-independent cellular immunity that may characterize the RSTR phenotype (6); however, these TNF responses were more likely among LTBI compared to RSTR subjects in Uganda (1). An important distinction from these prior studies is that TNF was not a DEG in either population. Accordingly, it is unlikely the protective RSTR response is explained solely by elevated TNF expression by myeloid cells following M. tuberculosis infection; rather, we hypothesize that other TNFα hallmark genes contribute to resistance to TST/IGRA conversion.

Several DEGs were leading-edge members in the Uganda TNFα hallmark GSEA, both under unstimulated (RSTR expression down) and under M. tuberculosis-stimulated (RSTR expression up) conditions. Accordingly, the cellular pathways in which DUSP2, ABCA1, and NR4A2 participate may distinguish RSTR and LTBI functions by modulating intracellular M. tuberculosis infection. Dual-specificity phosphatase 2 (PAC-1), encoded by DUSP2, is induced by Toll-like receptor stimulation and associated with proinflammatory (29) and anti-inflammatory phenotypes (30) that are mediated through its negative regulation of mitogen-activated protein kinase (MAPK) and STAT3 pathways in macrophages and T cells. Transcriptional analysis of M. tuberculosis-infected macrophages (31, 32) previously identified ABCA1 as one of the most transcriptionally upregulated genes, which may relate to its role as a cholesterol efflux pump located at the plasma and endosomal membranes, where it could regulate M. tuberculosis cholesterol uptake (32) and is linked to the IFN-γ response (31). Finally, NR4A2 is a nuclear receptor and early response gene that is upregulated by lipopolysaccharide (LPS), cytokines, and oxidized lipids (33), and as a transcription factor, it has pleotropic activities, including the regulation of macrophage polarity (34) and fatty acid metabolism (35). Further studies are required to establish a role for these gene pathways in modulating intracellular M. tuberculosis replication or resistance to TST/IGRA conversion, and if validated, they may provide potential targets for therapeutic intervention or the development of biomarkers for M. tuberculosis exposure.

In contrast to the increase in M. tuberculosis-stimulated expression of TNFα hallmark (and inflammatory response) gene sets among RSTR subjects, expression of the IFNγ response and IFNα response gene sets remained lower in RSTR subjects in unstimulated and M. tuberculosis-infected monocytes. In Gambian subjects who were followed for 3 months with TST/IGRA testing after household M. tuberculosis exposure, an increased type I IFN response among IGRA nonconverters was detected using transcriptomic profiling of unstimulated whole blood (3). The discordance in the type I IFN enrichment in these findings may relate to the difference in specimens (whole blood versus monocytes) and differences in study design, such as the timing of household exposure, which was recent among Gambian contacts and either remote (Uganda) or undefined (South Africa) in our analysis relative to blood collection. Alternatively, considering most genes in the IFNα response core enrichment are also members of the IFNγ hallmark gene set (see Fig. S2 in the supplemental material), it is likely the enrichment of these pathways among LTBI transcriptomes primarily reflects the IFN-γ response.

Both the global monocyte IFN-γ gene set enrichment and the identification of DEGs related to IFN-γ or IL-2 signaling (e.g., CIITA, IFNG, CXCL9, and IL2) were identified among LTBI monocytes. CXCL9 expression in whole blood was also found to be increased at baseline in Gambian subjects who ultimately converted TST after household exposure (3). As T lymphocyte cytokines, increases in IFNG and IL2 transcripts in M. tuberculosis-stimulated LTBI monocyte cultures may result from lymphocyte impurity in our monocyte cultures, but several alternative considerations should be made. First, the short stimulation (6 h) may restrict the breadth or magnitude of any paracrine effects due to IFN-γ secretion by lymphocytes. Second, since the IFN-γ gene set was significantly enriched in unstimulated LTBI monocytes, differences between LTBI and RSTR transcriptional programs are not solely a function of the IGRA-based categorical phenotype definitions that require the presence of M. tuberculosis antigens. Third, a distinction among these top DEGs should be made with respect to CIITA, which encodes a polymorphism (rs6498130) that independently associated with the RSTR phenotype in Uganda and recently was linked to M. tuberculosis susceptibility (36). Unlike other top DEGs involved in the T cell response (IFNG, IL2, and CXCL9) where expression increased following M. tuberculosis stimulation largely in LTBI monocytes (Fig. S5A), CIITA expression decreased following M. tuberculosis infection specifically in RSTR monocytes. As a key transcriptional cofactor for major histocompatibility complex class II (MHC-II) expression, the decreased expression of CIITA in RSTR monocytes following M. tuberculosis infection is of particular interest in that it may result in reduced MHC-II antigen presentation and, consequently, limit downstream antigen-specific IFN-γ⁺ CD4 T cell responses that fundamentally distinguish the RSTR and LTBI phenotypes. CIITA also regulates transcription of other genes that may regulate RSTR pathways, including BTN2A2, which also was suppressed following M. tuberculosis infection among RSTRs (Fig. S5B). Interestingly, deficiency of BTN2A2 in antigen-presenting cells is linked to increased CD4/8 T cell proliferation and decreased frequency of T regulatory cells (37, 38), suggesting relative BTN2A2 deficiency in RSTR myeloid cells could also enhance T cell functions. Taken together, we identified several top DEGs, including IFNG, that are expected based on LTBI phenotype classification, but the relative expression patterns of CIITA and other DEGs suggest additional pathways that contribute to the RSTR phenotype.

Our study has several limitations. We did not find any individual genes with significant differential expression in both Uganda and South Africa, which likely reflects either the smaller cohort of gold miners or heterogeneity of these cohorts. Accordingly, transcriptomic analyses from additional RSTR cohorts are needed to validate our findings. Enrichment of the TNFα hallmark genes and other gene sets in both Uganda and South Africa cohorts suggests common RSTR mechanisms despite their significant heterogeneity; however, narrowing subsequent analyses to shared genes that drive these enrichments may be required to uncover RSTR mechanisms or to identify molecular targets for host-directed therapies. As discussed previously, the degree to which impurities in our monocyte cultures influenced our differential gene expression due to lymphocyte cytokine paracrine signaling (e.g., IFNG and IL2) is unknown, but must be considered for genes downstream of these pathways. Another limitation is that resistance to TST/IGRA reactivity is an imperfect surrogate for resistance to infection that cannot be defined microbiologically, and whether RSTR subjects have lower risk of TB progression following exposure is unknown. Finally, to accommodate limitations in participant cell availability, we used a single laboratory M. tuberculosis strain (H37Rv) and chose a single intermediate RNA time point (6 h postinfection) to compare transcriptional profiles from readily accessible peripheral blood. It is likely additional phenotypes would be uncovered according to M. tuberculosis strain and at other time points. Furthermore, similar studies in alveolar macrophages isolated from these cohorts (39), which are the earliest targets for M. tuberculosis infection in vivo, may yield different results.

In summary, our study highlights genes and gene pathways that distinguish RSTR and LTBI monocyte responses to M. tuberculosis. These pathways suggest differences in inflammatory responses that may play important roles in mediating resistance to M. tuberculosis infection in vivo. Further studies are required to evaluate whether these responses mediate the early clearance of M. tuberculosis from infected myeloid cells or, alternatively, modulate antigen presentation or otherwise impact cell-mediated protective immunity. Elucidating these RSTR pathways may identify targets that augment antimicrobial therapy and vaccine design or may inform diagnostics that better identify M. tuberculosis-exposed individuals to help curb the epidemic.

MATERIALS AND METHODS

Subject recruitment.

Recruitment of HIV-negative Uganda household contacts was previously described (8). HIV-negative gold miner participants from South Africa who had worked in the industry for ≥15 years were recruited as part of the Highly Exposed TB Uninfected (HETU) study and underwent baseline TST and IGRA testing; the majority were followed for 12 months for repeat TST/IGRA (10, 40). Informed consent for participation was obtained from each subject at the time of enrollment as detailed in each study protocol that was approved by the respective institutional review boards. RSTR and LTBI phenotype classifications were reported previously with details available in Text S1 in the supplemental material (13).

Culture of M. tuberculosis and monocyte infections.

Log-phase cultures of virulent laboratory Mycobacterium tuberculosis strain H37Rv (a gift from David Sherman) were cultured in supplemented 7H9 medium (Middlebrook) and stored in aliquots at −80°C until monocyte infections to avoid heterogeneity between batches (supplemental material). Monocytes were isolated from cryopreserved peripheral blood mononuclear cells (PBMCs) from each participant and cultured in monocyte colony-stimulated factor (M-CSF) to establish adherence as previously described (13). To minimize artificial homogenization of transcriptional responses due to ex vivo culture conditions, experiments were completed within 48 h of PBMC thaw and M-CSF treatment. In a biosafety level 3 (BSL3) laboratory, monocytes were stimulated either with H37Rv (multiplicity of infection [MOI], 1.0) or medium, and RNA was isolated in TRIzol (Invitrogen) after 6 h. These samples were obtained in singlet for donors from South Africa, whereas for Uganda subjects, two independent experiments (Exp A and B in Fig. 1) with partially overlapping subjects were each performed in singlet.

RNA sequencing and data processing.

Preparation of cDNA libraries, RNA sequencing and alignments using STAR2.6.0a (41), filtering, and analysis were performed in R v4.0.2 (42, 43) (supplemental material). ComBat-seq (44) was used to correct batch effects for Uganda Exp A and B, while controlling for stimulation, phenotype, sex, and age. Counts were converted to log₂ counts per million using voom (45).

Differential gene expression, gene set enrichment analyses, and STRING network analysis.

We previously analyzed transcriptomes from unstimulated monocytes from South Africa and a subset of unstimulated Uganda monocytes (Exp A in Fig. 1) (13). To instead identify expression patterns that distinguish RSTR and LTBI phenotypes upon M. tuberculosis stimulation, we selected an expression model that incorporated an interaction term in addition to the main effects: Expression ~ phenotype + stimulation + phenotype:stimulation + covariates, with patient and genetic kinship (when available) included as random effects using R packages coxme (46) or lme4 (47) when kinship data were not available. We adjusted for age, sex, and sequencing batch (Uganda) or age alone for the South Africa samples, which were sequenced simultaneously and given the gold miner participants were all male. Our model sigma and outputs were stable with and without these adjustments (data not shown). Differentially expressed genes (DEGs) in the interaction term (FDR, <0.2) were then assessed using pairwise contrasts of the four phenotype:stimulation groups (e.g., contrasting phenotypes for each stimulation condition or contrasting stimulations for each phenotype) corrected for the same covariates and random effects as the interaction model (supplemental material). To understand biologic connectivity between each of the 260 Uganda DEGs, we used STRING v11 network analysis (48) followed by enrichment analysis using topGO (49). Gene set enrichment analysis (GSEA) was performed using the Molecular Signatures Database (MSigDB v7.2) (16) Hallmark and Gene Ontology (GO) collections. Fast gene set enrichment analysis (FGSEA) (14) was used to compare fold changes of all genes in phenotype:stimulation pairwise contrasts.

Data availability.

Access to raw transcriptomic data is available through the NCBI database of Genotypes and Phenotypes (dbGaP) Data Browser (https://www.ncbi.nlm.nih.gov/gap/) under accession no. phs002445.v1.p1 (for Uganda) and phs002746.v1.p1 (for South Africa), but first must be approved by data access committees (DACs) for each study site (supplemental material). All R code is available at https://github.com/hawn-lab/RSTR_RNAseq_Mtb_public.