Novel Biomarker Panel of Let-7d-5p and MiR-140-5p Can Distinguish Latent Tuberculosis Infection from Active Tuberculosis Patients

Introduction

Tuberculosis is one of the leading causes of death worldwide. In 2021, 10.60 million people were newly infected with Mycobacterium tuberculosis (Mtb), according to the data released by World Health Organization (WHO).¹ China ranks third among countries with a high burden of tuberculosis in the world, with 780,000 new tuberculosis patients. Mtb can survive in macrophages of tuberculosis granuloma in host.² The dominant host immune response can prevent bacterial growth and eliminate Mtb. On the contrary, Mtb can remain in the host for a long time and become latent tuberculosis infection (LTBI). Only 5–15% of population with LTBI will progress to active tuberculosis (ATB).³ Tuberculin skin test (TST) and interferon-γ release assay (IGRA) are common methods for screening Mtb infection, but they cannot distinguish those people with LTBI from ATB patients.⁴

The lack of effective diagnosis of LTBI has largely hampered the development of clinical diagnosis and therapy of tuberculosis. The T cell-based interferon γ release assay(IGRA) has contributed to the diagnosis of tuberculosis, but it still cannot distinguish those people with LTBI from ATB patients. Wang et al developed a method to diagnose LTBI, adding IL-2 as an additional biomarker, together with IFN-γ to distinguish those people with LTBI from those with ATB to enhance the diagnostic performance (area under the curve [AUC] of receiver operating characteristics [ROC] analysis is 0.7494 [95% CI: 0.6419–0.8569, P < 0.0001]).⁵ There are also lots of studies on diagnostic markers based on blood samples, miRNA, cytokines and proteins in peripheral blood have the potential to be used as diagnostic biomarkers.^6–8

MicroRNA is an endogenous single stranded noncoding small RNA with length of 18–22 nucleotides. It participates in the pathogenesis of various diseases by triggering the translation inhibition or degradation of its target mRNA.^9–11 MicroRNA plays a key role in the progress of tuberculosis. After tuberculosis infection, a large number of miRNA expression are up-regulated or down-regulated to regulate apoptosis,¹² autophagy, inflammation and innate immune responses.^13–15 More and more evidences show that miRNA may become a new tool for the diagnosis of tuberculosis.¹⁶

The stability of miRNA makes it as a good candidate with high potential to be used as biomarker. More and more evidences show that miRNAs play an important role in regulating tuberculosis infection. Zhu et al found that the down-regulation of miR-18b leads to the expression increasing of its target gene HIF-1α and promotes the inflammatory response, clearing Mtb inside of the macrophages.¹⁴ Zhu et al found that miR-378d was down regulated during the infection with Mtb, and the intracellular survival rate of Mtb was reduced by targeting Rab10.¹⁷ Yuan et al released that the overexpression of mir-18a down regulated its target gene ATM, and then promoted the survival of intracellular Mtb via the ATM pathway and autophagy process.¹⁸

In this study, we screened the data sets published in GEO, conducted joint differential analysis and target gene enrichment analysis, after filtering the disease-related database, we screened the differential miRNA related to TB, verified expression of real-time fluorescence quantitative PCR (qPCR), and validated the biomarkers for their performance for distinguishing populations with LTBI from ATB patients.

Methods and Materials

Data and Samples Collection

‘Tuberculosis’ was used as a keyword to search data sets in GEO database. Inclusion criterion of these data are as follows:

1. Species should be human.
2. Sample type should be peripheral blood mononuclear cell (PBMC) or whole blood or serum or plasma.
3. Each data set contains ATB, LTBI and Healthy Control (HC) three groups.
4. Data set should be normalized.

Three data sets, including two miRNA data sets and one mRNA data set, were downloaded from GEO databases and information in detail about each data set can be found in Table 1. All miRNA data sets in this study were openly available in GEO database (NCBI) and downloaded by searching GEO ID. Two miRNA target public databases including mirtarbase (https://mirtarbase.cuhk.edu.cn/) and mirDB (http://mirdb.org/) were also downloaded for analysis (Figure 1).

Table 1 Public Data sets

Figure 1 The study design and data analysis workflow.

The study was performed using samples stored at the First Affiliated Hospital of University of Science and Technology of China (Anhui Provincial Hospital Infection Hospital).The diagnosis of ATB was based on clinical presentation, chest radiography, and bacterial plating culture results. Furthermore, IGRA results used to distinguish Mtb group and HC group. Inclusion criteria: 1) Active pulmonary tuberculosis, the diagnostic criteria for patients with latent tuberculosis were referred to China 2001 Guidelines for the diagnosis and treatment of tuberculosis developed by the Tuberculosis Branch of the Medical Association; 2) Patients with tuberculosis through sputum smear acid-fast bacillus or tuberculosis culture positive diagnosis; 3) The healthy subjects underwent routine physical examination and PPD test antigen specificity test and Mtb test were negative; 4) Patients with latent infection were specific for PPD tests and Mtb antigens test is positive, but the patient has no TB-related symptoms; 5) Age range from 18- to 69-years-old; 6) Informed consent of patients was obtained in this study. Exclusion criteria: 1) lung tumor disease; 2) Pulmonary fibrosis, silicosis, pneumoconiosis; 3) Patients infected with human immunodeficiency virus (HIV); and 4) along with other types of infectious diseases.The basic clinical characteristics of the participants were summarized in Table 2. A total of 84 serum samples were collected from 29 ATB patients, 25 LTBI patients and 30 HCs.

Table 2 Characteristics of Study Participants with ATB, LTBI and HCs

MiRNA Differential Analysis

Limma (3.46.0) was used to perform differential analysis for each miRNA data set. P value (p<0.01) and fold change (|log2foldchange|≥1) were used to determine the differentially expressed miRNAs (DEMs) among the three groups of interest. DEMs results from two miRNA data sets were combined as total DEMs in each comparison group including ATB vs. LTBI and LTBI vs. HC.

Differential miRNA Target Genes Prediction

The mirtarbase contains miRNA-target interactions (MTIs) information. The collected MTIs are experimentally validated by reporters’ studies, such as Western blotting or microarray assays with overexpression or knockdown of miRNAs. And mirDB is an online database for miRNA target prediction and functional annotations. After selecting DEMs, this two databases were used to search target genes of each DEM. Furthermore, another mRNA public data set was used to filter DEMs’ target genes based on the theory that most miRNAs could down-regulated target genes.

Limma package was also used to analyze the differential expression of mRNA with p value (p<0.01) and fold change (|log2foldchange|≥1). Combined with MTI information, interactions with opposite expression trend were selected, which meant miRNA was up-regulated while mRNA down-regulated, vice versa.

Enrichment of Target Genes

In this study, metascape, a web enrichment tool was used to perform enrichment analysis through default databases and default parameters. After that, we intersected ATB vs. LTBI and LTBI vs. HC and determined the top 20 of enrichment results to find vital enrichment results in tuberculosis infection process.

MiRNA Filtering According to Enrichment Result

In order to research the roles of miRNAs in enrichment functions, single miRNA impact score (S_impact) was defined to estimate the importance of miRNAs in specific functional term (see Equation 1.1). MiRNAs with single impact score higher than 0.01 were selected. Furthermore, according to HMDD (the human miRNA disease database), only miRNAs related to tuberculosis were select to further analysis.

(1.1)

Note that nmt represents the number of miRNA target genes, n represents total number of miRNAs in specific enrichment term.

Serum Sample PCR Validation

MiRNA Extraction

In order to perform quantitative analysis of the target miRNAs, a final concentration of 10⁻⁴ pM commercial external control miRNAs (Cat: CR100-01, TIANGEN, Beijing, China) was equally added to each serum sample from ATB, LTBI and HC before miRNA extraction, which is a reference gene applied in the subsequent test.

Total miRNAs were extracted from patients’ serum using an miRNA isolation Kit (Cat: DP503, TIANGEN, Beijing, China). Briefly, 900μL lysis solution was add to 200μL serum and mixed given its two solutions sufficiently, then mixture sample was incubated for 5 min at room temperature to permit the complete disruption of nucleoprotein complexes. The tube containing the homogenate and 0.2mL chloroform was shaken vigorously to facilitate the denaturation and precipitation of proteins. After centrifugation at 12,000×g for 15 min at 4°C, the raw RNA in supernatant was transferred in aqueous to a fresh tube and purified following the manufacturer’s instruction.

Candidate miRNA Quantification by Real-Time qRT-PCR

TIANGEN miRcute Plus miRNA qPCR Kit (Cat: FP411-01, TIANGEN, Beijing, China) was used to determine the expression levels of miRNAs according to the manufacturer’s instruction in a CFX ConnectTm Real-Time System (Bio-Rad, Hercules, CA, USA).

The first-strand cDNA was synthesized from the total RNA using a miRcute Plus miRNA First-Strand cDNA Kit (Cat: KR211-02, TIANGEN, Beijing, China). Subsequently, quantitative real-time PCR was conducted using gene-specific primers to measure the expression levels of miRNAs. Gene-specific primers used in this assay were listed in Table 3. Following initial denaturation at 95°C for 30s, 40 cycles of PCR amplification were performed at 95°C for 5s and 60°C for 30s with a dissociation curve at the end of the amplification reaction.

Table 3 Primers Used in the qPCR Assay

The obtained qRT-PCR data were calculated using the 2^−ΔΔCT method¹³ while the relative amount of the target gene miRNAs were normalized according to the difference in CT values of commercial external control miRNA.

The results of qRT-PCR were analyzed using GraphPad Prism (version 7.0). Differences between groups were statistically analyzed by analysis of Mann–Whitney-Wilcoxon. A P-value <0.05 was considered statistically significant.

Model Construction to Identify miRNA Marker

We evaluated different combinations of biomarkers by logistic regression to obtain a biomarker panel with good performance in diagnosing LTBI. LTBI samples as a whole was labeled as “1” and all HC+ATB samples as a control group with label “0”. Data were randomly divided into two sets for each of the three comparisons, with one set serving as a training set and the other as a test set. The training set was defined as the model used for training and test set was an independent data set to verify the validity of the model. A stepwise logistic regression model was used to select diagnostic microRNA markers based on the training data set and then checked by leave-one-out cross-validation.¹⁴ A test set of samples was used to verify the first set using stepwise logistic regression and the model established was then evaluated using a ROC curve.

Results

Differentially Expressed miRNAs Analysis

To identify DEMs in ATB vs. LTBI and LTBI vs. HC, P-value < 0.05 and |log2FC|≥1 were set as cut-off criteria. A total of 174 DEMs were identified among the two groups in GSE25435 data set including 51 upregulated miRNAs and 39 downregulated miRNAs in ATB vs. LTBI. However, only 24 DEMs were identified in GSE29190 data set (Table 4). The difference in the number of DEM screening between the two data sets may be due to the inconsistency of their experimental platform versions, where GES2435 is Agilent miRNA microarray (version 2.0) and GSE29190 is Agilent miRNA microarray (version 3). DEMs results were demonstrated by using volcano plots (Figure 2). 101 DEMs including 57 upregulated miRNAs and 44 downregulated miRNAs were found between ATB and LTBI after combined two data set results.97 DEMs including 43 upregulated miRNAs and 54 downregulated miRNAs were also identified between LTBI and HC by the same way.

Table 4 Two Data sets miRNA Differential Analysis Results

Figure 2 Volcano plot of differentially expressed miRNAs. (a) GSE25435 miRNA data set ATB vs. LTBI differential analysis; (b) GSE25435 miRNA data set LTBI vs. HC differential analysis; (c) GSE29190 miRNA data set ATB vs. LTBI differential analysis; (d) GSE29190 miRNA data set LTBI vs. HC differential analysis. Note that only top 10 up and down differential miRNAs were marked in the plots.

Differentially Expressed miRNAs’ Target Genes Prediction

MiRNAs were important in different biological functions and pathways by regulating their target mRNAs. In this section, according to the MTI information and difference analysis results of mRNA data set, the most important miRNA-mRNA association was selected to further analyze which functions would affected by miRNAs. Totally, 86,781 paired miRNA-target gene including 2459 miRNAs and 8849 target genes were obtained from both mirtarbase and mirDB. MiRNAs could downregulate target genes due to the biological process that miRNAs bind to the 3’UTRs of their target genes and interfere with translation.¹⁹ 2156 target genes regulated by 89 miRNAs between ATB and LTBI, and 2139 target genes regulated by 81 miRNAs between LTBI and HC were retained after removing miRNAs which upregulate mRNAs. Finally, the selected target genes were used to next function enrichment analysis.

Metascape Enrichment Analysis

The top 20 enriched functional terms of differentially expressed mRNA in ATB vs. LTBI were shown in Figure 3. Among that, response to growth factor (GO:0070848), pathways in cancer (hsa05200), blood vessel development (GO:0001568) and regulation of kinase activity (GO:0043549) were most significant. For LTBI vs. HC, Figure 4 was also shown the top 20 enriched functional terms. Some enriched terms such as response to growth factor, signaling by receptor tyrosine kinases (R-HSA-9006934), diseases of signal transduction by growth factor receptors and second messengers (R-HSA-5663202), pathways in cancer were identified. Some enriched terms were also found in previous researches.²⁰ Pathways in cancer enriched term was found in both ATB vs. LTBI and LTBI vs. HC. Microbes caused diseases like Hepatitis B virus could lead to liver cancer. Previous research found that chronic tuberculosis infection could cause DNA damage in their vicinity.²¹

Figure 3 Top 20 enrichment terms in ATB vs. LTBI differentially expressed miRNA target genes enrichment terms. The value at the end of the bar indicates the number of genes enriched to this function.

Figure 4 Top 20 enrichment terms in LTBI vs. HC differentially expressed miRNA target gene enrichment terms. The value at the end of the bar indicates the number of genes enriched to this function.

Response to growth factor was most enriched term in both ATB vs. LTBI and LTBI vs. HC. This term may play an important role in the development of tuberculosis infection. Some growth factors like vascular endothelial growth factor(VEGF) are increased in ATB patients.²² Previous study also showed that deletion of TGF-β1 from the systems resulted in improved bacterial clearance and lesion sterilization.²³ Furthermore, 181 genes regulated by 57 miRNAs were enriched in response to growth factor in ATB vs. LTBI. And 170 genes regulated by 60 miRNAs were also enriched in same term.

Vital miRNAs Filtered by Impact Score

In order to find vital miRNAs which can play important roles in most enriched terms. We defined miRNAs’ impact score to estimate miRNA impact value in specific term. Two comparisons’ (ATB vs. LTBI and LTBI vs. HC) most enriched term shared 11 miRNAs including hsa-let-7d-5p, hsa-miR-140-5p, hsa-miR-155-5p, hsa-miR-194-5p, hsa-miR-19a-3p, hsa-miR-29b-3p, hsa-miR-29c-3p, hsa-miR-30e-5p, hsa-miR-32-5p, hsa-miR-374a-5p and hsa-miR-92a-3p after removing miRNAs that impact score less than 0.01. In order to further screen out more reliable biomarker combinations, miRNAs were required to have a known relationship to tuberculosis. Thus, only 4 miRNAs including hsa-let-7d-5p, hsa-miR-140-5p, hsa-miR-155-5p, hsa-miR-32-5p were left after removing miRNAs not related to tuberculosis infection from HMDD database. The expression of four miRNAs in the three groups were validated by qPCR in serum samples from the validation/test cohort. (Figure 5).

Figure 5 PCR results for four miRNAs. Notes: (a) serum hsa-let-7d-5p expression levels in three groups; (b) serum hsa-miR-155-5p expression levels in three groups; (c) serum hsa-miR140-5p expression levels in three groups; (d) serum hsa-miR-32-5p expression levels in three groups. The expression level of the validated miRNAs were normalized to CR100-01 (Log 2 relative level). Differences between groups were statistically analyzed by analysis of Mann–Whitney-Wilcoxon. *P<0.05, **P<0.01

MiRNAs Classification Model Building

Before model building, the assumptions of Logistic Regression was checked. First, as mentioned in method, the output label was binary and all samples were independent. The Pearson correlation coefficient between hsa-let-7d-5p and hsa-miR-140-5p was 0.2418 (P-value = 0.034). This indicates a low correlation between the two independent variables. Result of Box-Tidwell test shown that hsa-let-7d-5p (P = 0.9301) and hsa-miR-140-5p (P = 0.6186) have the linearity with log(odd) and that the assumption is satisfied. Models for the training sets were built using the logistic regression models: logit(P = LTBI) = −1.065+0.138*(hsa-let-7d-5p) −0.131 *(hsa-miR-140-5p).The diagnostic value of the data models derived from the training sets was validated with data from the test sets by examining ROC curves. The AUC value of the ROC curve was 0.930 (sensitivity = 100%, specificity = 88.5%) and 0.923 (sensitivity = 100%, specificity = 92.3%) with the biomarker panel for the training set and test set, respectively (Figure 6 and Table 5).

Table 5 The AUC, Specificity and Sensitivity of Two miRNA Panel

Figure 6 ROC curves for the biomarker panel of data sets. ROC curves of training set (a) showed AUC value was 0.93, the test set (b) showed AUC value was 0.923.

Discussion

According to WHO estimation, the global TB latent infection population is close to 2 billion. The identification of LTBI is very important for TB prevention and control worldwide.¹ Therefore, early detection and early treatment of LTBI is very necessary.

Biomarkers from peripheral blood are the best choice for clinical application due to their ease of collection and subsequent testing, which is also the main consideration in our selecting data sets. Here, we used a new filtration screening method to find serological diagnostic markers that can well distinguish the population with LTBI from those with pulmonary TB. Firstly, the miRNA expression of LTBI, ATB and HC from GEO data set were analyzed to select the DEMs. After filtering by miRNA impact score and disease database, 4 miRNAs were regarded as candidate biomarkers finally. After that, the combine panel of let-7d-5p and miR-140-5p obtained the lowest AIC in logistic regression model. ROC results of LOOCV shown that our model trained by let-7d-5p and miR-140-5p expression was good at distinguishing LTBI from ATB or HC with AUC higher than 0.9. As a result, we successfully discovered two reliable LTBI miRNA biomarkers.

Moreover, MiR-140-5p and let-7d-5p are not only important diagnostic markers for tuberculosis in this study. They are also widely studied in other diseases. In other studies, the levels of miR-140-5p were decreased in peripheral blood of children with Necrotizing pneumonia (NP), and miR-140-5p regulated the expressions of TLR4 in A549 cells.²⁴ Moreover, Yang et al has showed that miR-140-5p regulates the inflammation of Acute Lung Injury (ALI) via TLR4/MyD88/NF-κB signaling pathway.²⁵ Zhuo et al assessed that miR-140-5p inhibited cell proliferation and induced apoptosis of lung cancer cells and hindered the migration and invasion of lung carcinoma cells.²⁶ In addition, Dong et al reported that miR-140-5p inhibited pulmonary fibrosis by blocking the activity of the Wnt1/β-catenin signaling pathway.²⁷

As shown in study of Li et al, let-7d-5p was up-regulated in A549 cell exosomes, and participated in the generation and maintenance of cancerous bone pain in vivo.²⁸ Let-7d-5p was also up-regulated in isolated pulmonary presentation (PSP-BHD), and inhibited the repair response of lung lesions.²⁹ There was study shown that let-7d-5p may play its inhibitory role in the inflammatory response and apoptosis of intestinal epithelial cells in neonatal rats with necrotizing enterocolitis through the TLR4/NF-κB signaling pathway.^30,31 Let-7d-5p expression was decreased in the plasma of patients with acute exacerbation of idiopathic pulmonary fibrosis (AE-IPF), indicating that let-7d-5p has potential as a marker for diagnosis and prediction of disease progression.¹⁰

However, the roles of miR-140-5p and let-7d-5p in the pathogenesis of tuberculosis have not been elucidated. In the present study, the expression of them were down-regulated in the serum of patients compared with HCs. And their target genes were enriched in terms of growth factor response, vascular development and kinase activity regulation. Whether they were involved in the pathogenesis of tuberculosis by regulating inflammation, apoptosis and other processes remains to be further explored.

For most marker identification methods, differentially expressed miRNAs or genes were used to identify markers. In this research, we selected miRNA markers to distinct different groups by combining miRNA functions. A new method was used to estimate miRNAs’ importance in specific pathway or enriched term. This method could be applied to other researches which aim to identify most important miRNAs in specific function terms.

Our research aimed at exploring biomarkers from peripheral blood, which are considered as highly valuable biomarkers and are conducive to large-scale further application. However, the sample size used in our research is insufficient, it is more necessary to use a larger number of samples to verify these targets. In addition, more subgroups such as patients with anti-TB therapies and cured TB need to be included in the study for further evaluation.

Conclusion

Response to growth factor was the most enriched term in both ATB vs. LTBI and LTBI vs. HC, which indicated that this term may have essential impact on the whole progression of tuberculosis infection. In our study, we coined S_impact to estimate miRNA values in specific function terms. This could help researchers to filter out important miRNAs. We believe that miRNAs’ S_impact which larger than 0.01 were important miRNAs. In this study, HMDD database was also considered to remove unrelated to tuberculosis infection miRNAs. It showed that the logistic regression model built by let-7d-5p and miR-140-5p have the ability to distinguish LTBI from active TB patients. However, our model still has some limitations. At present, our research results have only been verified in a small-scale sample, and its reliability has not yet been verified by a large-scale sample size, so it cannot be generalized for the time being. In addition, since the exact clinical features of LTBI have not yet been well-defined, there may be more influencing factors in the selection of LTBI samples that lead to changes in expression data. Although there are some limitations, our research can also provide strong data and theoretical support for future in-depth research and the development of tuberculosis early diagnosis products. The screening of biomarkers and their combinations is a huge and continuous research work. Although better miRNA combinations are screened at present, we believe that the combination of miRNA and other markers (such as mRNA, protein, etc.) will bring new breakthroughs in this field in the future. Finally, stable and reliable biomarker panels can also provide quantitative standards in the clinical definition of LTBI.