Deregulated lncRNAs in B Cells from Patients with Active Tuberculosis

Role of lncRNAs in human adaptive immune response to TB infection is largely unexplored. To address this issue, here we characterized lncRNA expression profile in primary human B cell response to TB infection using microarray assay. Several lncRNAs and mRNAs were chosen for RT-qPCR validation. Bioinformatics prediction was applied to delineate function of the deregulated mRNAs. We found that 844 lncRNAs and 597 mRNAs were differentially expressed between B cell samples from individuals with or without TB. KEGG pathway analysis for the deregulated mRNAs indicated a number of pathways, such as TB, TLR signaling pathway and antigen processing and presentation. Moreover, corresponding to the dysregulation of many lncRNAs, we also found that their adjacent protein-coding genes were also deregulated. Functional annotation for the corresponding mRNAs showed that these lncRNAs were mainly associated with TLR signaling, TGF-β signaling. Interestingly, SOCS3, which is a critical negative regulator of cytokine response to TB infection and its nearby lncRNA XLOC_012582, were highly expressed in active TB B cells. Subsequent RT-qPCR results confirmed the changes. Whether upregulated XLOC_012582 causes SOCS3 overexpression and is eventually involved in the context of exacerbations of active TB represents an interesting issue that deserves to be further explored. Taken together, for the first time, we identified a set of deregulated lncRNAs in active TB B cells and their functions were predicted. Such findings provided novel insight into the pathogenesis of TB and further studies should focus on the function and pathogenic mechanisms of the lncRNAs involved in active TB.


Introduction
Tuberculosis (TB), caused by Mycobacterium tuberculosis (Mtb), remains a major challenge to human public health [1]. According to the World Health Organization, there are reported to be more than 30% of population infected with Mtb. Host immune response against Mtb is complex and multifaceted. The underlying mechanisms of TB pathogenesis are still poorly a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 understood, especially the mechanisms explaining how the immune system dysfunction in TB development [2]. Accordingly, the elucidation of molecular mechanisms in TB has been the subject of extensive research over past decades. Better understanding of the interplay between Mtb and host immune response is critical for TB control and prevention [3].
Effective cell-mediated immune response plays an essential role in the host defense against intracellular pathogen Mtb [4,5]. Immunity to TB is still understood to be driven and maintained by T-cell derived immune response [6]. However, B cells, as one of the major adaptive immune cells, are thought to play a limited role [7]. Recent evidence shows that B cells can regulate immune response to Mtb via the modulation of cytokine production and macrophage activation [8]. Moreover, one most recent study indicates that human B cells can phagocytose Mtb, which can in turn regulate the immune activation of B cells [9]. The data suggest that B cells, as effectors in both the innate and adaptive immune response, can modulate host defense against Mtb infection and play a significant role in determining the clinical outcome of TB infection [10]. Despite this, the contributory role of B cells in the protection against Mtb infection is less defined than the role of T cells.
Increasing evidence reveals that long non-coding RNAs (lncRNAs) have key functions in regulating diverse biological processes, such as transcription activation and inhibition, mRNA translation, organellar biogenesis as well as cell development [11]. LncRNAs are responsible for at least 80% of all genome transcripts and have been shown to be involved in various pathophysiologic processes and human diseases. Accumulated evidence also indicates that lncRNA plays an important role in host immune responses against invading pathogens and its misregulation has been shown in different types of infectious diseases [12]. For instance, lncRNA NRON can modulate HIV-1 replication in a NFAT-dependent manner [13]. We previously found that many lncRNAs were differentially expressed in CD4 + T cells from patients with active TB [14]. Recently, one study shows that downregulated lncRNA MEG3 eliminates mycobacteria in macrophages via autophagy [15] and another study indicates that lncRNA BC050410 inhibits CD8 + T-cell immune response in TB infection [16]. A more recent study reveals significantly altered lncRNA expression profiles in plasma from patients with TB disease [17] and another study shows that two lncRNAs, MIR3945HG V1 and MIR3945HG V2, are identified as novel candidate diagnostic markers for TB [18]. The data suggest that similar to its important role in other infectious diseases, lncRNA is demonstrated to be closely associated with Mtb infection. What remains to be seen is whether there exists an altered lncRNA profile in the active TB B cells.
Hence, the present study aimed to explore the expression profiles of lncRNA and mRNA in the B cells obtained from subjects with or without active TB, and deregulated lncRNAs and mRNAs were also evaluated in independent patient and control samples.

Human subjects
Patients with active pulmonary TB (male/female = 13:18) with a mean age of 42.1±15.5 yr (ATB group) were recruited from Weifang No. 2 People's Hospital of Shandong Province, China. Diagnosis of active pulmonary TB was based on both sputum smear and culture positive or at least one sputum culture positive, as well as typical pulmonary TB clinical symptoms. Healthy subjects (male/female = 15:20) with a mean age of 40.3±13.9 yr (control group) were enrolled from the staff of Weifang Medical University, China, were free of clinical symptoms of any infectious disease, and had no close contact of a TB patient. All involved subjects had no history of TB, diabetes, tumor or other infectious disease. Three samples randomly selected from each group (male/female = 1:2) were used in the microarray assay and all of the samples were used for further RT-qPCR confirmation.
The study was conducted in accordance with the Declaration of Helsinki. Written informed consent was obtained from each participant prior to enrollment. The study was approved by the Research Ethics Committee of Weifang Medical University, China (Consent NO.: 2013-054).

Isolation of B cells
For preparation of B cells, fasting venous blood was drawn from each subject into a coded EDTA-anticoagulant tube. Peripheral blood mononuclear cells (PBMCs) were collected from venous blood using density gradient centrifugation. B cells were then isolated from obtained PBMCs by negative selection using B Cell Isolation Kit (R&D Systems, Minneapolis, MN, USA) according to the manufacturer's instructions. Briefly, B cells were negatively selected by depletion of unwanted cells. Flow cytometric analysis clearly showed that the purity of isolated B cells was over 90% (S1 Fig). Qualified samples were immediately stored in liquid nitrogen until further use.

RNA extraction and RNA quality control
For RNA preparation, total RNA was isolated from purified B cells using TRIzol 1 Reagent (Invitrogen, Carlsbad, CA, USA) and further purified with RNeasy mini kit (Qiagen, Hilden, Germany). Quantification and quality check were performed with Nanodrop ND-1000 (Nano-Drop Technologies,Wilmington, USA) and Agilent 2100 Bioanalyzer (Agilent Technologies Europe, Waldbroon, Germany), respectively. RNA integrity was assessed by standard denaturing gel electrophoresis. Only RNA sample with good quality was used for further downstream processing.

Microarray analysis of lncRNAs and mRNAs
For gene chip hybridization, RNA sample labeling and array hybridization were performed using Quick Amp Labeling Kit, One-Color (Agilent, Santa Clara, CA, USA) according to the manufacturer's instructions. Briefly, mRNA sample was purified from total RNA after removal of rRNA with mRNA-ONLY™ Eukaryotic mRNA Isolation Kit (Epicentre, Madison, WI, USA). Each sample was then transcribed to double-stranded cDNA, synthesized into cRNA and subsequently labeled with Cyanine-3-CTP. The labeled samples were then purified using RNAeasy Mini Kit (Qiagen, Hilden, Germany). The yield and specific activity of labeled cRNAs were then measured by NanoDrop ND-1000. Only if the yield is over 1.65 μg and the specific activity is more than 9.0 pmol Cy3 per μg cRNA, can labelled cRNAs proceed to next hybridization step. After passing quality test, 1 μg of each labeled cRNAs in hybridization solution was used for hybridization on Human LncRNA Microarray v3.0 (Arraystar, Rockville, MD, USA), which contains probes of 30,586 human lncRNAs and 26,109 human protein-coding transcripts. All the lncRNAs were obtained from authoritative databases, GENCODE, UCSC Knowngene, RefSeq, UCR and many related literatures. Positive control probes for 28 housekeeping genes (NM_003753, NM_005022, NM_002046, NM_002539, NM_001861, NM_001101, NM_001614, NM_000841, NM_006098, NM_022551, NM_001536, NM_000291, NM_002107, NM_021009, EIF3D, PFN1, GAPDH, ODC1, COX4I1, ACTB, ACTG1, GRM4, GNB2L1, RPS18, PRMT1, PGK1, H3F3A and UBC) and negative control probes were used for hybridization quality control. After hybridization, each microarray slide was washed and immediately scanned using an Agilent Microarray Scanner (Agilent G2505C, Santa Clara, CA, USA).

Data analysis
Agilent Feature Extraction software (version 11.0.1.1) was used to analyze acquired array images. Quantile normalization and subsequent data processing were performed using Gene-Spring GX v12.1 software package (Agilent Technologies). After quantile normalization of the raw data, lncRNAs and mRNAs that at least 3 out of 6 samples have flags in Present or Marginal (all targets value) were chosen for further data analysis. Differentially expressed lncRNAs and mRNAs with statistical significance between active TB group and healthy controls were identified through P-value corrected using the false discovery rate method. The threshold set for differentially expressed genes was fold-change >2.0 (P<0.05). A positive fold-change value indicates upregulation and a negative fold-change value indicates downregulation. Moreover, hierarchical clustering was performed to display the distinguishable lncRNA and mRNA expression patterns among the samples.
In addition, subgroups of the deregulated lncRNAs, including large intervening noncoding RNAs (lincRNAs), lncRNAs with enhancer-like function, antisense lncRNAs as well as their paired differentially expressed mRNAs, were also identified,

GO and KEGG pathway analysis
Gene Ontology (GO) project provides a controlled vocabulary to describe gene and gene product attributes (http://www.geneontology.org) [19]. For GO enrichment analysis of the differentially expressed mRNAs, GO categories are considered as significantly enriched only if Pvalue < 0.05. Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis is used to analyze involved biological pathways of the deregulated mRNAs. The P-value denotes the significance of the pathway correlated to the conditions. The lower the P-value, the more significant considered the pathway.

RT-qPCR analysis
The expression pattern of selected lncRNAs and mRNAs in all samples was analyzed with a SYBRGreen PCR kit (TaKaRa, Dalian, China). Primers were designed and synthesized by Generay Biotech (Shanghai, China). The thermal cycling conditions for PCR reaction were as follows: initial denaturation at 95˚C for 10 min, followed by 40 cycles at 95˚C for 10 s, 60˚C for 60 s and 72˚C for 15 s. All RT-qPCR experiments included no-template controls. Each sample was detected in triplicate. Glyceraldehyde 3-phosphate dehydrogenase (GAPDH) mRNA was used as an internal control and expression level of each RNA was normalized to that of GAPDH. Fold change of RNA expression in ATB versus control group was calculated using 2 -44Ct method and a P-value of less than 0.05 was considered statistically significant.

Statistical analysis
All data were presented as mean ± standard deviation (SD) or proportions where appropriate. Student's t-test or chi-square test was used for statistical analysis where appropriate. P < 0.05 was considered statistically significant.

Profile of microarray data
Expression profiling studies were performed on RNA from three independent B samples from each group. The microarray expression data discussed in the article have been deposited into National Center for Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) and are accessible through (GEO) Series accession number GSE89552. In total, 844 lncRNAs were identified with differential expression between samples from subjects with or without active TB (fold change > 2.0, P<0.05), of which 345 lncRNAs were increased and 499 lncRNAs were decreased in samples from active TB patients, respectively. Among them, ENST00000505706 was the top increased lncRNA and ENST00000562027 was the most decreased ones.
Further data analysis showed that 361 upregulated and 236 downregulated mRNAs were also identified in ATB group compared with healthy controls, of which 118 lncRNAs and 72 mRNAs exhibited fold changes > 3.0 (Fig 1). Among them, ANKRD22 was the top upregulated mRNA and HOXC4 was the most downregulated ones. Numerous studies indicate that a variety of cellular events are disturbed in the progression of TB, ranging from matrix synthesis to cytokine expression [20]. Underlying these alterations there is the dysregulated gene expression of particular molecules. Our data suggested that molecular events in peripheral blood B cells, such as lncRNAs and mRNAs, were altered during active TB infection.

Expression signatures of the deregulated lncRNAs
As lncRNA expression is cell specific [21], to further study the lncRNA expression pattern in B cells associated with active TB, general signatures of the deregulated lncRNAs were investigated, including lncRNA classification, length and chromosome distribution. According to the positional relationship between lncRNA and its adjacent protein-coding genes in the same chromosome, lncRNAs can be roughly classified as bidirectional, intron sense-overlapping, exon sense-overlapping, intergenic, intronic antisense and natural antisense. The majority of deregulated lncRNAs in the study were related to intergenic (~63%) (Fig 2A) and distributed 401-800 nt (~38%) in length (Fig 2B). Chromosome distribution of the deregulated lncRNAs showed that chromosome 1 and 2 were the most frequent ones and numbers on these two chromosomes accounted for approximately 11% and 9% of total deregulated lncRNAs, which were higher than the expected numbers based on chromosome total lncRNA numbers, respectively ( Fig 2C). The variation of lncRNAs expression in human B cells indicated that these deregulated lncRNAs may be involved in the onset and development of active TB.

GO analysis and pathway analysis of deregulated mRNAs
GO and KEGG pathway analysis were used to analyze the potential roles that the differentially expressed mRNAs played in GO biological process and pathways. GO analysis of the deregulated mRNAs in the study identified numerous biological processes with significantly altered expression of gene products involved. Many of these processes which are upregulated in active TB B cells are mainly involved in single organism process, immune system process, immune response, response to bacterium and molecule of bacterial origin. In contrast, other processes which are downregulated in active TB B cells are related to immune activation, such as T cell activation, T cell differentiation, positive regulation of lymphocyte mediated immunity and cytotoxicity, positive regulation of leukocyte mediated immunity and natural killer cell mediated cytotoxicity (Table 1). KEGG pathway analysis showed that 23 pathways corresponded to increased mRNAs and 8 pathways corresponded to decreased mRNAs, respectively. The top 10 pathways associated with overexpressed mRNAs were primarily enriched in TB, phagosome, toll like receptor signaling pathway and TNF signaling pathway. However, enriched pathways of underexpressed mRNAs were mainly related to antigen processing and presentation, T cell receptor signaling pathway as well as natural killer cell mediated cytotoxicity ( Table 1). The data demonstrated that B cell responses to active TB infection were differentially modulated and suppressed.

Fig 1. Differentially expressed lncRNAs (A) and mRNAs (B) between ATB group and control group.
Red indicated high relative expression and green indicated low relative expression. LncRNA or mRNA with expression fold change > 3 and with FDR adjusted P value < 0.05 was considered statistically significant. ATB group: BT1, BT2, BT3; Control group: CB1, CB2, CB3.

Confirmation of the microarray results by RT-qPCR
Four deregulated lncRNAs RP11-99H8.1, ENST00000507373, RP1-90G24.6, TCONS_00024847 as well as four mRNAs CH25H, NRG1, IL-32 and HOXC4 were randomly selected to confirm the microarray data in all samples from ATB group and control group (Fig  3). The results showed that RP11-99H8.1, ENST00000507373,CH25H and NRG1 were increased, while RP1-90G24.6, TCONS_00024847, IL-32 and HOXC4 were decreased in ATB versus control samples (P<0.05). The RT-qPCR results matched well with the microarray data, which demonstrated high credibility for the microarray analysis.

Subgroup analysis of the differentially expressed lncRNAs and their adjacent mRNA pairs
In the study, the differentially expressed lncRNAs and their neighboured protein-coding genes were focused. It was noteworthy that 9 antisense lncRNAs with their adjacent mRNA pairs were identified as coregulated transcripts, of which 8 pairs showed similar expression direction (upregulated or downregulated) ( Table 2). Moreover, there were 11 pairs of aberrantly expressed enhancer-like lncRNAs and their adjacent mRNAs (distance < 300 kb) (fold- change > 2, P < 0.05), of which 6 pairs were differentially expressed in similar direction (increased or decreased) ( Table 3). In addition, profiling data based on Rinn's lincRNAs indicated that there were 21 aberrantly expressed lincRNAs with their deregulated associated protein-coding genes (fold-change > 2, P < 0.05), of which 16 pairs were deregulated in similar direction (overexpressed or underexpressed) ( Table 4). Of these, expression of XLOC_012582 lncRNA and its paired SOCS3 mRNA was further confirmed by RT-qPCR in the study (Fig 3). Detailed relationship on genomic location between these lncRNA and the adjacent proteincoding genes was shown in S1 Table. It has been shown that transcription of lncRNAs can affect expression of their nearby coding genes at the level of chromatin modification, transcription and post-transcriptional processing, and their dysregulation was involved in many diseases [22][23][24]. In the study, not all pairs of IncRNAs and mRNAs are changed in the same direction. LncRNAs have the ability to regulate transcription by affecting gene promoters through interacting with initiation complex [25]. Moreover, some lncRNAs, like antisense lncRNAs, can influence mRNAs expression through splicing, editing and translation in the post-transcriptional processing [26]. Generally, an equidirectional transcriptive target gene is for promoting expression in the promoter region, otherwise it is for suppression. In some conditions, a reversed direction is possible to promote expression in the 3 0 -UTR region [27]. Our data suggested that these deregulated lncRNAs may positively or negatively regulate their adjacent mRNAs expression and by which, they may affect B cell function and so contribute to the pathogenesis of active TB. The differential expression of lncRNA and its relationship with the protein-coding genes are of great significance in active TB. Further studies were needed to confirm these lncRNAs functions with knockdown and over-expression experiments. Italic: enriched GO Term and KEGG pathway of downregulated mRNAs. Each P-value denoted the significance of the GO Term or KEGG pathway. The lower the P-value, the more significant GO Term or the pathway was.

Discussion
The molecular determinants of B cell immune response to TB are largely unknown. Only a handful of studies to-date have examined lncRNA expression in TB disease. In the current study, we demonstrated that there was a significantly altered lncRNA and mRNA expression profile in the active TB B cells in vivo. We identified 844 lncRNAs and 597 mRNAs with differential expression between TB and non-TB B cells and we confirmed a selection of these differentially expressed transcripts by RT-qPCR. GO and KEGG pathway analysis of the deregulated mRNAs showed that biological function and signaling pathway were altered in B cells and positive regulation of B cell response against TB infection was also changed. The aberrantly expressed lncRNAs observed in the study may provide clues to the dysfunction of B cells and so to the pathophysiological properties of active TB. However, their corresponding functions remain poorly understood. It has been shown that the transcription of lncRNAs can affect the expression of their nearby coding gene [28]. In the study, we analyzed lncRNAs and their associated protein-coding genes, which could help to predict and reveal the function and mechanism of lncRNAs in active TB. Corresponding to the dysregulation of many lncRNAs, we also found that their adjacent protein-coding genes were also deregulated. The following discussion was mainly focused on these paired lncRNA and mRNAs. Functions of most associated aberrantly expressed mRNAs in active TB were also little known, hence, we discussed these lncRNAs based on the function of their associated mRNAs reported in other studies.
Our results showed that expression of 7 antisense lncRNAs was correlated with that of corresponding nearby mRNAs. Although the corresponding functions of these mRNAs in active TB remain largely unknown, based on other data, we can find that they are mainly involved in cell proliferation (CHRM3), protein ubiquitination (RNF175), genetic susceptibility (HLA-DQB1), antibody function (FCGR1A, FCGR1B), T-cell homeostasis (FAM198B) and promoter activity (TPM1) [29,30]. Moreover, nine enhancer-like lncRNAs are associated with 11 significantly differentially expressed mRNAs in the study. These mRNAs are related to TLR signaling (PTGER3), TGF-β signaling (BAMBI), cell division, growth or death (EBF3, INPP4A, RFPL3, SYN3), oxidative stress (AHSP), cell cycle (MGAT4A), lipid metabolism and inflammation (ALOX15B) [31,32]. However, the mechanisms by which they are connected to the designated lncRNAs and so to the pathogenesis of active TB remain unknown. The lncRNAs and related gene pathways detected in our study suggest the complicated molecular mechanism of active TB.
In addition, 17 differentially expressed lincRNAs and associated coding gene pairs were identified in the study. These associated protein-coding genes are linked to cell apoptosis (TACSTD2), notch and hedgehog signaling (PTCH2), cell survival or division (MUC4, IFT122), inflammasome response (NLRP3), tumor suppressor (CLDN22), TGF-β signaling (AZGP1), phosphotransferase activity (FAM69B), WNT/β-catenin signaling (ZNF488), autophagy (ZBTB16, LRRK2) and cell differentiation (KCNK10). SOCS3, a critical negative regulator of STAT3-dependent cytokine response, is one major controller of the outcome of TB infection [33]. Many reports have shown that TB progression is associated with increased SOCS3 in T cells [34], macrophages [35] and bronchoalveolar lavage [36]. We demonstrated here for the first time that SOCS3 was also upregulated in active TB B cells. Interestingly, we found that lincRNA XLOC_012582, located on the upstream of SOCS3, was also significantly increased in the same samples. Expression of SOCS3 mRNA and XLOC_012582 lncRNA was further confirmed by RT-qPCR in the study. This may imply a partial increase of regulation of SOCS3 expression, and although very little is known about the functions of this lncRNA, this transcript may hold relevance in the context of exacerbations of active TB, which represents an interesting issue that deserved to be further explored.
In summary, our result for the first time showed that many lncRNAs were differentially expressed in active TB B cells, and their functions could be predicted based on their positional relation with protein-coding genes. Although available datasets in the study were limited and these lncRNA signatures need further identification and validation, worth noting is the fact that the aberrant lncRNAs may play important role in misregulation of B cell response against TB infection, which in turn may affect TB development. These findings shed a novel light on the pathogenesis of TB and provide a basis for the diagnosis and therapy of TB.

S1 Fig. Flow cytometric analysis of the purity of isolated B cells. PBMCs before (A) and after (B) isolation of B cells.
(TIF) S1 Table. Positional relationship between lncRNA and the adjacent protein-coding genes. (XLS)