Advertisement

Comprehensive identification of immuno-related transcriptional signature for active pulmonary tuberculosis by integrated analysis of array and single cell RNA-seq

  • Yuzhong Xu
    Affiliations
    Department of Clinical Laboratory, Shenzhen Baoan Hospital, The Second Affiliated Hospital of Shenzhen University, Shenzhen 518101, China
    Search for articles by this author
  • Yaoju Tan
    Affiliations
    Department of Clinical Laboratory, Guangzhou Chest Hospital, Guangzhou/State Key Laboratory of Respiratory Disease, Hengzhigang Road 1066, Guangzhou 510095, China
    Search for articles by this author
  • Xianyi Zhang
    Affiliations
    Department of Clinical Laboratory, Shenzhen Baoan Hospital, The Second Affiliated Hospital of Shenzhen University, Shenzhen 518101, China
    Search for articles by this author
  • Minggang Cheng
    Affiliations
    Department of Clinical Laboratory, Shenzhen Baoan Hospital, The Second Affiliated Hospital of Shenzhen University, Shenzhen 518101, China
    Search for articles by this author
  • Jinxing Hu
    Affiliations
    Department of Clinical Laboratory, Guangzhou Chest Hospital, Guangzhou/State Key Laboratory of Respiratory Disease, Hengzhigang Road 1066, Guangzhou 510095, China
    Search for articles by this author
  • Jianxiong Liu
    Affiliations
    Department of Clinical Laboratory, Guangzhou Chest Hospital, Guangzhou/State Key Laboratory of Respiratory Disease, Hengzhigang Road 1066, Guangzhou 510095, China
    Search for articles by this author
  • Xinchun Chen
    Affiliations
    Guangdong Provincial Key Laboratory of Regional Immunity and Diseases, Department of Pathogen Biology, Shenzhen University School of Medicine, Shenzhen 518060, China
    Search for articles by this author
  • Jialou Zhu
    Correspondence
    Corresponding author.
    Affiliations
    Department of Clinical Laboratory, Guangzhou Chest Hospital, Guangzhou/State Key Laboratory of Respiratory Disease, Hengzhigang Road 1066, Guangzhou 510095, China

    Guangdong Provincial Key Laboratory of Regional Immunity and Diseases, Department of Pathogen Biology, Shenzhen University School of Medicine, Shenzhen 518060, China
    Search for articles by this author
Open AccessPublished:August 21, 2022DOI:https://doi.org/10.1016/j.jinf.2022.08.017

      Highlights

      • Comprehensive illustration of immuno-related transcriptional signature for active pulmonary TB.
      • Immune transcriptional profiling could effectively differ TB with LTBI and HC.
      • First scRNA-seq-based description of ADM expression abundance in TB Myeloid cells.
      • ADM could be used as a novel biomarker for distinguishing patients with TB from LTBI and HC.

      Abstract

      Background

      Tuberculosis (TB) continues to be a major cause of morbidity and mortality worldwide. However, the molecular mechanism underlying immune response to human infection with Mycobacterium tuberculosis (Mtb) remains unclear. Assessing changes in transcript abundance in blood between health and disease on a genome-wide scale affords a comprehensive view of the impact of Mtb infection on the host defense and a reliable way to identify novel TB biomarkers.

      Methods

      We combined expression profiling by array and single cell RNA-sequencing (scRNA-seq) via 10X Genomics platform to better illustrate the immuno-related transcriptional signature of TB and explore potential diagnostic markers for differentiating TB from latent tuberculosis infection (LTBI) and healthy control (HC).

      Findings

      Pathway analysis based on differential expressed genes (DEGs) revealed that immune transcriptional profiling could effectively differ TB with LTBI and HC. Following WGCNA and PPI network analysis based on DEGs, we screened out three key immuno-related hub genes (ADM, IFIT3 and SERPING1) highly associated with TB. Further validation found only ADM expression significantly increased in TB patients in both adult and children's datasets. By comparing the scRNA-seq datasets from TB, LTBI and HC, we observed a remarkable elevated expression level and proportion of ADM in TB Myeloid cells, further supporting that ADM expression changes could distinguish patients with TB from LTBI and HC. Besides, the hsa-miR-24–3p-NEAT1-ADM-CEBPB regulation pathway might be one of the critical networks regulating the pathogenesis of TB. Although further investigation in a larger cohort is warranted, we provide useful and novel insight to explore the potential candidate genes for TB diagnosis and intervention.

      Interpretation

      We propose that the expression of ADM in peripheral blood could be used as a novel biomarker for differentiating TB with LTBI and HC.

      Keywords

      Introduction

      Tuberculosis (TB), a chronic infectious disease caused by Mycobacterium tuberculosis (Mtb), is a major global health threat responsible for ∼10 million new cases and 1.3 million deaths in 2020

      World Health O. Global tuberculosis report 2021. Geneva: World Health Organization; 2021 2021.

      . Until the COVID-19 pandemic, TB was the leading cause of death from a single infectious agent, ranking above HIV/AIDS. About one-quarter of the world's population has been infected with the infectious pathogen Mtb. Approximately a small proportion (5–15%) of infected individuals progress to active TB due to failed host immune control that contains the bacterial spread, while the vast majority remain clinically asymptomatic in a state known as latent TB infection (LTBI)
      • Luo Y.
      • Suliman S.
      • Asgari S.
      • Amariuta T.
      • Baglaenko Y.
      • Martinez-Bonet M.
      • et al.
      Early progression to active tuberculosis is a highly heritable trait driven by 3q23 in Peruvians.
      . Specially, active TB can develop in certain conditions after many years of latency.
      At present, the “bacteriologically confirmed” diagnostic tests for TB, mainly based on culture, rapid molecular tests recommended by WHO, or sputum smear microscopy, face major challenges that high proportions of active TB with negative bacteriology and relatively low sensitivity of such test. Immunological tests, such as tuberculin skin test (TST) or interferon-gamma release assays (IGRAs), established by demonstrating an immune response to Mtb infection, do not require sputum or pulmonary samples. Unfortunately, both the two tests can neither distinguish LTBI from active TB nor the different stages of Mtb infection
      • Dheda K.
      • Barry 3rd, C.E.
      • Maartens G.
      Tuberculosis.
      . Therefore, the development of sputum-free, easy-to-operate, high sensitivity and specificity strategies to assist in identifying active and latent tuberculosis is crucial for TB management and control.
      Host immunological resistance to persistent Mtb infection, including innate and adaptive immune response, has been reported to play vital roles in controlling Mtb infection or its clearance, which determines whether the individual is asymptomatic or ill
      • Simmons J.D.
      • Stein C.M.
      • Seshadri C.
      • Campo M.
      • Alter G.
      • Fortune S.
      • et al.
      Immunological mechanisms of human resistance to persistent Mycobacterium tuberculosis infection.
      . Blood, which stores a history of the immune system, is arguably the most readily accessible biofluid in humans. Assessing changes in transcript abundance in blood between health and disease on a genome-wide scale affords a comprehensive view of the impact of Mtb infection on the host defense and a reliable way to identify novel TB biomarkers. Indeed, transcriptomic analysis of whole blood has provided a powerful tool to discover key factors in TB pathogenesis
      • Kaforou M.
      • Wright V.J.
      • Oni T.
      • French N.
      • Anderson S.T.
      • Bangani N.
      • et al.
      Detection of tuberculosis in HIV-infected and -uninfected African adults using whole blood RNA expression signatures: a case-control study.
      • Sweeney T.E.
      • Braviak L.
      • Tato C.M.
      • Khatri P.
      Genome-wide expression for diagnosis of pulmonary tuberculosis: a multicohort analysis.
      • Berry M.P.
      • Graham C.M.
      • McNab F.W.
      • Xu Z.
      • Bloch S.A.
      • Oni T.
      • et al.
      An interferon-inducible neutrophil-driven blood transcriptional signature in human tuberculosis.
      .
      While bulk RNA-seq or microarray reflects only the divergence of the most abundant cell populations, or moderate changes as a result of a cell medley, scRNA-seq is considered a powerful tool for resolving cellular heterogeneity in gene expression patterns in peripheral blood cell populations in health and disease
      • Papalexi E.
      • Satija R.
      Single-cell RNA sequencing to explore immune cell heterogeneity.
      • Deng Y.
      • Zheng Y.
      • Li D.
      • Hong Q.
      • Zhang M.
      • Li Q.
      • et al.
      Expression characteristics of interferon-stimulated genes and possible regulatory mechanisms in lupus patients using transcriptomics analyses.
      • Cai Y.
      • Dai Y.
      • Wang Y.
      • Yang Q.
      • Guo J.
      • Wei C.
      • et al.
      Single-cell transcriptomics of blood reveals a natural killer cell subset depletion in tuberculosis.
      . Here, we integrated expression profiling by array with scRNA-seq from TB, LTBI and HC to comprehensively delineate the immuno-related transcriptional signature of TB and explore potential diagnostic markers for differentiating TB from LTBI and HC.

      Material and methods

      Dataset collection and processing

      The gene expression profiles of GSE28623
      • Maertzdorf J.
      • Ota M.
      • Repsilber D.
      • Mollenkopf H.J.
      • Weiner J.
      • Hill P.C.
      • et al.
      Functional correlations of pathogenesis-driven gene expression signatures in tuberculosis.
      , GSE41055
      • Verhagen L.M.
      • Zomer A.
      • Maes M.
      • Villalba J.A.
      • Del Nogal B.
      • Eleveld M.
      • et al.
      A predictive signature gene set for discriminating active from latent tuberculosis in Warao Amerindian children.
      and GSE19439
      • Berry M.P.
      • Graham C.M.
      • McNab F.W.
      • Xu Z.
      • Bloch S.A.
      • Oni T.
      • et al.
      An interferon-inducible neutrophil-driven blood transcriptional signature in human tuberculosis.
      were downloaded from Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo). In brief, 46 TB, 25 LTBI and 37 HC whole blood samples from adults were included in GSE28623, and 9 TB, 9 LTBI and 9 HC whole blood samples from children were analyzed in GSE41055, and 13 TB, 17 LTBI and 12 HC whole blood samples from adults were included in GSE28623. For 10x genomics scRNA-seq, 3 TB, 2 LTBI and 2 HC peripheral blood mononuclear cell (PBMC) samples were included. General characteristics of the study subjects are given in Table S1. The overall research design is shown in Fig. 1. This study was approved by the Institutional Ethics Board and informed consent was obtained from each participant.

      Identification of DEGs and gene ontology enrichment analysis

      The R package limma
      • Ritchie M.E.
      • Phipson B.
      • Wu D.
      • Hu Y.
      • Law C.W.
      • Shi W.
      • et al.
      limma powers differential expression analyses for RNA-sequencing and microarray studies.
      was applied to identify the DEGs among TB, LTBI and HC in GSE28623, with the filtering threshold setting as adjusted p-value ≤0.05 and |log2 fold change (FC) | ≥1. Then the R package ClusterProfiler
      • Yu G.
      • Wang L.G.
      • Han Y.
      • He Q.Y.
      clusterProfiler: an R package for comparing biological themes among gene clusters.
      was exploited to perform Gene Ontology (GO) enrichment analysis based on the DEGs between TB and LTBI or HC. Enriched GO pathways were determined according to the cut-off criterion of adjusted P-value ≤ 0.05.

      Consensus clustering analysis

      Human immuno-related gene sets were obtained after intersecting the Immprot and InnateDB databases. The R package ConsensusClusterPlus
      • Wilkerson M.D.
      • Hayes D.N.
      ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking.
      was then used to perform consensus clustering analysis to investigate the association between immuno-related genes and the disease state, with the maximum number of clusters set at 8.

      Weighted gene co-expression network analysis

      Weighted Gene Co-Expression Network Analysis (WGCNA) is a well-established method for studying biological networks and diseases. Thus, WGCNA was performed with the R package WGCNA
      • Langfelder P.
      • Horvath S.
      WGCNA: an R package for weighted correlation network analysis.
      to identify key modules related to TB by using immuno-related genes within GSE28623 dataset. First, various β parameters were assessed to construct a scale-free network, and we selected the power of β=10 as the soft thresholding value. Second, the dynamic tree cut algorithm was applied with the minimal gene number of 30 within each module. Third, similar modules were merged through dynamic hybrid cutting at a threshold of 0.05. Finally, trait-related modules were selected based on the Pearson correlation coefficient between ME and sample characteristics of each module.

      Protein-protein interaction network construction

      The Metascape (metascape.org) database was used to construct a protein-protein interaction (PPI) network for co-expressed genes within the modules. Then, the Cytoscape software (3.6.0) (http://cytoscape.org) was applied to visualize the network interactions. Hub genes were defined as genes with top 10 ranked connectivity degrees. In brief, the node fill color was set based on the hub gene. Yellow represents hub genes. The node height and width were set based on the interaction degree with other genes. The interaction network layout was set to circular.

      Single-cell data preprocessing, gene expression quantification and cell-type determination

      Single cell raw sequencing data (10X Genomics) was obtained from our previous study under the accession number (SRP247583) of NCBI Short Read Archive (SRA)
      • Cai Y.
      • Dai Y.
      • Wang Y.
      • Yang Q.
      • Guo J.
      • Wei C.
      • et al.
      Single-cell transcriptomics of blood reveals a natural killer cell subset depletion in tuberculosis.
      . Raw sequencing reads were processed by Cell Ranger (6.1.2) and aligned to the human reference genome (GRCh38). The unique molecular identifier (UMI) count matrices were then imported into the R package Seurat (4.1.0). Cells expressing <200 or >6000 genes or a high mitochondrial transcript ratio (> 0.2), and genes expressed in less than 3 cells were all excluded. Doublets were also checked and removed with the R package scDblFinder (1.8.0)
      • Germain P.L.
      • Sonrel A.
      • Robinson M.D.
      pipeComp, a general framework for the evaluation of computational pipelines, reveals performant single cell RNA-seq preprocessing tools.
      . For the remaining cells, gene expression matrices were normalized using NormalizeData function in Seurat package
      • Stuart T.
      • Butler A.
      • Hoffman P.
      • Hafemeister C.
      • Papalexi E.
      • Mauck 3rd, W.M.
      • et al.
      Comprehensive Integration of Single-Cell Data.
      , and Seurat FindVariableFeatures were applied to select the top 3000 genes exhibiting the highest cell-to-cell variation. Gene expression matrices from different samples were then integrated and the batch effects were removed by canonical correlation analysis and mutual nearest neighbors-anchors using Seurat package. Subsequently, the integrated, batch-corrected expression matrix for all cells was scaled with the Seurat ScaleData function to apply a linear transformation. Principle component analysis (PCA) was performed and Uniform Manifold Approximation and Projection (UMAP) were used for dimensionality reduction. Before clustering the cells, a shared nearest neighbor graph based on the Euclidean distance in PCA space was conducted using Seurat FindNeighbors. Clustering was then performed with Seurat FindClusters based on Louvain algorithm. Marker genes for each cluster were determined with the Wilcoxon rank-sum test using Seurat FindAllMarkers.

      Construction of miRNA-lncRNA-mRNA-TF gene regulatory network

      In brief, miRNA-mRNA and miRNA-lncRNA targeting regulation information were retrieved from starBase v3.0 (https://starbase.sysu.edu.cn) and the miRNA-mRNA interaction pairs should be included in at least four software or databases including PITA, RNA22, miRmap, DIANA-microT, miRanda, PicTar and TargetScan. The TF-mRNA interaction pairs were obtained from TRRUST v2.0 (https://www.grnpedia.org/trrust). After importing all screened mRNAs, miRNAs, lncRNAs and transcription factors (TFs), Cytoscape was used to construct the miRNA-lncRNA-mRNA-TF regulatory network.

      Statistical analysis

      All statistical analysis was performed using R software (4.1.0). A difference with a p-value or adjusted p-value ≤0.05 was considered statistical significant. The predictive sensitivity and specificity of selected marker genes for TB progression were assessed with receiver operating characteristic (ROC) curves by R package
      • Robin X.
      • Turck N.
      • Hainard A.
      • Tiberti N.
      • Lisacek F.
      • Sanchez J.C.
      • et al.
      pROC: an open-source package for R and S+ to analyze and compare ROC curves.
      .

      Results

      Identification of DEGs in the TB patients

      According to the data retrieving and processing procedure mentioned in methods, the DEGs among 37 HC, 25 LTBI and 46 TB samples in the GSE28623 dataset were identified, respectively. We observed 71 DEGs, including 55 up-and 16 down-regulated genes, in TB compared with HC groups (Fig. 2A, Table S2). We also found 163 DEGs in TB compared with LTBI groups, with 141 genes up-regulated and 22 genes down-regulated, respectively (Fig. 2B, Table S3). Then, we sought the overlapping DEGs between the two parts. After intersecting the two sets of genes, a total of 64 DEGs associated with TB were obtained (Fig. 2C).
      Fig 2
      Fig. 2Identification of differentially expressed genes (DEGs) among TB, HC and LTBI in the GSE28623 dataset. (A) Volcano plot of DEGs between TB and HC; (B) Volcano plot of DEGs between TB and LTBI; (C) Venn diagram shows the comparison of deregulated genes between TB and HC or LTBI; (D) GO analysis of the DEGs among TB, LTBI and HC; (E) Consensus clustering matrix of immuno-related genes within GSE28623.
      To explore the potential roles of DEGs among TB, LTBI and HC, GO enrichment analysis was performed. GO analysis showed that the up-regulated genes between TB and HC were mainly involved in biological processes (BP) associated with response to molecule of bacterial origin, defense response to bacterium, humoral immune response and myeloid leukocyte activation (Fig. 2D). Down-regulated genes were linked to leukocyte activation involved in immune response, lymphocyte differentiation and T cell activation (Fig. 2D). Among the DEGs between TB and LTBI, GO analysis indicated that the genes primarily enriched in defense response to bacterium, leukocyte mediated immunity and regulation of lymphocyte differentiation (Fig. 2D). Of interest, the DEGs between TB and HC and between TB and LTBI showed quite a similar GO enrichment pattern.
      Since the DEGs among TB, LTBI and HC were mainly enriched in pathways involved in immune response, thus we focused on the immuno-related genes from the intersecting parts of Immprot and InnateDB database in TB. The association between immuno-related genes and the disease state was further investigated. Consensus clustering analysis showed that all the samples could be grouped into four clusters. Of these, cluster 3 and 4 were primarily related to TB, indicating that immuno-related genes could differ TB from HC or LTBI (Fig. 2E).

      Identification of key immuno-related hub genes related to TB by WGCNA and PPI network

      WGCNA has been demonstrated to be an effective method for discovering key co-expressed genes linked to the clinicopathological characteristics of interest
      • Langfelder P.
      • Horvath S.
      WGCNA: an R package for weighted correlation network analysis.
      . Thus, WGCNA was then performed to identify key modules related to TB by using immuno-related genes within GSE28623 dataset. The power of β=10 (scale-free R2 = 0.92) was selected as the soft thresholding value to construct a scale-free network (Fig. 3A, B). After constructing similar module clustering through dynamic hybrid cutting at a threshold of 0.05, a total of 13 modules were identified (Fig. 3C). The Pearson correlation coefficient between ME and sample characteristics of each module was calculated. The results showed that the blue module had the highest negative correlation to TB, and the pink and black modules were highly positively correlated to TB (Fig. 3D).
      Fig 3
      Fig. 3WGCNA and PPI networks to identify key immuno-related hub genes related to TB. (A) The effect of different soft-thresholding powers on scale Independence (left) and mean connectivity (right); (B) Evaluation of scale free topology at selected β-values; Left: a histogram of network connectivities; Right: a log-log plot of the same histogram; (C) The dendrogram clustered by Dynamic Tree Cut algorithm; (D) The heatmap profiling the correlations between module eigengenes and the disease state; (E) PPI network from blue module; (F) PPI network from black module; (G) PPI network from pink module. Yellow represents hub genes. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)
      Subsequently, we constructed PPI networks with all the genes in blue, pink and black modules identified to be highly correlated to TB in WGCNA analysis, respectively. Hub genes were defined as genes with top 10 ranked connectivity degrees and are represented by yellow nodes in Fig. 3. The results showed that FYN, PRKCQ, CARD11, LCK, MAP3K14, ITPR3, ZAP70, SIGIRR, SYTL1 and LAT were the top 10 hub genes in blue module (Fig. 3E); ADM, FCGR2A, MAPK14, PLXNC1, TLR4, NAMPT, AQP9, FPR2, NFIL3 and S100A9 were the top 10 hub genes in black module (Fig. 3F); TAP1, IFITM2, SECTM1, IFITM1, IFITM3, TRAFD1, IFIT3, SERPING1, TRIM21, IRF7 were the top 10 hub genes in pink module (Fig. 3G). Taken together, 30 key hub genes related to TB were identified.
      Next, we intersected the key hub genes in TB and 64 DEGs generated from expression difference analysis, and finally extracted 3 key immuno-related genes highly associated with TB for the following analysis, including ADM, IFIT3 and SERPING1.

      Validation of key immuno-related genes expression in children with TB

      To better explore the expression pattern of selected genes among different types of samples, the expression levels of ADM, IFIT3 and SERPING1 were exhibited in the adult dataset (GSE28623) and further validated in a children's dataset (GSE41055) and an adult's dataset (GSE19439) (Fig. 4). Compared with healthy controls (HC) or individuals with latent TB infection (LTBI), the mRNA expression level of ADM significantly increased in TB patients in both the adult and children's datasets (Fig. 4), while IFIT3 and SERPING1 only showed a difference in adult, but not in children (Fig. 4), advising that the ADM expression changes could distinguish patients with TB from LTBI and HC. Thus, ADM was chosen for the following analysis.
      Fig 4
      Fig. 4Expression levels of ADM, IFIT3 and SERPING1 in the GSE28623 dataset (A-B), GSE41055 dataset (C-D) and GSE19439 dataset (E-F).

      The expression pattern of the key immuno-related genes in scRNA-seq

      As single cell RNA-seq (scRNA-seq) provides an unbiased high-throughput approach to illustrate the transcriptomic profile of individual single cells or immune-cell subsets, we then utilized scRNA-seq on PBMCs derived from seven individuals, including two HC, two LTBI and three TB. In total, 63,685 cells (19,248 from HC, 15,999 from LTBI, and28438 from TB) were taken into analysis (Table S4). Three major cell clusters across all three groups were identified, including myeloid cells expressing CD14 and LYZ, T cells expressing CD3D and CD3E, and B cells expressing CD79A and MS4A1 (Fig. 5A, B; Figure S1A, B). Subsequently, the expression levels and locations of the three key immuno-related genes ADM, IFIT3 and SERPING1 were determined. ADM and SERPING1 were dominantly expressed in Myeloid cells that nearly none of them were expressed in T and B cells, while IFIT3 was primarily expressed in Myeloid cells, but also expressed in T and B cells, though in a much lower percentage compared to Myeloid cells (Fig. 5C, D). Furthermore, when considering the disease state of each sample, a significantly elevated expression level and percentage of ADM was observed in TB compared to HC or LTBI in Myeloid cells (Fig. 5E, F; Figure S1C, D). In contrast, the expression level or percentage of IFIT3 significantly decreased in TB myeloid cells (Fig. 5E, Figure S1D). Besides, there was no marked difference regarding the expression status of SERPING1 between TB and LTBI or HC (Fig. 5E, Figure S1D). Taken together, our scRNA-seq data revealed ADM dominantly expressed in TB Myeloid cells at single cell resolution, which was in agreement with the above findings. ADM might serve as a useful diagnostic marker for differentiating TB from LTBI and HC.
      Fig 5
      Fig. 5Single-cell transcriptomic characterization of PBMCs from HC, LTBI and TB. (A) UMAP visualization of major cell types across all three groups. (B) Violin plots of classical markers defining each cell type. (C) Expression of the key immuno-related genes. (D) Dot plot showing the key immuno-related genes in each cell type. Dot size represents the proportion of cells within the group expressing each gene, and dot color represents its expression level. (E) Dot plot showing the key immuno-related genes from HC, LTBI and TB with each cell type. (F) Expression of ADM in each sample type.

      Construction of miRNA-lncRNA-mRNA-TF interaction network

      Recent studies have revealed the vital function of miRNA-lncRNA-mRNA-TF regulatory network in discovery of crucial markers and regulatory pathways
      • Matboli M.
      • Gadallah S.H.
      • Rashed W.M.
      • Hasanin A.H.
      • Essawy N.
      • Ghanem H.M.
      • et al.
      mRNA-miRNA-lncRNA regulatory network in nonalcoholic fatty liver disease.
      ,
      • Song Q.
      • Bian Q.
      • Liang T.
      • Zhang Y.
      • Zhang K.
      Identification of immune-related genes and susceptible population of pulmonary tuberculosis by constructing TF-miRNA-mRNA regulatory network.
      . As such, a miRNA-lncRNA-mRNA-TF competitive regulation network was constructed. According to the procedure mentioned in methods, miRNA-mRNA interaction pairs were firstly predicted. A total of 14 human miRNAs targeting the 3 key immuno-related genes (ADM, IFIT3 and SERPING1) were obtained. miRNA-lncRNA interaction pairs were then predicted that 62 lncRNAs and 163 miRNA-lncRNA interaction pairs were identified. Subsequently, these 3 key immuno-related genes were used to predict transcription factors (TFs). 8 TFs including CEBPB, GTF3A, TFAP2A, STAT2, SPI1, IRF1, IRF9 and STAT1 were extracted. The TF-miRNA-lncRNA-mRNA competitive regulation network was finally established by the identified 3 key immuno-related genes, 14 miRNAs, 62 lncRNAs and 8 TFs (Fig. 6). According to the regulation relationship between nodes in the network and previously published reports, hsa-miR-24–3p-NEAT1-CEBPB-ADM might be one of the possible TF-miRNA-lncRNA-mRNA regulation pathways regulating the incidence of TB.
      Fig 6
      Fig. 6Construction of miRNA-lncRNA-mRNA-TF regulatory network based on 3 key immuno-related genes (ADM, IFIT3 and SERPING1). Blue arrow represents miRNA, red circle represents lncRNA, green diamond represents mRNA, and purple rectangle represents TF. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

      Diagnostic potential of ADM

      Compared to HC and LTBI, the ADM expression level and percentage significantly increased in active TB. Thus, the diagnostic potential of ADM was evaluated by employing ROC analysis in datasets GSE28623, GSE41055 and GSE19439, respectively. The areas under the ROC curves for ADM were 0.786, 0.796 and 0.899 in the group of TB vs. HC and 0.858, 0.818 and 0.862 in the group of TB vs. LTBI, respectively (Fig. 7). ADM could distinguish TB from HC or LTBI not only in the test dataset but also in the validation dataset. It might be served as a promising marker in diagnosing TB.
      Fig 7
      Fig. 7The diagnostic potential of ADM in the GSE28623 dataset (A-B), GSE41055 dataset (C-D) and GSE19439 dataset (E-F).

      Discussion

      The severe threat of active TB to human health prompted us to seek early diagnosis and treatment strategies to effectively control TB globally. Due to the complex immune mechanisms that underlie resistance to Mtb infection and cellular heterogeneity, an in-depth understanding of transcript abundance and their roles in TB pathogenesis is urgently needed. Here, we combined expression profiling by array and scRNA-seq to comprehensively illustrate the immuno-related transcriptional signature of TB and explore potential diagnostic markers for differentiating TB from LTBI and HC. In total, 71 and 163 DEGs in TB compared with HC and LTBI were identified, respectively. Of these, 64 genes associated with TB overlapped between the two parts. The DGEs were mainly enriched in pathways involved in innate and adaptive immune response, consistent with previous reports that host immune defense plays a crucial role in controlling Mtb infection
      • Simmons J.D.
      • Stein C.M.
      • Seshadri C.
      • Campo M.
      • Alter G.
      • Fortune S.
      • et al.
      Immunological mechanisms of human resistance to persistent Mycobacterium tuberculosis infection.
      . Of interest, a similar pathway enrichment pattern was observed in the DEGs between TB and HC or LTBI, as well as the high overlapping rate of DEGs, suggesting that immune transcriptional profiling could effectively differ TB with LTBI and HC. Consequently, TB was visibly separated from LTBI and HC regarding all the immuno-related genes in the consensus clustering analysis.
      Due to the crucial role of host immune response in TB management and control, we tried to find out the key immuno-related modules highly related to TB. We identified that three immune modules were greatly correlated to TB, and PPI network analysis was performed to explore the hub genes in the pathogenesis and progression of TB. Eventually, we determined 30 hub genes related to TB. Of these, three key immuno-related hub genes highly associated with TB, including ADM, IFIT3 and SERPING1, were differentially expressed among TB, LTBI and HC and were finally extracted. We further verified that the ADM expression significantly increased in TB patients in both adult and children's datasets, which suggested that this change could distinguish patients with active TB from LTBI and HC.
      Taking advantage of scRNA-seq in providing an unbiased high-throughput approach to assess changes in transcript abundance of individual single cells or cell subsets, we found that ADM and SERPING1 were dominantly expressed in Myeloid cells at single cell resolution, while IFIT3 was also expressed in T and B cells. Specifically, we observed a significantly elevated expression level and proportion of ADM in TB compared to HC or LTBI in Myeloid cells, in line with the result of expression profiling by array. Nevertheless, the IFIT3 expression or percentage, contrary to the findings in array, significantly decreased in TB. This might be due to the vast expression of IFIT3 in multinuclear leukocytes which are absent in PBMC
      • Javed S.
      • Marsay L.
      • Wareham A.
      • Lewandowski K.S.
      • Williams A.
      • Dennis M.J.
      • et al.
      Temporal expression of peripheral blood leukocyte biomarkers in a macaca fascicularis infection model of tuberculosis; comparison with human datasets and analysis with parametric/non-parametric tools for improved diagnostic biomarker identification.
      . Taken together, these findings demonstrate the real power of scRNA-seq in determining the whole transcriptomic abundance of individual single cells or certain cell subsets through one-step operation, superior to traditional bulk RNA-seq or flow cytometry combined with real-time qPCR. We locate the origin of ADM expression that ADM was dominantly expressed by Myeloid cells at single cell level, suggesting the potential role of ADM in TB Myeloid cells. More importantly, ADM was found highly expressed in TB in both adult and children's datasets, not only from array that averages out underlying differences, but also in single cell resolution. We propose that the expression of ADM in peripheral blood could be used as a novel diagnostic biomarker for differing TB from LTBI and HC.
      Adrenomedullin (ADM), a 52 amino acid multifunctional peptide belonging to the calcitonin gene-related peptide family, is widely expressed in different tissues, especially in highly vascularized organs such as lungs and cardiac atria
      • Voors A.A.
      • Kremer D.
      • Geven C.
      • Ter Maaten J.M.
      • Struck J.
      • Bergmann A.
      • et al.
      Adrenomedullin in heart failure: pathophysiology and therapeutic application.
      . Apart from its vasodilatory and hypotensive effect, ADM has been demonstrated to exert anti-inflammatory effects through cAMP-mediated PKA activation, which prevents NF-κB from entering the nucleus or activates cAMP response element-binding protein (CREB), resulting in reduced transcription of pro-inflammatory genes or augmented transcription of anti-inflammatory cytokines
      • Geven C.
      • Kox M.
      • Pickkers P.
      Adrenomedullin and adrenomedullin-targeted therapy as treatment strategies relevant for sepsis.
      . Meanwhile, ADM was reported to induce the differentiation of IFN-γ-producing T cells via the phosphorylation of AKT and STAT3 or stimulating macrophages to produce IL-12, thus exerting a pro-inflammatory effect, contributing to the development of gastritis during H. pylori infection
      • Kong H.
      • You N.
      • Chen H.
      • Teng Y.S.
      • Liu Y.G.
      • Lv Y.P.
      • et al.
      Helicobacter pylori-induced adrenomedullin modulates IFN-gamma-producing T-cell responses and contributes to gastritis.
      . Although ADM has been identified to play multiple roles in inflammatory response and antimicrobial effect in various models, our data represent the first association of ADM with active TB as shown above. ADM might play a pivotal role in the occurrence and development of TB. Still, further studies are warranted to gain insights into the function and mechanism of ADM in TB pathogenesis.
      By exploiting the vital role of miRNA-lncRNA-mRNA-TF regulatory network in recognizing crucial regulatory pathways, a miRNA-lncRNA-mRNA-TF competitive regulation network was constructed through combining 3 key immuno-related hub genes, 14 miRNAs, 62 lncRNAs and 8 TFs. Specifically, the hsa-miR-24–3p-NEAT1-ADM-CEBPB interaction axis might be one of the critical networks regulating the pathogenesis of TB. LncRNA NEAT1 was reported highly expressed during Mtb infection, and it might be associated with the antituberculosis effect
      • Huang S.
      • Huang Z.
      • Luo Q.
      • Qing C.
      The expression of lncRNA NEAT1 in human tuberculosis and its antituberculosis effect.
      ,
      • Zheng J.
      • Wang X.
      • Shi J.
      • Tian J.
      • Chang X.
      • Wang X.
      • et al.
      Expression and clinical significance of lncRNA NEAT1 in patients with spinal tuberculosis.
      . Meanwhile, miR-24–3p was identified as a target of NEAT1
      • Luan L.
      • Hu Q.
      • Wang Y.
      • Lu L.
      • Ling J.
      Knockdown of lncRNA NEAT1 expression inhibits cell migration, invasion and EMT by regulating the miR-24-3p/LRG1 axis in retinoblastoma cells.
      . CEBPB, a transcription factor, is known to be implicated in differentiation and activation of macrophages that regulates inflammatory responses
      • Cain D.W.
      • O'Koren E.G.
      • Kan M.J.
      • Womble M.
      • Sempowski G.D.
      • Hopper K.
      • et al.
      Identification of a tissue-specific, C/EBPbeta-dependent pathway of differentiation for murine peritoneal macrophages.
      ,
      • Ruffell D.
      • Mourkioti F.
      • Gambardella A.
      • Kirstetter P.
      • Lopez R.G.
      • Rosenthal N.
      • et al.
      A CREB-C/EBPbeta cascade induces M2 macrophage-specific gene expression and promotes muscle injury repair.
      . It was predicted to interact with ADM and miR-26a-5p in a TF-miRNA-gene co-regulatory network in TB
      • Lin Y.
      • Duan Z.
      • Xu F.
      • Zhang J.
      • Shulgina M.V.
      • Li F.
      Construction and analysis of the transcription factor-microRNA co-regulatory network response to Mycobacterium tuberculosis: a view from the blood.
      . A better dissecting of the hsa-miR-24–3p-NEAT1-CEBPB-ADM co-regulatory axis will help us understand the molecular mechanism of immune response to TB further.
      Although numerous biomarker studies regarding TB diagnosis have been carried out, current active TB detection test still relies on sputum-based bacteriological confirmation, benefiting little from the existing laboratory work. However, the bacteriological approach does not make any sense while quite a few active TB patients produce no phlegm. Here, we found that the change in ADM expression is sensitive and specific to discriminate active TB from LTBI and HC in both adult and children's datasets from the perspective of bulk and single cell resolution. In considering the most readily accessibility of blood in human samples, we propose that the peripheral blood ADM expression could be served as a potential diagnostic biomarker for differing TB from LTBI and HC. There is still a high demand for further work to optimize and expand ADM application for active TB diagnosis in a larger cohort.
      In summary, our study comprehensively delineated the immuno-related transcriptional signature and explored the key molecular differences among TB, LTBI and HC through a potent brew of array-based expression profiling and scRNA-seq. To the best of our knowledge, our data represent the first scRNA-seq-based description of ADM expression abundance in TB Myeloid cells, as well as the diagnostic potential for differentiating TB from LTBI and HC. Further deciphering the mechanism behind dysregulation of key immuno-related hub genes and signaling networks, such as ADM, might provide crucial clues for establishing new therapeutic strategies for TB treatment.

      Funding

      This work was supported by Natural Science Foundation of China ( 82002109 and 82000017 ) and Natural Science Foundation of Guangdong Province ( 2022A1515012438 and 2018A0303130227 ) and Guangzhou ( 201904010071 ).

      Declaration of Competing Interests

      The authors declare that they have no competing interests.

      Acknowledgments

      The authors gratefully acknowledge all participants and staff for their contribution to the study.

      References

      1. World Health O. Global tuberculosis report 2021. Geneva: World Health Organization; 2021 2021.

        • Luo Y.
        • Suliman S.
        • Asgari S.
        • Amariuta T.
        • Baglaenko Y.
        • Martinez-Bonet M.
        • et al.
        Early progression to active tuberculosis is a highly heritable trait driven by 3q23 in Peruvians.
        Nat Commun. 2019 Aug 21; 10 (PubMed PMID: 31434886. Pubmed Central PMCID: PMC6704092. Epub 2019/08/23): 3765
        • Dheda K.
        • Barry 3rd, C.E.
        • Maartens G.
        Tuberculosis.
        Lancet. 2016; 387 (Mar 19PubMed PMID: 26377143. Epub 2015/09/18): 1211-1226
        • Simmons J.D.
        • Stein C.M.
        • Seshadri C.
        • Campo M.
        • Alter G.
        • Fortune S.
        • et al.
        Immunological mechanisms of human resistance to persistent Mycobacterium tuberculosis infection.
        Nat Rev Immunol. 2018; 18 (SepPubMed PMID: 29895826. Pubmed Central PMCID: PMC6278832. Epub 2018/06/14): 575-589
        • Kaforou M.
        • Wright V.J.
        • Oni T.
        • French N.
        • Anderson S.T.
        • Bangani N.
        • et al.
        Detection of tuberculosis in HIV-infected and -uninfected African adults using whole blood RNA expression signatures: a case-control study.
        PLoS Med. 2013; 10 (OctPubMed PMID: 24167453. Pubmed Central PMCID: PMC3805485. Epub 2013/10/30)e1001538
        • Sweeney T.E.
        • Braviak L.
        • Tato C.M.
        • Khatri P.
        Genome-wide expression for diagnosis of pulmonary tuberculosis: a multicohort analysis.
        Lancet Respir Med. 2016; 4 (MarPubMed PMID: 26907218. Pubmed Central PMCID: PMC4838193. Epub 2016/02/26): 213-224
        • Berry M.P.
        • Graham C.M.
        • McNab F.W.
        • Xu Z.
        • Bloch S.A.
        • Oni T.
        • et al.
        An interferon-inducible neutrophil-driven blood transcriptional signature in human tuberculosis.
        Nature. 2010; 466 (Aug 19PubMed PMID: 20725040. Pubmed Central PMCID: PMC3492754. Epub 2010/08/21): 973-977
        • Papalexi E.
        • Satija R.
        Single-cell RNA sequencing to explore immune cell heterogeneity.
        Nat Rev Immunol. 2018; 18 (JanPubMed PMID: 28787399. Epub 2017/08/09): 35-45
        • Deng Y.
        • Zheng Y.
        • Li D.
        • Hong Q.
        • Zhang M.
        • Li Q.
        • et al.
        Expression characteristics of interferon-stimulated genes and possible regulatory mechanisms in lupus patients using transcriptomics analyses.
        EBioMedicine. 2021; 70 (AugPubMed PMID: 34284174. Pubmed Central PMCID: PMC8318865. Epub 2021/07/21)103477
        • Cai Y.
        • Dai Y.
        • Wang Y.
        • Yang Q.
        • Guo J.
        • Wei C.
        • et al.
        Single-cell transcriptomics of blood reveals a natural killer cell subset depletion in tuberculosis.
        EBioMedicine. 2020; 53 (MarPubMed PMID: 32114394. Pubmed Central PMCID: PMC7047188. Epub 2020/03/03)102686
        • Maertzdorf J.
        • Ota M.
        • Repsilber D.
        • Mollenkopf H.J.
        • Weiner J.
        • Hill P.C.
        • et al.
        Functional correlations of pathogenesis-driven gene expression signatures in tuberculosis.
        PLoS ONE. 2011; 6 (PubMed PMID: 22046420. Pubmed Central PMCID: PMC3203931. Epub 2011/11/03): e26938
        • Verhagen L.M.
        • Zomer A.
        • Maes M.
        • Villalba J.A.
        • Del Nogal B.
        • Eleveld M.
        • et al.
        A predictive signature gene set for discriminating active from latent tuberculosis in Warao Amerindian children.
        BMC Genomics. 2013; 14 (Feb 1PubMed PMID: 23375113. Pubmed Central PMCID: PMC3600014. Epub 2013/02/05): 74
        • Ritchie M.E.
        • Phipson B.
        • Wu D.
        • Hu Y.
        • Law C.W.
        • Shi W.
        • et al.
        limma powers differential expression analyses for RNA-sequencing and microarray studies.
        Nucleic Acids Res. 2015; 43 (Apr 20PubMed PMID: 25605792. Pubmed Central PMCID: PMC4402510. Epub 2015/01/22): e47
        • Yu G.
        • Wang L.G.
        • Han Y.
        • He Q.Y.
        clusterProfiler: an R package for comparing biological themes among gene clusters.
        OMICS. 2012; 16 (MayPubMed PMID: 22455463. Pubmed Central PMCID: PMC3339379. Epub 2012/03/30): 284-287
        • Wilkerson M.D.
        • Hayes D.N.
        ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking.
        Bioinformatics. 2010; 26 (Jun 15PubMed PMID: 20427518. Pubmed Central PMCID: PMC2881355. Epub 2010/04/30): 1572-1573
        • Langfelder P.
        • Horvath S.
        WGCNA: an R package for weighted correlation network analysis.
        BMC Bioinformatics. 2008; 9 (Dec 29PubMed PMID: 19114008. Pubmed Central PMCID: PMC2631488. Epub 2008/12/31): 559
        • Germain P.L.
        • Sonrel A.
        • Robinson M.D.
        pipeComp, a general framework for the evaluation of computational pipelines, reveals performant single cell RNA-seq preprocessing tools.
        Genome Biol. 2020; 21 (Sep 1PubMed PMID: 32873325. Pubmed Central PMCID: PMC7465801. Epub 2020/09/03): 227
        • Stuart T.
        • Butler A.
        • Hoffman P.
        • Hafemeister C.
        • Papalexi E.
        • Mauck 3rd, W.M.
        • et al.
        Comprehensive Integration of Single-Cell Data.
        Cell. 2019; 177 (Jun 13e21. PubMed PMID: 31178118. Pubmed Central PMCID: PMC6687398. Epub 2019/06/11): 1888-1902
        • Robin X.
        • Turck N.
        • Hainard A.
        • Tiberti N.
        • Lisacek F.
        • Sanchez J.C.
        • et al.
        pROC: an open-source package for R and S+ to analyze and compare ROC curves.
        BMC Bioinformatics. 2011; 12 (Mar 17PubMed PMID: 21414208. Pubmed Central PMCID: PMC3068975. Epub 2011/03/19): 77
        • Matboli M.
        • Gadallah S.H.
        • Rashed W.M.
        • Hasanin A.H.
        • Essawy N.
        • Ghanem H.M.
        • et al.
        mRNA-miRNA-lncRNA regulatory network in nonalcoholic fatty liver disease.
        Int J Mol Sci. 2021; 22 (Jun 24PubMed PMID: 34202571. Pubmed Central PMCID: PMC8269036. Epub 2021/07/03)
        • Song Q.
        • Bian Q.
        • Liang T.
        • Zhang Y.
        • Zhang K.
        Identification of immune-related genes and susceptible population of pulmonary tuberculosis by constructing TF-miRNA-mRNA regulatory network.
        Tuberculosis (Edinb). 2021; 131 (DecPubMed PMID: 34740018. Epub 2021/11/06)102139
        • Javed S.
        • Marsay L.
        • Wareham A.
        • Lewandowski K.S.
        • Williams A.
        • Dennis M.J.
        • et al.
        Temporal expression of peripheral blood leukocyte biomarkers in a macaca fascicularis infection model of tuberculosis; comparison with human datasets and analysis with parametric/non-parametric tools for improved diagnostic biomarker identification.
        PLoS ONE. 2016; 11 (PubMed PMID: 27228113. Pubmed Central PMCID: PMC4882019. Epub 2016/05/27)e0154320
        • Voors A.A.
        • Kremer D.
        • Geven C.
        • Ter Maaten J.M.
        • Struck J.
        • Bergmann A.
        • et al.
        Adrenomedullin in heart failure: pathophysiology and therapeutic application.
        Eur J Heart Fail. 2019; 21 (FebPubMed PMID: 30592365. Pubmed Central PMCID: PMC6607488. Epub 2018/12/29): 163-171
        • Geven C.
        • Kox M.
        • Pickkers P.
        Adrenomedullin and adrenomedullin-targeted therapy as treatment strategies relevant for sepsis.
        Front Immunol. 2018; 9 (PubMed PMID: 29520277. Pubmed Central PMCID: PMC5827550. Epub 2018/03/10): 292
        • Kong H.
        • You N.
        • Chen H.
        • Teng Y.S.
        • Liu Y.G.
        • Lv Y.P.
        • et al.
        Helicobacter pylori-induced adrenomedullin modulates IFN-gamma-producing T-cell responses and contributes to gastritis.
        Cell Death Dis. 2020; 11 (Mar 17PubMed PMID: 32184393. Pubmed Central PMCID: PMC7078296. Epub 2020/03/19): 189
        • Huang S.
        • Huang Z.
        • Luo Q.
        • Qing C.
        The expression of lncRNA NEAT1 in human tuberculosis and its antituberculosis effect.
        Biomed Res Int. 2018; 2018 (PubMed PMID: 30534569. Pubmed Central PMCID: PMC6252192. Epub 2018/12/12)9529072
        • Zheng J.
        • Wang X.
        • Shi J.
        • Tian J.
        • Chang X.
        • Wang X.
        • et al.
        Expression and clinical significance of lncRNA NEAT1 in patients with spinal tuberculosis.
        Dis Markers. 2022; 2022 (PubMed PMID: 35465262. Pubmed Central PMCID: PMC9023229 publication of this paper. Epub 2022/04/26)5748756
        • Luan L.
        • Hu Q.
        • Wang Y.
        • Lu L.
        • Ling J.
        Knockdown of lncRNA NEAT1 expression inhibits cell migration, invasion and EMT by regulating the miR-24-3p/LRG1 axis in retinoblastoma cells.
        Exp Ther Med. 2021; 21 (AprPubMed PMID: 33732340. Pubmed Central PMCID: PMC7903428. Epub 2021/03/19): 367
        • Cain D.W.
        • O'Koren E.G.
        • Kan M.J.
        • Womble M.
        • Sempowski G.D.
        • Hopper K.
        • et al.
        Identification of a tissue-specific, C/EBPbeta-dependent pathway of differentiation for murine peritoneal macrophages.
        J Immunol. 2013; 191 (Nov 1PubMed PMID: 24078688. Pubmed Central PMCID: PMC3808250. Epub 2013/10/01): 4665-4675
        • Ruffell D.
        • Mourkioti F.
        • Gambardella A.
        • Kirstetter P.
        • Lopez R.G.
        • Rosenthal N.
        • et al.
        A CREB-C/EBPbeta cascade induces M2 macrophage-specific gene expression and promotes muscle injury repair.
        Proc Natl Acad Sci U S A. 2009; 106 (Oct 13PubMed PMID: 19805133. Pubmed Central PMCID: PMC2762675. Epub 2009/10/07): 17475-17480
        • Lin Y.
        • Duan Z.
        • Xu F.
        • Zhang J.
        • Shulgina M.V.
        • Li F.
        Construction and analysis of the transcription factor-microRNA co-regulatory network response to Mycobacterium tuberculosis: a view from the blood.
        Am J Transl Res. 2017; 9 (PubMed PMID: 28469803. Pubmed Central PMCID: PMC5411946. Epub 2017/05/05): 1962-1976