The human TB lesion signature shows a distinct and heterogeneous transcriptional profile as compared with non-lesional lung tissue

In our study, outlined in Fig. 1, we collected 48 samples from 14 individuals and analysed 44 paired samples from 13 individuals (6 DS-TB and 7 MDR/XDR-TB) to evaluate the human TB lung granuloma transcriptomic changes by RNA sequencing.

Although the patients included in this study exhibited normal to high BMI, low CRP levels, relatively low SGRQ scores and were considered microbiologically cured, they nonetheless required lung resection surgery due to the persistence of TB cavities.

We analysed total RNA from three different sections: Central Lesion (C; n = 6), Internal Wall (I; n = 12) and External Wall (E; n = 13) collected from each patient’s lesion biopsy. Fewer C- samples could be analysed compared to I and E, due to poorer RNA recovery. Additionally, surrounding non-lesional (NL) tissue from the involved lung was collected as a comparator (n = 13) (Fig. 2a). Patients were matched according to their sex and Mtb drug-sensitivity classification to avoid potential confounding factors (Supplementary Table 1). Moreover, clinical and demographic data, and resected TB lesion characteristics and pathology were assessed at the time of surgery and are reported for each participant (Supplementary Table 1 and Supplementary Data 2).

TB lesion samples were collected from each patient included in the SH-TBL cohort: central lesion (C), internal wall (I) and external wall (E) and, altogether, samples from each patient represent the human TB lesion. An additional sample from surrounding non-lesional lung tissue (NL) was also collected from the same patient as control (a). 48 samples from 14 patients (6 DS-TB and 8 MDR/XDR-TB) were RNA sequenced to evaluate the human TB lung lesion transcriptomic changes. A set of 4630 DEGs was identified after comparing the human TB lesion counts with NL lung tissue expression, using DESeq2 with adjusted p < 0.05. b heatmap depicts the top 40 DEGs ranked by the adjusted p-value comparing the human TB lesion versus NL lung tissue expression profiles (44 paired samples from 13 patients). The intensity of each colour denotes the standardized ratio between each value and the average expression of each gene across all samples. Red pixels correspond to an increased abundance of mRNA in the indicated sample, whereas blue pixels indicate decreased mRNA levels. Source data are provided as a Source Data file. Image in (a) was created in BioRender. Vilaplana, C. (2025) https://BioRender.com/x16o926.

We found a total of 4630 significantly differentially expressed genes (DEGs), using DESeq2 with adjusted p ≤ 0.05 (Supplementary Fig. 1a). Of these, 2496 genes were over-expressed in lesion tissues, whereas 2134 were under-expressed, as compared to NL lung tissue (Supplementary Fig. 1a). The top 40 ranked DEGs clearly separated lesion samples from NL lung samples (Fig. 2b), showing distinct transcriptional profiles for the two tissues. Among them, genes involved in immune system/cytokine signalling (IRF4, CCL19, LTB, JAK3, INPP5D, FCER2, MMP1) and B cell activation and differentiation (CD22, BLNK, CARD11) were over-expressed, suggesting an inflammatory signature in the TB lesion.

Seven TB lesion samples clustered together with NL tissue samples, consisting of six samples from the external compartment and one from the internal compartment. This observation may suggest a transitional transcription profile across the lesions, particularly evident in the external tissue due to its proximity to the NL samples, but also not discarding the heterogeneity in the transcriptional profiles of the lesions (Fig. 2b).

Altogether, our data show a distinct segregation of the TB lesion when compared to the NL lung tissue with respect to an inflammatory profile, as previously proposed¹¹. Our findings also indicated a range of molecular diversity within the TB lesion samples, prompting our decision to delve deeper into the heterogeneity at a transcriptional level.

Compartments within the TB lesion reveal distinct gene expression profiles with an enriched inflammatory response across the lesion

To further explore the TB lesion heterogeneity and investigate the contribution of each compartment, we first performed an enrichment analysis derived from single sample Gene Set Enrichment analysis (ssGSEA) using the top 40 DEGs discriminating the TB lesion from the NL lung tissue. The expression of these genes in the different tissue compartments revealed a more pronounced enrichment score of these DEGs in central and internal lesion samples, suggesting that these two compartments might be the main contributors for the overall TB lesion transcriptional signature (Fig. 3a).

a show the enrichment score derived from single sample analysis GSEA using the top 40 genes discriminating TB lesion (G) from NL lung tissue. Data on the enrichment for each compartment (C, I, E and NL) are represented as medians with an interquartile range (IQR). Boxplots show minimum and maximum values, the interquartile range (IQR, 25th to 75th percentile), and the whiskers representing 1.5 times the interquartile range. Outliers are indicated as individual points outside the whiskers. Statistical analysis was performed by applying the two-sided t-test. Statistical differences refer to a p-value < 0.05. In (b) the heatmaps show differences in the top 40 ranked genes from DESeq2 with adjusted p < 0.05 by separately comparing the central (C), internal (I) and external (E) compartments with the NL lung tissue gene expression derived (). The intensity of each colour denotes the standardized ratio between each value and the average expression of each gene across all samples. Red pixels correspond to an increased abundance of mRNA in the indicated sample, whereas blue pixels indicate decreased mRNA levels. c pictures modular transcriptional of the seventeen modules of co-expressed genes derived from WGCNA for our TB lesion dataset separated by compartment. Fold enrichment scores derived using QuSAGE are depicted, with red and blue indicating modules over or under expressed compared to the control. Only modules with fold enrichment (FDR) < 0.1 were considered significant. Source data are provided as a Source Data file.

Next, we compared the expression profiles derived from each TB lesion compartment with the NL tissue. The list of DEGs (DESeq2 with adjusted p ≤ 0.05) for the C, I and E vs NL tissue comparisons respectively constituted 3228 (1539 genes were over-expressed, whereas 1689 were under-expressed); 5275 (2676 over-expressed and 2599 under-expressed); and 1045 genes (552 over-expressed and 493 under-expressed) (Supplementary Fig. 1b). For central and internal compartments, the hierarchical clustering of the 40 most significant DEGs showed an evident separation when compared each compartment against the NL lung tissue (Fig. 3b). Though less noticeable, the external compartment was still distinguishable from the NL tissue. Therefore, the magnitude of differential expression relative to NL decreased gradually towards the edge of the TB lesion structure, including between adjacent compartments (Fig. 3b and Supplementary Fig. 2).

Among the highly variable genes in central lesion, we found genes involved in the immune system/cytokine signalling (CCL19, CXCL10) to be upregulated in comparison to the NL tissue (Fig. 3b). On the other hand, we found extracellular matrix organization-related genes to be downregulated (LRP4, MUC15), while others were upregulated (ADAM12, CTSK). Moreover, collagen-encoding genes (COL1A1, COL3A1, COL11A1) were upregulated in the central compartment, which could reflect the fibrosis observed in all patients’ lesions (Supplementary Data 2); as well as of genes associated to immunoglobulin heavy and light chains (IGHV4 − 61, IGKV1 − 39, IGKV1D − 39, IGHV1 − 18, IGHV3 − 74, IGLV1 − 40, IGHV4 − 34, IGHV1 − 3, IGLV1 − 51, IGLV2 − 18, IGHV6 − 1, IGHV4 − 55), related to humoral immunity (Fig. 3b). Furthermore, genes involved in complement fixing (C1QA, C1QB, C1QC) were significantly upregulated, although not among the top 40 DEGs (Supplementary Data 1). For the internal compartment, genes involved in the immune system/cytokine signalling (LTB, FCMR, AIM2, CXCL10, IRF8, IRF4) were upregulated compared to NL (Fig. 3b). Furthermore, immune system/cytokine signalling genes (LTB, CCL19, CXCL9, TNFAIP3, TNFRSF13C, FCMR, AIM2, CXCL10) were over-expressed in the external compartment relative to NL (Fig. 3b), evidencing an inflammatory signature throughout the lesion.

We then applied weighted gene co-expression analysis (WGCNA) to perform a modular analysis of co-expressed genes in the TB lesions and in the three compartments separately, comparing all samples to NL control tissues. We identified 17 modules from co-expression networks related to the whole human TB lesion (Fig. 3c and Supplementary Data 2). The identified TB lesion modular signature showed that neutrophil degranulation, cell signalling, adaptive/humoral immunity, extracellular matrix, interferon/cytokine signalling, and innate/pathogen recognition receptors (PRR) modules were overabundant. These observations were consistent throughout the compartments, except for the neutrophil degranulation and innate/PRR modules, which were apparent in the total lesion and the internal lesion only, but not in the central or external lesions (Fig. 3c). Conversely, the Epithelial to Mesenchymal Transition (EMT), cholesterol biosynthesis and metabolism modules were found to be underabundant in the whole lesion and external and internal, but not in the central compartment. In addition to cholesterol and metabolism, the organelle biogenesis module was underrepresented, across all compartments, suggesting that some pathways present in the healthy lung are diminished in the lesion (Fig. 3c).

To broaden our understanding on the distribution of the immune response among the modular signature in the whole lesion and its compartments, we have used the LM22 signature matrix to profile the distinct human hematopoietic cell populations²¹. We found that the adaptive/humoral and the innate/PRR modules presented most of the genes related to the immune populations (Supplementary Fig. 3a), with said populations also varying in proportions across TB lesion compartments (Supplementary Fig. 3b). By expanding these modules, we found that most of the submodules composing them were significantly differently enriched when comparing the compartments, particularly for the adaptive/humoral submodules (Supplementary Fig. 4a).

In summary, our results showed a significant enrichment of modules related to inflammation, including pathways of innate immunity in the TB lesion, the central and internal compartments, and of adaptive/humoral immunity across all compartments. Meanwhile a decrease in modules related to extracellular matrix organisation and cholesterol biosynthesis and metabolism was observed in the lesion. Furthermore, profiling of the LM22 populations as well as the expansion of the adaptive/humoral and innate/PRR modules revealed a differential distribution of the immune cells in the different compartments, contributing to the identified modular signature.

Patients’ clinical status is associated with differential modular transcriptomic profiles in TB lesions

The heterogeneity in the host immune response to infection, considering the involvement and contribution of physically distinct compartments, together with the bacteria and the inflammatory environment, defines granuloma fate and disease manifestation^19,22. Hence, we next aimed to associate the modular signature changes in the TB lesion (considering the three compartments together) with clinical data (Supplementary Data 2), using surrogates of treatment response and disease severity (SGRQ symptoms sub-score; being a fast or slow sputum culture converter; DS vs MDR-TB case; being a relapse or new TB case and number of lesions present in the CXR). We quantitatively tested the association of each clinical parameter with each of the significant module’s eigengene (ME) expression patterns (Wilcoxon p ≤ 0.05).

Regarding the sputum culture conversion (SCC), the modular signature of the TB lesion revealed a significant association of DNA binding and interferon/cytokine signalling modules with SCC, with the enrichment of these modules being significantly higher in those individuals converting the sputum culture later (FDR < 0.1; Fig. 4a, b). No significant modular expression was found to be associated with Mtb drug sensitivity of the individuals, relapsed or new cases, or number of lesions (Supplementary Data 2).

Modular analysis of RNA-seq data from TB lesions of 14 patients. Patients were clinically defined accordingly to sputum culture conversion (SCC) and TB disease impact on lung function, measured using the Saint George’s Respiratory Questionnaire (SGRQ) symptom score, as surrogates of treatment response and TB severity. Heatmap represent the key TB lesion modules significantly associated to individual’s’ clinical surrogates of TB severity and treatment response (a). Fold enrichments were calculated for each WGCNA module using hypergeometric distribution to assess whether the number of genes associated with each clinical status is larger than expected. Fold enrichment scores derived using QuSAGE are depicted, with red and blue indicating modules over or under expressed compared to the control. The colour intensity represents degree of perturbation. Modules with fold enrichment scored FDR p-value < 0.1 are considered significant. b, c show TB individuals’ stratification according to SCC (fast n = 28 or slow converters n = 20) and SGRQ symptom score (low impact if SGRQ < 20 with n = 23 or high impact if SGRQ > 20 with n = 25), respectively, and the significant association using their corresponding derived WGCNA significant eigengene modules (ME) (p < 0.05). Data are represented as median with an interquartile range (IQR). Boxplots show minimum and maximum values, the interquartile range (IQR, 25th to 75th percentile), and the whiskers representing 1.5 times the interquartile range. Outliers are indicated as individual points outside the whiskers. Statistical analysis was performed by applying the two-sided Wilcoxon-rank sum test. Source data are provided as a Source Data file.

When considering the severity of TB disease, in terms of a higher presence and severity of symptoms, we found that DNA binding, neutrophil degranulation, interferon/cytokine signalling, cholesterol biosynthesis and myeloid activation modules were significantly overabundant and associated with higher SGRQ symptoms score (FDR < 0.1; Fig. 4a and c), pointing to higher inflammation status with more severe disease manifestation. In contrast, the EMT module was significantly underabundant in these individuals’ TB lesions (Fig. 4c). When stratifying against the clinical data, results showed that there was no clustering of the clinical surrogates with neither the central nor the internal compartments (Supplementary Fig. 5).

Further examination of the associations between submodules and severity correlates revealed that only three submodules exhibited statistically significant differences among clinical surrogates (Supplementary Fig. 4b), derived from the innate/PRR module: the submodule Response to inflammation, linked to neutrophils and granulocytes, with genes related specifically to response to IL-1 and type II IFN and neutrophil chemotaxis; and the Innate response regulation submodule, linked to neutrophil, monocytes and macrophages, with genes related specifically to immune response activation and regulation. We also observed an enrichment in the CD4 + T helper lymphocyte response submodule. derived from the Adaptive/humoral module, linked to antigen-presenting cells and CD4 T cells, with genes specifically related to the regulation of T cell activation, lymphocyte differentiation and regulation of the adaptive immune response.

To gain further insights into the differences found between clinical surrogates, we next identified a set of seven transcription factors differentially expressed (ETV7, STAT1, AR, SOX5, ERG, ASCL2 and PRDM5) between fast and slow SCC and patients with less severe or more severe symptoms. Transcription factors corresponding to the IFN/cytokine signalling module were overexpressed in slow converters and/or more severe patients, whereas transcription factors belonging to the EMT module had higher expression in patients with less severe symptoms, complementing the modular analysis (Supplementary Fig. 6a, b and Supplementary Table 2).

To support RNA sequencing data, we validated by immunohistochemistry the protein products of three genes significantly expressed between TB lesion and NL tissue. In our data set, CXCL9, GBP5 and STAT1 are representative genes from the module with the highest enrichment in the whole TB lesion associated with TB surrogates of severity. We quantified the presence of the respective proteins in TB patient lesions compared to non-TB controls and found a significantly higher expression of these proteins in the TB patient lesions (Fig. 5a–d).

a shows representative immunohistochemistry staining for CXCL9, GBP5 and STAT1 from the TB lesion of a representative patient (TB-05) compared to a patient presenting bullous emphysema (TB-42), as non-TB control. The top row corresponds to whole sections of the TB lesion (at the left of the images) and of non-lesional tissue (at the right of the images). Scale bars correspond to 1000 µm. NC necrotic core, M macrophage region, F fibrotic region, L lymphocyte-enriched region, AS alveolar space. b–d show the quantification of CXCL9, GBP5 and STAT1 protein levels respectively in lesion sections of all TB patient (n = 14) compared to the non-TB control tissue sections (n = 3). n refers to biologically independent tissue sections from different individuals. Data on the percentage of stained area are represented as median with an interquartile range (IQR) Boxplots show minimum and maximum values, the interquartile range (IQR, 25th to 75th percentile), and the whiskers representing 1.5 times the interquartile range. Outliers are indicated as individual points outside the whiskers. Statistical analysis was performed by applying the two-sided Wilcoxon-rank sum test. Statistical differences refer to a p-value < 0.05. Source data are provided as a Source Data file.

All these data reinforce our findings from the transcriptional comparison between TB lesion and NL tissue, showing that a slower SCC, reflecting a poorer response to treatment and slower clearance of Mtb, and severer TB cases are associated to an increased inflammatory response at site.