Relationships of Non-coding RNA with diabetes and depression

In order to study the molecular differences between type 2 diabetes mellitus (T2DM) and T2DM with depression (DD), we aimed to screen the differential expression of lncRNA, mRNA, and circRNA in the blood of patients with T2DM and DD. Based on the self-rating depression scale (SDS), patient health questionnaire 9 (PHQ9), blood glucose and HbA1c, we divided the patients into T2DM and DD group. Peripheral blood was collected from the two groups of patients to perform lncRNA, mRNA, and circRNA expression profiling and screening DD-related specific molecules. Subsequently, bioinformatics analysis was performed to investigate the functions of differentially expressed genes (DEgenes). Finally, RT-PCR and lncRNA-mRNA regulatory network was performed to verify the expressions of lncRNAs and mRNAs related to the occurrence and development of DD. 28 lncRNAs, 107 circRNAs, and 89 mRNAs were identified in DD differential expression profiles. GO and pathway analysis found that 20 biological process (BP) related entities and 20 pathways associated with DD. The analysis shows that the genes that are differentially expressed in the DD group involved in the development of the neuropsychiatric system, immunity, and inflammation. Then, we screening for the important DElncRNA and mRNA associated with DD were verified by RT-PCR experiments and the results of RT-PCR were consistent with the sequencing results. LncRNA, circRNA, and mRNA differential expression profiles exist in DD patients compared with T2DM. The lncRNA-mRNA regulatory network analysis confirmed the crosslinking and complex regulation patterns of lncRNA and mRNA expression and verified the authenticity of the regulatory network.

Type 2 diabetes mellitus (T2DM) and depression have become chronic psychosomatic diseases that pose the most serious threat to human health following cancer, cardiovascular and cerebrovascular diseases and AIDS. One out of every four patients with T2DM is accompanied by significant depressive symptoms 1 . The effects of T2DM with depression (DD) are extremely high and have a higher risk of death than patients with T2DM or depression alone 2 . Depression increases the development of T2DM and the subsequent risk of complications such as hyperglycemia, insulin resistance, microvascular and macrovascular; in contrast, T2DM also increases the risk of depression in the patient and may lead to more severe depression 3,4 . This association reflects a common etiology consisting of complex bidirectional interactions between multiple variables, including autonomic disorders, neurohormonal disorders, weight gain, inflammation, and hippocampal structural changes 5 .
However, the specific molecular mechanisms of the correlation between DD and T2DM are remain unclear. Long non-coding RNAs (lncRNAs) are a group of RNAs between 200 nt and 100 kb in length and have limited protein-coding potential, functioning as regulatory genes at the epigenetic, transcriptional and post-transcriptional levels 6,7 . Recent studies have shown that lncRNAs can participate in a variety of important biological control processes and the development of a variety of diseases, including a variety of neuropsychiatric and metabolic diseases including diabetes and depression [8][9][10] . CircRNA, a type of closed circular RNA molecule, can play a competitive endogenous RNA (ceRNA) role by regulating transcription and blocking miRNAs' inhibition of their target genes 11 . CircRNA is involved in the development of many diseases, and the differential expression of circRNAs plays an important role in the disease development 12 . Therefore, understanding the role of lncRNAs and circRNAs in diseases and their mechanisms is an important topic in current research. Therefore, we studied the nature of DD from LncRNA, mRNA and mRNA, and explored the pathological mechanism of DD from the transcriptomic level.

Results
Clinical characteristics of the participants. This study enrolled 5 T2DM patients and 5 DD patients.
Baseline comparisons were made between this 10 subjects, including gender, age, self-rating depression scale (SDS), and patient health questionnaire 9 (PHQ9) scores (Table 1 and Supplementary Table 1). There were no significant differences between the groups in terms of age. The mean SDS and PHQ-9 were significantly different between the two groups. The average of SDS total coarse points and standard score of DD patients were 45.6 and 56.8, and the average of PHQ-9 was 9.2.
Differentially expressed mRNA (DEmRNA) between T2DM and DD. Compared with the T2DM group, we found a total of 87 up-regulated and 2 down-regulated mRNAs in the DD group (Fold Change > 1.5, P < 0.05, Supplementary data 1). The up-regulated DEmRNA-OR2W3 had the highest difference (Fold Change = 2.179), and the Fold Change of the two down-regulated mRNAs (UBC and ORMDL3) were 0.653 and 0.343, respectively (Table 2). Among them, the up-regulated expression of CD27 in the mRNA correlates with depression 13 . Therefore, CD27 can be used as an important molecular target for the study of DD. Cluster analysis of DEmRNA showed that these DEmRNAs clearly divided the samples into T2DM and DD groups; in addition, volcano and scatter plots also visually displayed DEmRNA between the two groups ( Fig. 1). The expression levels of CD27 and CCR7 selected, as determined by quantitative real-time PCR, were consistent with the sequencing results ( Fig. 2), verifying the accuracy and reliability of the sequencing data.
Differentially expressed circRNAs (DEcircRNAs) between T2DM and DD. Compared with the T2DM group, 75 down-regulated and 32 up-regulated circRNAs were detected in the DD group (Fig. 3, Supplementary data 2). As shown in Table 3, we listed the top 10 circRNAs with up-regulation and down-regulation, respectively, based on log2 (Fold Change). Compared with the T2DM group, circRNA-TFRC and TXLNG2P were the most significantly up-regulated and down-regulated expressions in the DD group, respectively. Therefore, we speculate that these DEcircRNAs may affect depression in T2DM patients through many molecular mechanisms.
Differentially expressed lncRNAs (DElncRNAs) between T2DM and DD. Compared with T2DM group, there were 25 up-regulated lncRNAs and 3 down-regulated lncRNAs detected in DD group (Fold Change > 1.5 and P < 0.05, Fig. 4, Supplementary data 3). XIST and AP000350.5 were up-regulated and down-regulated expressions of the most significant lncRNAs, respectively, which were 6.7 and 0.59 fold change higher than those in the T2DM group, respectively (Table 4). We found that XIST belongs to intergenic lncRNA and is closely related to the pathogenesis of T2DM with depression. We will use it as an important research target for the development of DD. The expression levels of three selected lncRNAs (XIST, RP11-706O15.3, and www.nature.com/scientificreports www.nature.com/scientificreports/ RP11-415F23.2), as determined by quantitative real-time PCR, were consistent with the sequencing results ( Fig. 2), verifying the accuracy and reliability of the m sequencing data.
Classification of differentially expressed lncRNAs. DElncRNAs are divided into six classes, including exon sense-overlapping, intron sense-overlapping, intergenic, natural antisense, intronic antisense and bidirectional, based on the relative positions of the DElncRNAs and their protein-coding genes on the genome (Fig. 5). In the two groups of DElncRNAs, the number of lncRNAs of the exon sense-overlapping type was the most, there were 477 up-regulated and 247 down-regulated expressions; followed by intron sense-overlapping lncRNAs, a total of 66 DElncRNAs were belonged to this type, 52 of them are up-regulated and 14 is down-regulated. There are 50 and 42 lncRNA belong to intronic and natural antisense, respectively. In summary, we found that the DElncRNAs in the DD group were mainly located in the exon sense-overlapping part, followed by intron sense-overlapping part. , the x-axis and y-axis represent the mean FPKM values (log2 transformation) of each group of sample genes. The two oblique lines divide the upper and lower genes (1.5 fold difference) and the unmodified genes. Volcano-Plot (B), the x-axis represents the log2Fold_Change value and the y-axis represents -log10p_value. The vertical two green lines are up (right) and down (left), and the green parallel lines correspond to p-value. Hierarchical clustering (C) was performed using FPKM values of significant expression genes obtained by comparison between groups, each row representing one gene and each column representing one sample. The red dots represent up-regulated DEgenes, the green dots represent downregulated DEgenes, and the gray dots represent indistinguishable genes. PCA analysis (D), each color represents each grouping, and each axis represents a principal component.
www.nature.com/scientificreports www.nature.com/scientificreports/ Co-analysis of differentially expressed antisense lncRNAs and mRNAs. Antisense lncRNAs can induce epigenetic changes in chromatin and DNA, thereby affecting the corresponding sense mRNA expression 14 . Therefore, we integrated the antisense lncRNAs and their sense mRNAs to further infer the function of differential lncRNAs. We found that in DD patients, one antisense lncRNA may have multiple sense mRNAs corresponding to it. In addition, there are two types antisense (intronic and natural antisense) of the relationship between antisense lncRNAs and sense mRNA. In this study, we found a total of 34 sense mRNAs corresponding to antisense lncRNAs, of which 20 were intronic antisense and 14 were natural antisense (Table 5).
GO analysis. GO analysis found that DEmRNAs in the blood of DD patients were mainly enriched in 20 BP entries, 22 CC entries, and 13 MF entries, respectively. In the BP category, DEmRNAs are mainly enriched in homeostasis of number of cells、cell surface receptor signaling pathway and hypothalamus development; in CC entries, DEmRNAs are mainly enriched in the external side of plasma membrane、 hemidesmosome and neuronal cell body, etc. In addition, DEmRNAs are mainly concentrated in the cysteine-type endopeptidase inhibitor activity involved in apoptotic process、RNA polymerase II activity and C-C chemokine receptor activity in the MF entry. Therefore, we speculate that the occurrence and development of DD may be related to GO items such as brain and neuropsychiatric system development and chemokines (Fig. 6).

KEGG Pathway analysis.
Pathway analysis showed that DEmRNAs were involved in 20 pathways. Among them, there were 11 statistically significant pathways (P < 0.05): Hematopoietic cell lineage, RNA polymerase, Cytokine-cytokine receptor interaction, T cell receptor signaling pathway, Thyroid hormone signaling pathway, ras signaling pathway, Cell adhesion molecules, CAMs, HTLV-I infection and Ribosome (Fig. 7). The T cell receptor and Ribosome signaling pathway are closely related to the occurrence and development of DD. A total of three DEmRNAs participate in the T cell receptor pathway and five DEmRNAs participate in the C-C chemokine receptor activity pathway. CD27 is an important target we screened and it is involved in Cytokine-cytokine receptor interaction. Therefore, we speculate that Cytokine-cytokine receptor interaction plays an important role in the development of DD. Fig. 8, differential lncRNAs with the same biological process and the same pathway are clustered together. Compared with the T2DM group, we found that the DElncRNAs in the DD group were mainly related to mitotic spindle organization, microtubule based movement, sister chromatid, positive regulation of cytokinesis and other BP and natural killer cell mediated cytotoxicity, leukocytes transendothelial migration, oxidative phosphorylation, T cell receptor signaling pathway, cysteine, methionine metabolism and other pathways. GSEA analysis makes up for the insufficiency of single lncRNAs in the analysis, providing insights for the next in-depth study of the role of differential lncRNAs in the development of DD.

Differentially lncRNAs Gene Set Enrichment Analysis (GSEA Prerank). As shown in the
CNC network analysis. We selected two significantly up-regulated lncRNAs (ENSG00000271964.1_2 and ENSG00000279995.1_2) to construct the lncRNA-mRNA interaction network using Cytoscape v2.8.3 (Fig. 9). According to the network, we found a total of 8 mRNAs related to these two differentially expressed lncRNAs. Specifically, ENSG00000271964.1_2 was associated with seven mRNAs, including the important differential mRNA-CD27, CCR7, NELL2, and EPPK1 that we screened; whereas ENSG00000279995.1_2 was positively correlated with DDX11. These DEmRNA and lncRNA will be the focus of our further research.

Discussion
By annotating the DElncRNAs, we found that up-regulated expression of lncRNA-XIST was associated with the development of DD. Genetic studies have found that inactive X chromosomes often lead to psychiatric symptoms in patients 15 . XIST (X inactive specific transcript) is a unique length of 19 kb lncRNA, which is the major lncRNA regulating X chromosome inactivation (XCI) 16,17 . XIST was significantly overexpressed in lymphoblasts www.nature.com/scientificreports www.nature.com/scientificreports/ of patients with severe depression, which may be due to compensatory responses caused by XCI deficiency 18 . The unique pattern of neuronal expression of the XIST chromosomal gene may lead to the prevalence of certain neuropsychiatric disorders. Our study found that compared with the T2DM group, the expression level of XIST in the DD group was significantly upregulated. This result is consistent with previous reports, indicating that the expression level of lncRNA-XIST is indeed correlated with the onset of depression and is expected to be useful as a biomarker for diagnosing mental disorders in patients. However, the specific relationship between them will be our future research direction.
In addition, up-regulated expression of RP11-706O15.3 and RP11-415F23.2 also attracted our attention, with an up-regulation multiple of 2.527 and 1.969, respectively. Although these two lncRNAs were not found in previous studies, due to their significant differential expression, we suspect that they may play an important role in the development of DD.
To gain insight into the regulatory functions of lncRNAs, we focused on differentially expressed antisense lncRNAs. Antisense lncRNAs are RNA molecules that are transcribed from the antisense strand and regulate   www.nature.com/scientificreports www.nature.com/scientificreports/ their corresponding sense mRNAs at the transcriptional and post-transcriptional levels [19][20][21] . In this study, we found a total of 34 target genes associated with antisense lncRNAs. Among them, interferon γ receptor 2 (IFNGR2) corresponds to antisense lncRNA-ENST00000484377. Interferon γ is a pleiotropic protein secreted by immune cells, and its signal exerts an immune action through interferon γ receptors 22 . The gene corresponding to antisense lncRNA-ENST00000613908 is calcium-regulated heat-stable protein 1 (carhsp1), a biomarker of diabetic complications, which can inhibit the hepatic gluconeogenesis gene expression by inhibiting the expression level of PPARα 23 . Therefore, carhsp1 can be used as a molecular target for diagnosis and treatment of diabetic complications. In addition, we found that up-regulated expression of antisense lncRNA-ENST00000448942 in the DD group could down-regulate TNFAIP3 expression. Studies have shown that the expression level of mRNA-TNFAIP3 is inversely related to the degree of depression 24 . Therefore, we can use up-regulated antisense lncRNA-ENST00000448942 as a molecular target for diagnosing DD.
The function and regulation mechanism of most lncRNAs remains unclear. Therefore, we constructed mRNA expression profiles to further screen DEgenes associated with DD. Inflammation is considered to be an important pathogenesis of neuropsychiatric diseases, and inflammatory immune activation is closely related to the development of depression 25 . Chemokines are part of the immune system and have an impact on brain cognition and learning behaviors, including anxiety 26 . The chemokine receptor CCR7 (CC-chemokine receptor 7) is a GTP-protein-coupled transmembrane receptor for CC chemokines and is normally expressed in the cell membranes of immune cells and endothelial cells. The study found that the activation of the pro-inflammatory   www.nature.com/scientificreports www.nature.com/scientificreports/   Table 5. Analysis of differential antisense lncRNAs and associated sense mRNAs.
www.nature.com/scientificreports www.nature.com/scientificreports/ surface marker CCR7 is associated with depressive behavior, and that activation of CCR7 may aggravate anxiety and depression behavior in the male mice 27,28 . Compared with the T2DM group, the expression level of CCR7 in the DD group was significantly upregulated (Fold change = 2.05), which may be related to the occurrence of depression.
In our study, some immune-related DEmRNAs were also found. B cells play a key role in the inflammatory response by secreting proinflammatory and anti-inflammatory factors. It has been found that the up-regulated expression of CD27, a molecule on the surface of B cells, is associated with severe depression 13 . In our study, the expression level of CD27 was significantly higher in the DD group than in the T2DM group (Fold change = 1.89). This suggests that the occurrence of DD may be related to the activation of CD27, and they may be used as targets for preventing and improving DD in future studies.
We analyzed the function of DEmRNAs through the GO and KEGG pathways and found that they are involved in many aspects of the progression of T2DM and involve a number of related pathways of depression such as CC chemokine receptor activity, cytokine-cytokine receptor interactions and T cell receptor signaling pathways. Among them, the "T cell receptor signaling pathways" and the "CC chemokine receptor activity 29 " have an important correlation with the occurrence and development of DD. Therefore, we speculate that the occurrence and development of DD may be related to GO items and pathways such as brain, neuropsychiatric system development, and chemokines.
Depression is a complex disease that affects certain specific brain circuits and components of the nervous system and is induced by the interrelationship of epigenetics and the environment 30 . CircRNA-TFRC (transferrin receptor, TFRC) is associated with insulin resistance, and the overexpression of TFRC can aggravate the risk of type 2 diabetes and metabolic diseases 31,32 . Our study found that the expression level of TFRC in blood samples of DD patients was significantly higher than that of T2DM patients, indicating that DD patients had more severe insulin resistance and symptoms of diabetes than T2DM. Our results provide evidence from the molecular level that T2DM with depression is more dangerous than T2DM alone. In addition, another study found that TFRC  www.nature.com/scientificreports www.nature.com/scientificreports/ expression levels are associated with mental disorders, suggesting that we can use TFRC as the molecular target for the study of the occurrence and development mechanism of DD 33 .
Another up-regulated circRNA-TNIK (Traf2-and Nck-interacting kinase, TNIK) is associated with depression-related mental disorders 34 . TNIK is a highly expressed serine/threonine kinase in the brain that regulates synaptic composition and activity density by stabilizing critical postsynaptic levels and has a favorable role in dendritic development and synaptic transmission 35 . In this study, compared with the T2DM group, we found that the expression level of TNIK in the DD group was significantly increased (Fold Change = 3.505, P < 0.05), suggesting that we can focus on the function of TNIK in the next study.

Materials and Methods
Sample information description. This study was conducted in accordance with the guidelines of the Declaration of Helsinki, and the study protocol was approved by the Ethical Committee of Beijing University of Chinese Medicine (BUCM) (2016BZHYLL0105). After fully explaining the purpose and nature of all the procedures used, consent has been obtained from each patient or subject and informed consent has been signed. Enrolled patients were divided into T2DM group (patient ID: 302, 109, 215, 202, and 204) and DD group (patient ID : 206, 103, 205, 117, and 104). DD patients were screened in T2DM patients using the SDS and the PHQ9.  www.nature.com/scientificreports www.nature.com/scientificreports/ Library construction and Illumina sequencing. The total RNA was extracted according to our previous experimental method 36 . Briefly, using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) and purified with an RNeasy Mini Kit (Qiagen, Hilden, Germany) according to manufacturer protocol. Enrich mRNA with oligo (dT) beads (rRNA removal kit for RNA degradation or prokaryotic); RNA sequencing library are completed by the KAPA Stranded RNA-Seq Library Prep Kit (Illumina). Finally, the constructed library was checked with the Agilent 2100 (NanoDrop ND-1000). The sequencing library was denatured to a single-stranded DNA by 0.1 M NaOH, diluted to a concentration of 8 pM, and then amplified in situ on the TruSeq SR Cluster Kit v3-cBot-HS (#GD-401-3001, Illumina), and the resulting fragment at the end was used the Illumina HiSeq 4000 sequencer was sequenced for 150 cycles.
Quantitative analysis of DERNAs. Transcript abundance calculations were performed by software StringTie and expression unit expressed by FPKM (Fragments Per Kilobase of gene/transcript model per Million mapped fragments). Genes and lncRNA were expressed at a threshold of FPKM ≥ 0.5 and the average number of genes with FPKM exceeding 0.5 in each group was considered as expressed in the group and statistically analyzed. The expression of circRNA (Backsplice junction reads count) was quantified by the software STAR and then calculated by CIRCexplorer2 for Backsplice junction reads. CircRNA expression threshold is that CircRNAs with mean CPM exceeding 100 in each group are considered as expressed in the group and statistically analyzed.
Principal component analysis. Principal Component Analysis (PCA) is a classification of data presented by the expression data of the sample. The visual distribution of the sample between the experimental group and the control group was obtained by PCA analysis. We used PCA analysis of genes with significant differences (P value ≤ 0.05) in mean values across (ANOVA) all samples and mapped PCA (using all genes for PCA plots).
Bioinformatics analysis. We used topGO software and the newest KEGG database (Kyoto Encyclopedia of Genes and Genomes) for GO and Pathway analysis to obtain significant enriched GO terms and important biological functions involved in DEgenes, ultimately identifying and differentiating the most relevant pathway for DEgenes.
GSEA is a computational method used to determine whether a given gene set has significant differences among different groups 37 . Briefly speaking, we follow the Subramanian method 38 , using Java software, first calculate the Enrichment Score (ES), then estimate the importance of the ES, and finally evaluate their importance by adjusting multiple hypothesis tests.
Coding-noncoding network (CNC) analysis. An lncRNA-mRNA regulatory network was constructed based on the inter-regulatory association between DEgenes in the blood between DD and T2DM patients using Cytoscape v2.8.2 software (http://www.cytoscape.org/). Statistical analysis. SPSS software (Version 20.0) was used in this study for the statistical analyses. The results are expressed as mean ± SEM (standard error). Comparisons between multiple groups was analyzed using one-way ANOVA analysis. Two-sided P values < 0.05 were considered statistically significant.  www.nature.com/scientificreports www.nature.com/scientificreports/

Conclusions
The onset of DD is a complex and highly regulated biological process. In this experiment, we compared the lncRNA, circRNA and mRNA between DD and T2DM and finally screened out the DEgenes profile of DD. Through the annotation of DEgenes, function, some molecules closely related to the development of DD were discovered, such as TNIK, TFRC, XIST, RP11-706O15.3, RP11-415F23.2, CD27 and CCR7. We expect that these DERNA molecules will become attractive biomarkers to enhance the research and diagnostic potential associated with DD.

Data Availability
The datasets used and analyzed during the current study are available from the first or corresponding authors on reasonable request.