Background

Mycobacterium tuberculosis is aerobically transmitted, obligate aerobic bacteria that is widely considered a global health crisis that affects nearly a quarter of the population and is the world’s deadliest infectious disease. Mycobacterium tuberculosis establishes dormancy in alveolar macrophages, remaining asymptomatic for years. It reactivates upon immune suppression, causing severe outcomes in HIV patients. Early signs of progression, response variability, and specific mechanisms remain unknown. The goal of this study is to investigate and propose biomarkers for early signs of progression from latent to active tuberculosis. 

Methods

GEO2R, SR Plot, KEGG, and GO were used to find the top 50 differentially expressed genes (DEGs), visualize, and analyze them. GEO dataset GSE189996 was used to find DEGs. This dataset used human, alveolar macrophage samples from 28 donors to investigate DEGs associated with host response variability. Groups were defined as control versus infected from time periods of 2, 24, and 72 hours. 

Results  

50 DEGs were identified based on padj, p-value, and stat (Wald Statistic). Thirty were upregulated and 20 were downregulated (see figure 2). GO analysis indicated enrichment in the biological process ontology, and ER-stress-induced-apoptosis and RNA-transcription-dysregulation being the most enriched. KEGG analysis identified MAPK and TNF signaling pathways including TNF, FOS, SIRT1, KLF4, KLF5, HERPUD1, PMAIP1 (NOVA), and EIF2AK3 (PERK) were top genes and pathways. 

Discussion

We propose that the biomarkers of increased SIRT1 expression and high KLF4:KLF5 as the most easily testable biomarkers of progression from latent to active tuberculosis due to being the upstream controllers of the proposed mechanisms.

Tuberculosis (TB), caused by Mycobacterium Tuberculosis, continues to be a health crisis. Glaziou et al.1 states that nearly one quarter of the population worldwide has latent TB infection or will develop it, where the bacterium is dormant within the human body in a Ranke complex without symptoms and is not contagious. De Martino et al. [2] describe active tuberculosis as a curable, aerially transmittable disease which can create symptoms of coughing, fevers, pain in the chest, etc. They additionally state that the majority of people infected with Tuberculosis will have inactive tuberculosis, also known as latent or dormant tuberculosis; latent TB has no symptoms and can exist within a host for years before showing symptoms, making its progression unpredictable. In 2016, TB caused 1.3 million deaths and is the leading infectious disease and 10.4 million new cases arose [1]. 

Current tests like Tuberculosis Skin Test [3] and QuantiFERON [4] can detect exposure from TB, but the focus of this study is to identify progression risks from latent to active. In a clinical setting, QuantiFERON could identify latent TB and an intervention based on genes identified from this study could identify progression risks [4]. This fixes costly issues in public health systems where they are forced to either treat all latent TB individuals which risks overtreating those who won’t progress, side effects, and drug resistance or treat none which misses preventable deaths [5]. In this study, it is investigated which pathways or genes in early M. TB infection can be used to identify patients of high risk of progression to active disease from latent within a defined timeframe and if these pathways or genes be used for targeted therapies? Moreover, do these signatures reveal specific pathways or genes that could be intervened to prevent progression? Because of its obligate aerobic nature it infects the lungs which allows it to spread aerobically [2,6]. The lipid cell wall (mycolic acid) causes it to be acid-fast, which requires the Ziehl-Neelsen gram stain rather than normal techniques. Then, M. TB is a facultative intracellular organism, meaning it can live in other cells, and more specifically, macrophages. Its cord factor allows it to evade immune response, and catalase peroxidase will enable it to resist host cell oxidation. M. TB creates a TH-1 response which eventually leaves a Ghon complex, a.k.a the caseating tissue (earlier stage of Ranke complex). If M. TB gets reactivated due to immunorepressors, then it will eventually lead to bronchopneumonia or spread through the vascular system to other organs in the body [2,6].

Finally, M. TB is found through screening using a TB Skin test (active) or QuantiFERON (latent), and if it’s positive, then a Chest x-ray into a sputum sample staining for acid-fast bacilli [2,7]. 

Bioinformatics studies compare gene expression profiles using transcriptomics (RNA-seq) and they have identified distinct transcriptional signatures that can identify progression from latent to active TB. For example, Burel et al. [8] have identified blood transcriptomic signatures that can identify high risk latent TB individuals based on reactions to chemprophylaxis. These signatures usually find upregulated immune and inflammatory pathways which align with the TB skin test checking for size of induration [3,8]. While transcriptional signatures can stratify latent TB progression risk months to years before symptom onset, validated biomarkers currently exist to predict progression during the earliest stages of infection. Moreover,  It is also unknown exactly how early these signatures can appear, how stable they are, or how influenced they are by other factors [2,6,9]. Additionally, distinguishing a real signature from unimpactful ones in longitudinal data; accounting for variation between patients, especially with TB; and translating complex multi-gene signatures into clinically available interventions can be extremely challenging.  The goal of this study is to propose genetic mechanisms of infection in order to identify molecular signatures in the host which can be used to predict the progression within 72 hours of a latent M. TB infection which can allow for more precise interventions for individuals found as high-risk. I hypothesize that latent TB individuals who later progress to active disease will exhibit distinct gene expression profiles compared to stable non-progressors. This research addresses the limitations of what exists currently for TB patients and can possibly fix the inefficiency of current latent to active TB treatments overall reducing deaths and optimizing the use of limited healthcare resources by focusing on what is most immediately needed.

In this study, the dataset Alveolar macrophage immunobiology and functional genomics: Unlocking human to human variation in host response to M. tuberculosis was collected from NCBI GEO2R10 by using the search details:

 (“mycobacterium tuberculosis”[MeSH Terms] OR “Mycobacterium tuberculosis”[Organism] OR mycobacterium tuberculosis[All Fields]) AND (“macrophages, alveolar”[MeSH Terms] OR alveolar macrophages[All Fields])

Following figure 1, searching this will provide 396 results and the focus of this study was on result 18 (where 17 was a follow-up, replication study). Then the 124 samples from dataset GSE18999611 was defined or categorized into groups “mtb-2hr” AND “mtb-24hr” AND “mtb-72hr” AND “control-2hr” AND “control-24hr” AND “control-72hr.” They were then analyzed using the no-code GEO2R bioinformatics tool that uses R programming language in order to generate the venn diagram, volcano plot, and list of differentially expressed genes and finally input into SRPlot.

Figure 1: Research Methodology of the  steps and bioinformatics used in this study. Research method workflow from data collection to functional analysis. Steps include dataset collection from NCBI GEO, differential gene expression analysis using GEO2R. Functional and pathway analysis of key genes was performed using SR Plot, integrating Gene Ontology (GO) and KEGG databases.

Identification of the Top Differentially Expressed Genes

To identify the most significant differentially expressed genes, statistical analysis was applied. This process used stat, p-value, and padj value to prioritize the most important genes out of 12381 based on their differential expression across samples. Prioritizing padj, stat, and p-value in that order. Stat being the Wald Statistic. 

Data Analysis Using SRPlot, KEGG, and GO Bioinformatics Tools

Then SRPlot12, KEGG13, and GO14,15 bioinformatics tools and databases were utilized to analyze the functions of these top genes. These tools helped uncover the potential roles of the genes in Mycobacterium Tuberculosis. SRPlot was used in tandem with KEGG and GO to generate the visuals shown in figure 3.

Identification of Differentially Expressed Genes 

GEO2R was initially used to identify the top 50 differentially expressed genes (DEGs) out of the primary 12381 genes identified from the samples. From the first two results, some genes were expressed differently between M. TB infected host cells and controls from all collected time periods (see figure 2). The volcano plot (Figure 2) is utilized to quickly identify the most meaningful changes in DEGs by combining a measure of statistical significance with the magnitude of change. The red dots represent genes which are upregulated, or more expressed; blue dots represent genes which are downregulated, or less expressed; and the black represents neither, or not significant. From the venn diagram (Figure 2), 1 DEG was common with all time periods; 3 were common with the earliest stage, 2 hours, and the latest stage, 72 hours; and 398 were common with the 24 hour stage and the 72 hour stage.

Figure 2: Differentially Expressed Genes: Volcano plot shows gene expression of patients infected with mtb versus the control after 72 hours while venn diagram shows the count of shared and separate DEGs between every recorded time period.

Identification of 50 Statistically Significant Differentially Expressed Genes  (DEGs)

To identify the top 50 DEGs, the 12381 genes were organized by padj, p-value, and stat. Padj is the adjusted p-value statistic, making it more significant.10 Stat, or the Wald Statistic, measures the statistical significance of data.10 Out of the top 50 genes selected, 30 of the genes were chosen with the highest stat value and were upregulated and 20 of the genes with the lowest stat value, downregulated, but maintained a padj value < 7.00E-15. 

Potential Functions and Enrichment of the Identified Genes and/or pathways

SR Plot was used to visualize and determine the function of genes. From the KEGG results, the TNF and MAPK signaling pathways were identified as statistically significant due to the high presence of top 50 DEGs. KEGG additionally identified TNF, FOS, FOSB, SIRT1, KLF4, KLF5, PMAIP1, HERPUD1, DUSP1, IRS2, ODC1, SLC16A1, DUSP2, JUNB, CXCL3, EIF2AK3, and GADD45B as significant genes of equal significance in discovered pathways. Out of these genes, we identified TNF, FOS, SIRT1, KLF4, KLF5, PAMIP1, HERPUD1, and EIF2AK3 as key genes. 

The GO results identified the regulation of endoplasmic reticulum stress-induced intrinsic apoptosis signaling pathway and positive regulation of pre-miRNA transcription by RNA polymerase II as the most major biological process pathways. GO results were analyzed in tandem with genes identified from KEGG to find related potential functions and genes.

Fig. 3

Functional pathway and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis of top 50 differentially expressed genes (DEGs) for Mycobacterium Tuberculosis. GO Bar plot of Biological Processes (BP), Molecular Functions (MF), and Cellular Components (CC) of upregulated DEGs (A), Bubble plot of upregulated DEG-enriched biological processes (B), KEGG Cnet plot of upregulated DEGs (C), and KEGG Pathway of the TNF signaling pathway (D) via SRPlot. 

Summary of Findings: 

The goal of this research was to identify molecular signatures in the host which can be used to predict the progression from latent to active TB which can allow for identifying high-risk individuals. GEO2R Analysis, gene ontology, and KEGG found the biological processes of the regulation of endoplasmic reticulum stress-induced intrinsic apoptotic signaling pathway and positive regulation of pri-mRNA transcription by RNA polymerase II as the most enriched and statistically significant processes. Gene ontology (GO) identified apoptosis and positive regulation of mRNA transcription as the two most enriched biological processes, and the CNET plot provides associated genes:

Positive regulation is associated with FOS, KLF5, KLF4, and TNF. 

ER stress-induced Apoptosis is associated with PMAIP1, EIF2AK3, HERPUD1, and SIRT1. 

Positive regulation of miRNA transcription is defined as any process that activates or increases the frequency, rate or extent of microRNA (miRNA) gene transcription [16]. The apoptotic signaling pathway in response to ER stress is the series of molecular signals in which an intracellular signal is conveyed to trigger the apoptotic death of a cell. The pathway is induced in response to a stimulus indicating endoplasmic reticulum (ER) stress, and ends when the execution phase of apoptosis is triggered. ER stress usually results from the accumulation of unfolded or misfolded proteins in the ER lumen [17]. FOS, KLF5, and KLF4 are all transcription factors. However, KLF5 and KLF4 are Kruppel-like factors which are associated with colorectal cancer. KLF4 acts as a tumor suppressor and KLF5 acts as an oncogene [18]. These genes are also in the MAPK signaling pathway [19]. TNF plays a crucial role in inflammation, the TNF signaling pathway, and the MAPK signaling pathway [20,21]. PMAIP1 or Noxa helps induce cell death.22,23 PERK (EIF2AK3) and HERPUD1 are used for ER stress response, inhibiting translation and removing misfolded proteins respectively [24]. SIRT1 is involved in a lot of processes as it deacetylates, but importantly it helps cells adapt to oxidative pressures, like HERPUD1, and is shown to have anti-inflammatory properties [25-27].

Interpretation of Results: 

In the RNA Polymerase II process, FOS and KLF5 seem to override [28,29] KLF4’s [30,31] tissue-protective signaling which causes immune polarization failure. FOS would continue unregulated and cause an influx of cytokine production. A high-risk signature would be a certain ratio of upregulated FOS and KLF5 to downregulated KLF4 (high KLF5:KLF4 ratio). In the ER stress-induced pathway, SIRT1 is likely deacetylating p5332, which induces PMAIP1, and suppressing apoptotic pathways, causing PERK hyperactivity. This leads to chronic stress memory in macrophages leading to a heightened inflammatory response upon subsequent stimulation, even without the original stressor [33]. Moreover, SIRT1 and PERK are likely making a pro-stress cellular environment [34] in latent TB, which prevents PMAIP1 driven apoptosis. In this case, latent TB patients with upregulated SIRT1 and PERK in monocytes can be tested for progression to validate this biomarker. Finally, TNF is possibly in a feedback loop where it induces FOS [35] to bind to the TNF promoter, which creates more TNF, overall fueling an apoptosis blockade via SIRT1. An important biomarker would be high TNF, FOS, and SIRT1 indicating a very high-risk individual.

Table 1: Summary of Results

Comparison with Previous Studies: 

ER Stress-induced Apoptosis Process

M TB infection inhibits apoptosis by upregulating factors like Mcl-1 which are anti-apoptotic and suppresses PMAIP1-mediated cell death. This aligns with Simper et al. [36] findings that M. TB induces PPARγ expression, leading to macrophage survival and bacterial persistence. PERK ER stress responses are critical in TB pathogenesis as studies [32,37] confirm that chronic ER stress in macrophages creates a “stress memory” which allows for bacterial survival and granuloma necrosis. SIRT1 deacetylating pro-apoptotic proteins and blocking mitochondrial apoptosis is consistent with TB models [37,38] showing SIRT1 overexpression promoting bacterial reservoirs by inhibiting infected cell death. 

RNA Polymerase II Process and TNF Feedback Loop Process

Herta et al. [31] support the link between KLF4 and KLF5 and the development of increased cytokines due to their opposing roles in inflammation and immune regulation. Wang et al. [39] support the increased attention in M. TB therapies to TNF receptors and FOS related genes suggesting and backing up their importance in the RNA Polymerase II process.

Implications: 

By identifying similarities between the infection of TB to the progression of latent to active TB, the genes and pathways identified in this study can serve as pre-72 hour biomarkers to stratify patients by risk of progression and be utilized in clinical tests. By doing so, clinical resources can be more effectively allocated to patients as excessive testing and repetitive screening can be avoided entirely and look to decrease global cases of TB as less patients will be infected. 

Limitations: 

In our case, since we use bioinformatics datasets from RNA-seq experiments conducted by other researchers, one limitation is that the identified genes will need to be further studied in the laboratory or clinical environment before developing any interventions. Additionally, the chosen dataset was of a low sample size in Ohio and Texas and only validated by a single follow up study in South Africa, which although proves statistical significance it loses its effectiveness due to the extremely high variability of TB infections on hosts. Of course, the proposed mechanisms are hypothetical and have limited clinical validation and require further testing to validate. Finally, the chosen dataset was of non-infected cells derived from donors versus infected cells. Therefore, the proposed mechanisms in this study are only valid if similar mechanisms can be found from infection of TB and progression of TB.  

Future Directions: 

Studies can be conducted to identify similar mechanisms in infection of TB and progression of TB, which when compared to the results of this study can find proposed mechanisms which work most effectively. The mechanisms of the proposed processes and biomarkers can be further researched to confirm the logistics of them. Then, clinical trials can be conducted to identify if the proposed biomarkers are conducive of the progression of TB. Specifically, we would recommend to focus on SIRT1 activity and the KLF5:KLF4 ratio as upstream controllers of the proposed mechanisms. They are detectable, mechanistically grounded, and most clinically actionable.