LncRNA and mRNA expression profile of peripheral blood mononuclear cells in primary Sjögren’s syndrome patients

The aim of this study was to elucidate the expression profile and the potential role of long non-coding RNA (LncRNA) in the peripheral blood mononuclear cells of primary Sjögren’s syndrome (pSS) patients. RNA-seq technology was used to detect the differentially expressed LncRNAs and mRNAs between five age-and sex-matched paired pSS patients and healthy control PBMCs. The selected LncRNAs were detected in the validation study by RT-qPCR in 16 paired pSS patients and healthy controls. The GO, KEGG, co-localization, and co-expression analysis were performed to enrich the potential gene functions and pathways. In this study, 44 out of 1772 LncRNAs and 1034 out of 15,424 mRNAs were expressed differentially in the PBMCs of pSS patients. LINC00426, TPTEP1-202, CYTOR, NRIR, and BISPR were validated as aberrantly expressed, and these LncRNAs strongly correlated with disease activity of pSS. GO and KEGG pathway analysis revealed the significant enrichment of biological processes, cellular components, and molecular function of the up and down-regulated mRNAs, which were mainly concentrated in the immune response and immune system processes. Co-localization and co-expression analysis also revealed that differentially expressed LncRNAs in the PBMCs of pSS were strongly correlated to the mRNA functioning associated with immune response and cell metastasis. Numerous LncRNAs and mRNAs were found differentially expressed in the PBMCs of pSS patients, especially NRIR and BISPR; they interacted with the co-localized and co-expressed mRNAs, which might participate in the pathogenesis of pSS through the NF-κB, JAK-STAT, and other signaling pathways that regulate cell metastasis.

www.nature.com/scientificreports/ stages of CD4 + T cell differentiation, most of which could bind the transcription factors related to T cell functions (T-bet, GATA-3, and STAT4). LincR-Ccr2-5' AS 12 and Lnc-MAF-4 13 could regulate the specific Th2 differentiation. LncRNAs were recently reported to have important regulatory roles in autoimmune diseases (AID). NETA1 (activated MAPK signaling pathway through Toll like receptor (TLR) 4) was found to be elevated in systemic lupus erythematosus (SLE) patients 14 . Some inflammatory LncRNAs involved in the pathogenesis of autoimmune diseases were also found elevated in RA 15 , MS 16 , and other diverse autoimmune diseases. LncRNAs has significant role in the pathogenesis of autoimmune diseases. Some previous literatures have reported that LncRNAs had important roles in the pathogenesis (Immune regulation and cell metabolism) of labial gland of pSS patients 17,18 . While the different expression and function have not been found in pSS patients' peripheral blood mononuclear cells (PBMCs) yet. Thus, in the present study, we used the RNA-seq technology to analyze the expression profile of LncRNAs and mRNAs in pSS patients' PBMCs to investigate and further improve the potential roles of LncRNAs in the pathogenesis of pSS.

Methods
Patients and healthy controls. Twenty-one pSS patients and 21 age-sex matched healthy controls were enrolled in this study. All the pSS patients were recruited from the Peking Union Medical College Hospital and met the American-European consensus criteria for diagnosis 19 and have not received any medication treatment. This study was approved by the Ethics Committee of Peking Union Medical College Hospital (Approval number JS-2038). This study was performed in accordance with the Declaration of Helsinki, all methods were performed in accordance with the relevant guidelines and regulations. Patients with other autoimmune diseases, active/severe infection or malignant diseases were excluded from the study. All the pSS patients and matched healthy controls in this study have signed the informed consent. The RNA-seq analysis was carried out on 5 pSS patient samples and 5 matched healthy controls while the following validation study consisted of 16 pSS patients and 16 matched healthy controls. After the patient enrolling stage, we enrolled 21 female pSS patients who met our inclusion criteria. Patient's disease activity was determined by EULAR Sjögren's syndrome disease activity index (ESSDAI) 20 , and serological examinations, including Immunoglobulin G (IgG), immunoglobulin A (IgA), immunoglobulin M (IgM), anti-SSA and SSB antibody, erythrocyte sedimentation rate (ESR), rheumatoid factor (RF), and so on, were performed for the study subjects. The detailed demographic and clinical features of patients and healthy controls are listed in Table 1.
Sample processing. Peripheral blood (10 mL) was obtained from each subject and all the samples were collected in the ethylene diamine tetraacetic acid (EDTA) tubes. Peripheral blood mononuclear cells (PBMCs) were isolated using the Ficoll density gradient centrifugation and were counted using Cellometer Auto T4 (Nexcelom Bioscience, USA). Total RNA was extracted using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) and the concentration of total RNA was measured by NanoDrop2000c spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA). Then, the total RNA was reverse transcribed into complementary DNA (cDNA) using PrimeScript RT Master Mix (TaKaRa, Dalian, China) and was stored at -80 ℃ for the following study. Quality control. Raw data (raw reads) of fastq format were firstly processed through in-house perl scripts.
In this step, clean data (clean reads) were obtained by removing reads containing adapter, reads on containing poly-N and low-quality reads from raw data. At the same time, Q20, Q30 and GC content of the clean data were calculated. And all the down stream analyses were based on the high-quality clean data.
Mapping to reference genome. Reference

Real-time quantitative polymerase chain reaction (RT-qPCR).
The stored cDNAs were processed by TB Green Premix Taq II kit (TaKaRa, Dalian, China) and 7900 HT Fast Real-Time PCR System (ABI, Foster City, CA, USA). SYBR Green Buffer, other reagents in this kit, and the cDNAs were in a total volume of 10 uL and the RT-qPCR was performed in 384-well microplates. The conditions we used to select LncRNAs were: (a) The fold change of the LncRNAs must be greater than 1.5 fold compared with healthy controls and also P-value must be less than 0.05; (b) Genes with LncRNA-mRNA repeated seruences and LncRNAs which could not find information in databases were not selected in the RT-qPCR validation study. And based on the results of the RNA-seq, we selected 11 differentially expressed LncRNAs meeting the above conditions for validation study. The selected LncRNAs' primer sequences were listed in supplementary Table 1. The expression of each LncRNA was represented as fold changes using the ΔΔ Ct method to obtain RT-qPCR results.
GO and KEGG enrichment analysis. Gene ontology (GO) enrichment analysis of differentially expressed genes or LncRNA target genes were implemented by the GOseq R package, with the gene length bias being corrected. GO terms with corrected P-value less than 0.05 were considered significantly enriched by differentially expressed genes. KEGG was a database resource for understanding high-level functions and utilities of the biological system, and KOBAS software was used to test the statistical enrichment of differential expression genes or LncRNA target genes in KEGG pathways 21,22 .

Co-expression and co-location analysis.
Cis role is LncRNA acting on neighboring target genes. We searched coding genes 10k/100k upstream and downstream of LncRNA, found the adjacent mRNAs and then analyzed their function next. Trans role is LncRNA to identify each other by the expression level. We calculated the expressed correlation between LncRNAs with R function to find out the co-expressed mRNAs and then predicted the potential functions of the LncRNAs in pSS.

Statistical analysis.
The results of the RT-qRCR were verified through t-test if satisfied normality distribution or through Nonparametric test if did not satisfy normality distribution. And the correlation between clinical features and validated LncRNAs were tested through Spearman's correlation tests. All the data were analyzed using R Statistical Package (version 3.6.3; https ://www.r-proje ct.org). A P-value (two-tailed) < 0.05 was thought to be statistically significant. www.nature.com/scientificreports/ Co-localization and co-expression between differentially expressed LncRNAs and mRNAs. Tables 3 and 4 show the results of co-localization and co-expression analysis. LncRNA might regulate the adjacent protein-coding genes. Among the 22 co-localized mRNAs, HMGB1, RSAD2, CMPK2, MVB12A, and BST2 were reported to be strongly related to the development and differentiation of immune cells and regulation of the release of the inflammatory cytokines. PGLS, MRPL34, SLC27A1, and ABHD8 were reported to be significantly related to the mitochondrial function, glycolysis, and fatty acid metabolism. Other co-localized mRNAs were related to the degradation and modification of DNAs and RNAs. Among the 21 coexpressed mRNAs, IFI44, RAPGEF5, HDHD2, CRIP2, MIB2, and NF1 were strongly related to the immunerelated signaling pathways and the release of the inflammatory cytokines. LOXL2, CALCOCO1, CRIP2, and HSD17B11 were associated with glycolysis, fatty acid metabolism, and cell metastasis. Other co-localized mRNAs were related to DNA repair, the regulation of tumor and cell morphogenesis. The gene network of the co-localized and co-expressed mRNAs are shown in Fig. 4. The differentially expressed LncRNAs might have been involved in the pathogenesis of pSS through interaction with co-expressed and co-localized mRNAs.

Differentially expressed LncRNAs and mRNAs in PBMCS of pSS patients. Five pSS patients and
The results of GO and KEGG analysis for the differentially expressed mRNAs. GO and KEGG pathway analysis was mainly used to discover the potential pathogenic mechanism in pSS patients. Figure 5 shows the results of the Go and KEGG analysis. We found that the top significantly enriched biological process (BP) for the down-regulated mRNAs were cellular macromolecule metabolic process, cellular macromolecule biosynthetic process, intracellular transport, cytoplasmic transport, and translation while the top significantly enriched BP for the up-regulated mRNAs were immune system process, response to stress, immune response, defense response, and immune effector process. The most significantly enriched cellular component (CC) for the down-regulated mRNAs were intracellular parts, intracellular organelle, membrane-bounded organelle, and intercellular organelle. In case of the up-regulated mRNAs, it was cell part, cell, intracellular part, intercellular,  www.nature.com/scientificreports/ and organelle. The significantly enriched molecular function (MF) for down-regulated mRNAs were molecular function, binding, protein binding, organic cyclic compound binding and heterocyclic compound binding. The same MFs were enriched for the up-regulated mRNAs. The results of KEGG pathway analysis for the differentially expressed mRNAs indicated that the down-regulated mRNAs in PBMCs of pSS patients were significantly enriched in ribosome, acute myeloid leukemia, legionellosis, chronic myeloid leukemia, and colorectal cancer; whereas the up-regulated mRNAs were significantly enriched in primary immunodeficiency, antigen processing and presentation, herpes simplex infection, RIG-I-like receptor signaling pathway, and influenza A, which were all closely related to the immunity.

Discussion
With the development of the sequencing technology, LncRNAs have increasingly being found to be involved in various biological processes. Through its interaction with DNA, RNA, and proteins, LncRNAs could influence the activity of immune cells at transcriptional or post-transcriptional levels 23 . Nowadays, an increasing number of aberrantly expressed LncRNAs and their functions are being reported in the autoimmune diseases and infectious diseases 24,25 . Although the potential role of LncRNAs in pSS has been speculated previously, the detailed role of LncRNAs in the pathogenesis of pSS remains unknown. In this study, we used the RNA-seq technology to identify the differentially expressed LncRNAs and mRNAs in PBMCs of pSS patients, and to find out the potential function of these genes.
In this study, we enrolled 21 paired pSS patients and healthy controls. Through sequencing analysis we found that 44 out of 1772 LncRNAs and 1034 out of 15,424 mRNAs were differentially expressed in the PBMCs of pSS patients. www.nature.com/scientificreports/ The up-regulated mRNAs were significantly enriched in the immune system process, response to stress, and immune response, indicating that a variety of immune-related genes were aberrantly up-regulated and participated in the pathogenesis of pSS. And the KEGG pathway analysis revealed that up-regulated mRNAs were enriched in primary immunodeficiency, antigen processing and presentation, herpes simplex infection, RIG-Ilike receptor signaling pathway, and influenza A. Similar to other autoimmune diseases, pSS patients also had a Table 3. The results of co-localization analysis between differentially expressed LncRNAs and mRNAs.  www.nature.com/scientificreports/ variety of autoantigens (mainly SSA and SSB antigen), which could be presented to the T cells or be recognized by B cells to trigger the release of autoantibodies to promote disease progression 26,27 , and this attached importance  www.nature.com/scientificreports/ indicating that LncRNAs were related to the regulation of the immunity and metastasis of PBMCs; thus, contributing to the pathogenesis of pSS. Eleven selected LncRNAs were involved in the RT-qPCR validation stage. Among them, LINC00426, TPTEP1-202, CYTOR, NRIR, and BISPR were conformed to be aberrantly expressed. These differentially expressed LncRNAs were strongly related to the ESSDAI and serum IgG, CRP, and C4 (The correlation of CRP and C4 seemed to be affected by the particularly high outlier samples and were still needed the larger sample size study), which were closely correlated with disease activity of pSS. Thus, these differentially expressed LncRNAs were thought to affect pSS disease activity, which would be further verified in the subsequent experimental study. It was noted that the NRIR and BISPR were significantly up-regulated in the pSS patients as compared to the healthy controls. LncRNA-NRIR was located closely to the protein-coding interferon stimulated genes (ISGs), and CMPK2 was strongly correlated to the JAK-STAT signaling pathway, which acts as a negative regulator of IFN response 29 . It was reported that the knockdown of the LncRNA resulted in an increase in the Type I IFNstimulated transcription. While BISPR was a BST2 (IFN-α stimulation) promotor-sharing LncRNA 30 that was dependent on the JAK-STAT signaling pathway, overexpression of BISPR could up-regulate BST2 and induce the release of IFN-α. The activation of Type I IFN pathway has been found to play a vital role in the pathogenesis of pSS patients before 31,32 . Various factors like viral infection or immunity disorders can promote the release of type I IFN in pSS patients and activate the immune cells to release the inflammatory cytokines, antigen presentation, and production of autoantibodies, which was consistent with our LncRNA sequencing results. The differentially expressed LncRNA-NRIR and BISPR also attached importance to the Type I IFN signaling pathway in the pathogenesis of pSS. Wang et al. in their study found that LncRNA-TMEVPG1, which could regulate the release of IFN-γ, elevated in the pSS patients' CD4+ T cells 33 and were related to patients Th1 cells, ESR, serum IgG and anti-SSA antibody, this results just similar to the result of our study that the LncRNAs, especially the IFN related LncRNAs, were significantly differently expressed in the pSS patients PBMCs. While Some novel LncRNAs were also found differentially expressed in the sequencing study, like XLOC_167004, XLOC_102167, and XLOC_078623, their fold changes showed that they were highly expressed and their potential function could be verified in the future experiments.
LncRNAs could regulate the function of adjacent or correlated mRNAs. Therefore, in our study, we used the co-expression and co-localization analysis between the conformed differentially expressed LncRNAs and mRNAs to elaborate the potential mechanism of the LncRNAs in PBMCs of pSS patients. The results revealed that many immune-related genes were involved in PBMCs of pSS. HMGB1, which is located downstream of LINC00426, is a DNA-binding nuclear protein that has been implicated in some inflammatory disorders 34,35 . It could activate the innate immune cells (macrophages/monocytes) through interaction with TLR-2/4 and induce the release of cytokines like IL-8 through the NF-κB or MAPK signaling pathways. RASD2 (located upstream of NRIR) and IFI44 (correlated with NRIR) were reported to be the IFN inducers dependent on the JAK-STAT signaling pathway 36 , which can activate the innate immune cells and promote the release of IFN. CMPK2 (located downstream of NRIR) was thought to be related with the differentiation of monocytes and the release of I-IFN through the JAK-STAT signaling pathway 37 . BST2 (located upstream of BISPR) was reported to promote the innate immune release of inflammatory cytokines through the NF-κB signaling pathways 38 . NF1 (correlated with NRIR) could promote the development and differentiation of the immune cells through mTOR signaling pathway 39 . These immune cells, inflammatory cytokines, and signaling pathways were reported to be important in the pathogenesis of pSS. LOXL2 (correlated with LINC00426), CALCOCO1 (correlated with LINC00426), CRIP2 (correlated with TPTEP1-202), HSD17B11 (correlated with CYTOR), and SLC27A1 (upstream of BISPR) could regulate the metastasis of the PBMCs to regulate the corresponding cell function and the release of cytokines [40][41][42][43][44] , which could provide new insights into the pathogenesis of pSS. LncRNA-PVT1 participated in the regulation of cell metabolism through reprogramming glycolysis in CD4+ T cells to regulate immunity have been reported recently 45 , and similar to our results, LncRNAs regulating immunity through cell metabolism provided new direction for the pathogenesis of pSS. The results of the co-localization and co-expression analysis indicated that LncRNAs which could regulate the immune response of PBMCs were involved in the pSS disease progression. However, the mechanism of LncRNA still needs further studies, and animal models should also be used to further confirm the function of LncRNAs.