Differential Profile of Plasma Circular RNAs in Type 1 Diabetes Mellitus

Article information

Diabetes Metab J. 2020;44(6):854-865
Publication date (electronic) : 2020 July 13
doi : https://doi.org/10.4093/dmj.2019.0151
1Department of Endocrinology, The Second Hospital of Jilin University, Changchun, China
2Department of Endocrinology, Sir Run Run Hospital, Nanjing Medical University, Nanjing, China
3Department of Endocrinology, Children's Hospital of Nanjing Medical University, Nanjing, China
4Division of Arthritis and Rheumatic Diseases, Oregon Health & Science University School of Medicine, Portland, OR, USA
5Section of Rheumatology, VA Portland Health Care System, Portland, OR, USA
Corresponding author: Yu Liu. Department of Endocrinology, Sir Run Run Hospital, Nanjing Medical University, Nanjing 211166, China. drliuyu@njmu.edu.cn
Corresponding author: Cong-Qiu Chu. Division of Arthritis and Rheumatic Diseases, Oregon Health & Science University School of Medicine, Portland, OR 97239, USA. chuc@ohsu.edu
Received 2019 July 31; Accepted 2019 November 30.

Abstract

Background

No currently available biomarkers or treatment regimens fully meet therapeutic needs of type 1 diabetes mellitus (T1DM). Circular RNA (circRNA) is a recently identified class of stable noncoding RNA that have been documented as potential biomarkers for various diseases. Our objective was to identify and analyze plasma circRNAs altered in T1DM.

Methods

We used microarray to screen differentially expressed plasma circRNAs in patients with new onset T1DM (n=3) and age-/gender-matched healthy controls (n=3). Then, we selected six candidates with highest fold-change and validated them by quantitative real-time polymerase chain reaction in independent human cohort samples (n=12). Bioinformatic tools were adopted to predict putative microRNAs (miRNAs) sponged by these validated circRNAs and their downstream messenger RNAs (mRNAs). Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were performed to gain further insights into T1DM pathogenesis.

Results

We identified 68 differentially expressed circRNAs, with 61 and seven being up- and downregulated respectively. Four of the six selected candidates were successfully validated. Curations of their predicted interacting miRNAs revealed critical roles in inflammation and pathogenesis of autoimmune disorders. Functional relations were visualized by a circRNA-miRNA-mRNA network. GO and KEGG analyses identified multiple inflammation-related processes that could be potentially associated with T1DM pathogenesis, including cytokine-cytokine receptor interaction, inflammatory mediator regulation of transient receptor potential channels and leukocyte activation involved in immune response.

Conclusion

Our study report, for the first time, a profile of differentially expressed plasma circRNAs in new onset T1DM. Further in silico annotations and bioinformatics analyses supported future application of circRNAs as novel biomarkers of T1DM.

INTRODUCTION

Type 1 diabetes mellitus (T1DM) is a chronic autoimmune disorder that is characterized by destruction of insulin-producing pancreatic islet β-cells. Although the exact number of individuals with T1DM worldwide is unknown, as many as approximately 78,000 youths are diagnosed with T1DM each year [1]. T1DM can cause serious complications including ketoacidosis, kidney failure, heart disease, stroke and blindness. Aside from its physical and social burdens, T1DM also incurs an estimated economic loss of $14.4 to 14.9 billion merely in the USA [2]. As a result, considerable efforts have been devoted to the discovery of biomarkers for T1DM, such as serum autoantibodies against β-cell antigens, HLA genes and several nongenetic markers, that are expected to either predict the onset/progression of T1DM or aid in diagnosing the asymptomatic stage of this disease [3]. However, all these biomarkers have major limitations in terms of specificity at risk identification, early diagnosis prior to full-blown autoimmunity and assessment of therapeutic outcomes [4]. Therefore, it is urgently needed to identify novel biomarkers that are associated with initiation and progression of T1DM.

Circular RNAs (circRNAs) are a class of non-coding RNA transcripts that is characterized by covalent bonding between 3′ and 5′ ends, generating an unique enclosed structure [5]. Although the majority of circRNAs still lack functional significance, some are known to regulate gene expression by acting as microRNA (miRNA) sponges, nuclear transcription regulators and messenger RNA (mRNA) splicing competitors [6]. The high biological stability of circRNAs is a prerequisite for their usage as clinical biomarkers. Indeed, certain circRNAs are involved in numerous physiological and pathological events, such as brain development [7], Alzheimer disease [8], epithelial-to-mesenchymal transition [9], and atherosclerosis [10]. Markedly, many circulating circRNAs are differentially expressed in patients with type 2 diabetes mellitus (T2DM) [11], wherein hsa_circ0054633 has been proposed as a diagnostic biomarker of pre-diabetes and T2DM [12]. However, the expression profile and clinical significance of circRNAs in T1DM have remained unknown.

Given the urgent need for novel biomarkers for T1DM and the potential diagnostic values of circRNAs, we performed the current study to identify differentially expressed plasma circRNAs and to explore their roles in new onset T1DM.

METHODS

Human subjects

This study was reviewed and approved by the Institutional Ethics Review Board of the Second Hospital of Jilin University (Approval ID: 2013-009). Written and informed consents were obtained from legal guardians of participants in accordance with the Declaration of Helsinki. Specifically, six subjects were recruited for the microarray study: three patients with newly diagnosed T1DM and three age-/gender-matched healthy controls (Table 1). For the validation study, we recruited another 12 patients with new-onset T1DM by the same criteria as used in microarray and another 12 age-/gender-matched healthy controls (Table 1). All participants were Han Chinese that aged between 5 to 18 years old. Patients' state of T1DM was determined by the diagnostic criteria from American Diabetes Association [13] and by the World Health Organization reports on the classification of diabetes [14]. Specifically, they were positive for at least one islet autoantibody (glutamic acid decarboxylase autoantibodies [GADA], insulinoma associated-2 autoantibodies [IA–2A], and zinc transporter-8 autoantibodies [ZnT8A]), had little fasting serum C-peptide (measured by radioimmunoassay), showed complete β-cell loss during follow-up visits and required exogenous insulin therapy. Diminished secretion capacity of β-cells was confirmed by a 2-hour oral glucose tolerance test, whereby patients were fasted for at least 10 hours before ingestion of glucose (1.75 g/kg or a maximum of 75 g of glucose). Blood samples were obtained at 0 (fasting level), 30, 60, and 120 minutes after glucose administration for C-peptide measurements. Participants receiving any types of drug interventions, including nonsteroidal anti-inflammatory drugs, disease modifying anti-rheumatic drugs, glucocorticoids and immune suppressants, or those with severe liver or kidney diseases, were excluded.

Characteristics and demographics of human participants

RNA extraction, sample preparation, and microarray hybridization

Total RNA was isolated from human plasma using TRIzol LS reagent according to manufacturer's instructions (Invitrogen, Carlsbad, CA, USA). Integrity and concentrations of total RNA samples were determined by NanoDrop ND-1000 before sending to KangChen Biotech (Shanghai, China) to perform Arraystar Human circRNA Microarray. All experimental procedures of microarray were conducted according to a standard protocol provided by the company.

Microarray data analysis

Scanned images were imported into Agilent Feature Extraction software for raw data extraction. Quantile normalization of raw data and subsequent data processing were performed using the R software limma package (R Foundation for Statistical Computing, Vienna, Austria). Fold-change and P value between T1DM and healthy controls were calculated by ratio of group means and t-test, respectively. Differentially expressed circRNAs were defined as having fold-change ≥2 and P<0.05, and were displayed in volcano plot. Hierarchical clustering was also created to visualize patterns of circRNA expression between T1DM and healthy controls. The raw data supporting the conclusions of this manuscript will be made available by the authors, without undue reservation, to any qualified researchers.

Prediction of miRNA response elements

The circRNA:miRNA interaction was predicted with Arraystar's home-made miRNA target prediction software based on TargetScan & miRanda. Differentially expressed circRNAs were annotated in detail with the information of their interacting miRNAs.

Quantitative real-time polymerase chain reaction validation

Total plasma RNA (free of RNase treatment) were reversely transcribed by SuperScript III reverse transcriptase according to the supplier's protocol (Invitrogen). Then quantitative real-time polymerase chain reaction (qRT-PCR) was performed on ViiA 7 Real-time PCR System machine (Applied Biosystems, Foster City, CA, USA) using HieffTM SYBR green assay kit (Yeasen, Shanghai, China). Specific divergent primers were designed to amplify the circular transcripts. PCR reaction settings included an initial stage at 95℃ and 60℃ for 60 seconds, and finally melting curve analysis. All the results were first normalized to β-actin, quantified by the 2−ΔΔCt method and expressed as fold-change over healthy controls. Details of all the primers are summarized in Supplementary Table 1.

Construction of the circRNA-miRNA-mRNA network

Downstream mRNAs targeted by circRNA-interacting miRNAs were predicted by TargetScan 7.2 software, wherein mRNAs with cumulative weighted context++ scores lower than −0.4 were selected as putative candidates. To further elucidate relations among the four validated circRNAs, their interacting miRNAs and the predicted downstream mRNAs, a circRNA-miRNA-mRNA network was constructed and displayed via Cytoscape 3.6.1. For the convenience of viewing, only mRNA candidates with higher stringency (cumulative weighted context++ scores lower than −0.7) were used to generate the network.

GO and KEGG analyses

To further explore functions of those four validated circRNAs, their predicted downstream target mRNAs were submitted and annotated by the Database for Annotation, Visualization and Integrated Discovery bioinformatic tool (DAVID; http://david.abcc.ncifcrf.gov/) for Gene Ontology (GO) enrichment analysis that includes three terms: biological processes (BP), cellular component (CC), and molecular function (MF), and for the Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis that identifies altered pathways potentially mediated by dysregulated circRNAs in T1DM.

Statistical analysis

Microarray data were first analyzed by t-test, then corrected for false discovery rate by Benjamini Hochberg. The qRT-PCR validation data were expressed as mean±standard error of the mean. Two-tail Student t-test was used for data analysis by GraphPad Prism. P<0.05 was considered to be statistically significant.

RESULTS

Identification of differentially expressed plasma circRNAs by new onset T1DM

To identify differentially expressed plasma circRNAs, we performed Arraystar Human circRNA Array analysis in three patients with new onset T1DM and three age-/gender-matched healthy controls (Table 1). A total of 8,124 circRNAs were detected (data available upon request). Box plots showed that distributions of circRNA expression intensities were virtually the same across all the six samples after normalization, indicating highly comparable datasets with similar global expression patterns (Fig. 1A). Scatter plot provided an overview of circRNAs with an absolute fold-change ≥2, which are indicated by data points lying beyond the upper and bottom threshold lines (Supplementary Fig. 1). Among those circRNAs, 68 were considered differentially expressed with P<0.05 and are depicted in the Volcano plot, wherein 61 and 7 were up- and down-regulated, respectively (Fig. 1B). These differentially expressed circRNAs are also detailed in Supplementary Table 2. A heatmap was constructed to visualize these circRNAs and showed distinguishable expression patterns between the two groups (Fig. 1C). Further classification of differential circRNAs based on their genomic origins showed that the upregulated circRNAs included 85% exonic, 5% intronic, 7% sense overlapping, and 3% antisense. In contrast, all the downregulated circRNAs were exonic (100%) (Fig. 2A). We also summarized distributions of differential circRNAs across human chromosomes (Fig. 2B). Together, our results demonstrated an altered expression profile of plasma circRNAs during new onset T1DM.

Fig. 1

Overview of differentially expressed plasma circular RNAs (circRNAs) identified in patients with newly diagnosed type 1 diabetes mellitus (T1DM) by microarray. (A) The box plot shows intensity distribution of expressed circRNA across all the samples after normalization. The central line within each box represents the median of the data, whereas the error bars represent the upper and lower quartiles. (B) Volcano plots show differentially expressed plasma circRNAs with fold-change greater than 2 and P value less than 0.05 between control and T1DM subjects. The upwards and downwards arrows indicate up- and down-regulated circRNA clusters, respectively. (C) Hierarchical cluster analysis (heat map) for visualizing differentially expressed circRNAs, wherein red and green colors denote high and low expression levels, respectively. C1–C3, healthy controls; D1–D3, T1DM patients.

Fig. 2

Classification of differentially expressed circular RNAs (circRNAs). (A) Genomic origins of differentially expressed circRNAs. (B) Chromosome distribution of differentially expressed circRNAs.

Validation of differentially expressed circRNAs by qRT-PCR

Most of our differentially expressed circRNAs were upregulated. Besides, changes of all the downregulated ones were less significant in terms of fold-change and P values. Therefore, we focused on upregulated circRNAs and randomly selected six candidates with highest magnitude of fold-change for further analysis. We first validated them by qRT-PCR using independent plasma samples from another 12 patients with new onset T1DM and 12 age-/gender-matched healthy controls (Table 1). Intriguingly, four circRNAs (hsa_circRNA_101062, hsa_circRNA_100332, hsa_circRNA_085129, and hsa_circRNA_103845) were significantly upregulated again, while changes of one circRNA (hsa_circRNA_005178) did not reach statistical significance (Fig. 3). The other circRNA (hsa_circRNA_059571) was under detection (data not shown). Specificity of qPCR products was confirmed by electrophoresis that showed a single distinct band (Supplementary Fig. 2). In summary, qRT-PCR validation results demonstrated a good consistency with those of microarray, confirming a high reliability of our microarray dataset.

Fig. 3

Verification of microarray data by quantitative real-time polymerase chain reaction. Using independent human samples, our results confirmed that hsa_circular RNA (circRNA)_085129 (A), hsa_circRNA_100332 (B), hsa_circRNA_101062 (C) and hsa_circRNA_103845 (D) were all upregulated in type 1 diabetes mellitus (T1DM). However, expression of hsa_circRNA_005178 (E) remained unchanged. circRNA, circular RNA. Data are expressed as fold-change over healthy controls and represented as the mean±standard error (n=12). NS, not significant. aStatistically significant difference (P<0.05), bStatistically significant difference (P<0.01).

Prediction of miRNA response elements in differentially expressed circRNAs

CircRNAs interact via intrinsic miRNA response elements (MREs) with complementary miRNAs to control the function of target mRNAs. To find these miRNAs, we used Arraystar's homemade miRNA target prediction software based on TargetScan & miRanda to predict circRNA/miRNA interaction. A total of 340 MREs (five per circRNA) were identified within all differential circRNAs (Supplementary Table 3). Then we focused on the four experimentally validated circRNAs to further explore their biological functions. Their most likely 20 complementary miRNAs are summarized in Supplementary Table 4. Detailed annotation of circRNA/miRNA interaction included depicted MRE sequence, the target miRNA seed type and 3′ pairing sequence in “2D structure” column. The “local AU” column presented compositions of A/U (red bars) and G/C (black bars) contents of 30 nucleotides upstream and downstream of the seed sequence. Height of the bars indicated extent of accessibility. The “Position” column showed position of MREs in the linear form of circRNA (Supplementary Figs 3 and 4). Remarkably, literature-based curation revealed a number of these predicted miRNAs to be implicated in normal immune response and/or pathogenesis of other autoimmune diseases (Supplementary Table 5), many of which share co-morbidity with T1DM. In summary, our results supported a tight association of these differential circRNAs with the pathogenesis of T1DM.

Construction of circRNA-miRNA-mRNA network

To better understand functions of the four validated circRNAs, using TargetScan 7.2 program we predicted the downstream targets of their complementary miRNAs and identified a total of 1,827 unique mRNAs (Data available upon request). Then, a circRNA-miRNA-mRNA network was constructed to reveal the relations among validated circRNAs, miRNAs and representative target mRNAs (Fig. 4). Markedly, hsa-miR-660-3p had the largest number of downstream target genes and hsa-miR-612 exhibited the most interactions with other circRNA clusters, suggesting that hsa_circRNA_101062 may serve as a hub regulator in the network by acting upstream of these two miRNAs.

Fig. 4

The circular RNA (circRNA)-microRNA (miRNA)-messenger RNA (mRNA) network for the validated four circRNAs. Each circRNA interacts with its five miRNA response elements (MREs) and representative downstream mRNAs. For the convenience of visualization, only mRNAs with cumulative weighted context++ score no more than −0.7 were included. Specifically, hsa-miR-660-3p had the largest number of downstream target genes and hsa-miR-5189-5p exhibited the most interactions with other circRNA clusters. Red diamonds, circRNAs. Blue squares, miRNAs. Green ovals, mRNAs.

GO and KEGG pathway analyses

To further explore the significance of our microarray data, we performed GO and KEGG analysis using the 1,827 predicted target mRNAs. For GO analysis, results were sorted by enrichment scores and top five terms of each category (BP, CC, MF) were displayed (Fig. 5A). Specifically, the top five enriched terms of BP analysis included positive regulation of protein insertion into mitochondrial membrane involved in apoptotic signaling pathway, exocytosis, axon extension involved in axon guidance, negative regulation of proteasomal protein catabolic process and leukocyte activation involved in immune response. The CC analysis suggested actomyosin, condensed chromosome kinetochore, cell junction, integral component of mitochondrial inner membrane and integral component of synaptic vesicle membrane as the top five enriched terms. In MF analysis, protein binding, phosphatidylinositol-4-phosphate binding, phosphatidylinositol-5-phosphate binding, beta-catenin binding and phosphatidylinositol-3-phosphate binding were identified as the top five enriched terms. For KEGG analysis, results were ranked according to their enrichment scores and the top 10 pathways were listed (Fig. 5B), among which cytokine-cytokine receptor interaction, GAB-Aergic synapse, retrograde endocannabinoid signaling, Wnt signaling pathway and influenza A were the top five pathways. Importantly, we identified cytokine-cytokine receptor interaction as the most significantly altered pathway that is critical for pathogenesis of T1DM [15]. Overall, functional annotations revealed possible mechanisms underlying the inflammatory state of T1DM and supported the use of our differentially expressed circRNAs as potential diagnostic biomarkers.

Fig. 5

Gene Ontology (GO) enrichment and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis using predicted target messenger RNAs (mRNAs) of validated four circular RNA (circRNA). (A) GO enrichment analysis comprises three categories: biological process (in blue), cellular component (in green) and molecular function (in orange). Top five terms of each category are displayed. (B) KEGG pathway analysis shows top 10 terms that may be involved in the regulatory network mediated by differentially expressed circRNAs in type 1 diabetes mellitus.

DISCUSSION

T1DM still remains as a dreadful threat to human health with a continually rising prevalence worldwide, despite enormous efforts dedicated for the search of therapeutic strategies. To better understand the precise mechanisms underlying T1DM and to promote discovery of novel biomarkers, here we studied the differential expression profile of plasma circRNAs via microarray in patients with new onset T1DM and healthy controls. Specifically, we detected a total of 8,214 distinct circRNAs, out of which 68 were considered significantly altered by T1DM, with 61 and seven being up- and down-regulated, respectively. Next, four upregulated circRNAs with highest magnitude of fold-change (hsa_circRNA_101062, hsa_circRNA_100332, hsa_circRNA_085129, and hsa_circRNA_103845) were confirmed by qRT-PCR in independent human plasma samples. Functional annotations by GO and KEGG pathway analyses provided novel insights to better understand the pathogenesis of T1DM. To our knowledge, this is the first report demonstrating a differential expression profile of plasma circRNAs that are potentially associated with the onset of T1DM.

Since circRNAs function as sponges to sequester complementary miRNAs [5], we performed literature-based curation on the 20 predicted complimentary miRNAs (Supplementary Table 5) and revealed many of them to be implicated in normal immune responses and pathogenesis of several autoimmune diseases [16171819]. For example, hsa-miR-892a targeted by hsa_circRNA_100332, is downregulated in CD4+ T cells during the inflammatory transition from latent to active tuberculosis [20]. Moreover, hsa-miR-145-5p targeted by hsa_circRNA_085129, ameliorates the inflammatory status of macrophages by inhibiting secretion of inflammatory factors [21]. Finally, hsa-miR-127-3p targeted by hsa_circRNA_103845, is downregulated by lipopolysaccharide (LPS) and may participate in LPS-triggered development of T1DM [22]. Therefore, our observed upregulation of all the upstream circRNAs predicts a pro-inflammatory state that potentially underlies etiology of T1DM.

Interestingly, many of the predicted complimentary miRNAs are potentially associated with other autoimmune disorders that share strong comorbidities with T1DM. One example of such disorders is celiac disease wherein hsa-miR-212-5p, hsa-miR-26b-3p, hsa-miR-145-5p, and hsa-miR-432-5p are all dysregulated. Furthermore, hsa-miR-200c-5p, hsa-miR-145-5p, hsa-miR-769-3p, and hsa-miR-127-3p are all implicated in the development of inflammatory bowel disease. Previous studies showed that hsa-miR-892a is dysregulated in rheumatoid arthritis [23], and that hsa-miR-493-5p and hsa-miR-145-5p are proposed as biomarkers for multiple sclerosis [24]. Finally, islet specific miRNAs such as hsa-miR-493-5p and hsa-miR-127-3p are tightly associated with insulin biosynthesis and secretion. In keep, our findings predicted these two miRNAs to be inhibited, possibly reflecting a destruction state of β-cells in the context of T1DM. In summary, our findings are consistent with the theory that T1DM is frequently associated with other organ-specific autoimmune diseases [25].

Two of the complementary miRNAs, hsa-miR-432-5p and hsa-miR-145-5p, were predicted to be sponged by hsa_circRNA_103845 and hsa_circRNA_085129, respectively. Consequent inhibition of these two miRNAs seems paradoxically to ameliorate the overt inflammatory status in T1DM by reducing production of inflammatory cytokines [2627]. This discrepancy may be partially explained by the fact that different sample sources were analyzed in each study. Indeed, we used human plasma samples whereas one of the contradictory studies used THP-1 monocyte cell line [27]. Alternatively, cytokines stimulated by these two miRNAs may regulate (i.e., suppress) autoimmune responses that would otherwise result in β-cell destruction. In all, the exact role of hsa-miR-432-5p and hsa-miR-145-5p in the pathogenesis of T1DM remains undetermined and warrants further investigation.

Our results do not exclude the possibility that circRNAs exert their diabetogenic effects in an miRNA sponge-independent manner. In fact, only a few reported candidate circRNAs might function as miRNA sponges, not to mention none of them has yet been validated. Some circRNAs can regulate gene expression by affecting transcription [5]. For example, some nuclear circRNAs interact with the elongating Pol II complex and promotes Pol II-dependent transcription of ankyrin repeat domain 52 (ANKRD52) and sirtuin 7 (SIRT7) genes. Another potential function of circRNAs could be to compete with mRNA splicing, as in the case of muscleblind (MBL) gene. Unfortunately, direct evidence for the regulatory effect of circRNA on T1DM still lacks.

The biological functions of our four validated circRNAs were predicted by GO and KEGG pathway analyses through interrogating downstream mRNAs (data available upon request), which were targeted by the 20 complementary miRNAs. Identified among the top GO terms were mitochondria-dependent apoptotic pathway (BP term) and integral component of mitochondrial inner membrane (CC term), supporting the hypothesis that mitochondria-mediated apoptosis is the main mode of β-cell death in the early stage of T1DM [28]. Besides, binding of multiple phosphoinositides, such as phosphatidylinositol-4-phosphate, phosphatidylinositol-5-phosphate and phosphatidylinositol-3-phosphate, were identified as top MF terms. These findings agree with phosphoinositides as required stimuli of novel protein kinase C isoforms (δ, ε, ζ, θ) to promote activation of immune cells and subsequent development of autoimmune diseases [29303132]. Notably, we identified beta-catenin binding among the top GO MF terms and concomitantly Wnt signaling pathway as one of the top dysregulated KEGG pathways. Therefore, our findings are consistent with the reported modulatory role of Wnt/β-catenin signaling in inflammatory responses and autoimmunity [3334].

Several significantly enriched GO terms and KEGG pathways are closely related to immune functions, such as cytokine-cytokine receptor interaction, inflammatory mediator regulation of transient receptor potential (TRP) channels and leukocyte activation involved in immune response. Indeed, a variety of cytokines are highly expressed in the insulitis lesion of autoimmune diabetes-prone rodent models, as well as in the pancreas of humans with T1DM [15]. Previous studies showed that some TRP channels regulate multiple biological functions such as insulin secretion out of pancreatic beta cells, neurotransmitter release and ion homeostasis, and are implicated in many inflammatory (autoimmune) diseases [3536]. Markedly, as the key regulator of autoimmune diseases, activated leukocytes are dominant in the initiation and progression of T1DM [37]. This stage of T1DM, characterized by β-cell destruction mediated by infiltrating T cells, can persist long before overt hyperglycemia occurs. Therefore, our differential circRNAs may represent such scenario and serve as potential biomarkers to diagnose this asymptomatic stage. Taken together, GO and KEGG annotations revealed possible mechanisms underlying the autoimmune nature of T1DM and corroborated circRNA-mediated regulatory network in the pathogenesis of T1DM.

In conclusion, we report for the first time a differential expression profile of plasma circRNA between patients with new onset T1DM and healthy controls. Bioinformatics analyses predicted downstream targets of validated circRNAs and revealed possible mechanisms associated with T1DM pathogenesis. Study of circRNAs in body fluids such as peripheral blood is crucial for discovering novel biomarkers noninvasively for T1DM diagnosis and treatment. Hence, our findings carry significant clinical benefits for the care of T1DM patients. However, one limitation of the current study is the lack of systematic investigation of their function experimentally. Also, we admit the association between our candidate circRNAs and stage 1 or 2 pre-symptomatic T1DM awaits exploration [38]. Thus, future studies, both in vivo and in vitro, are urgently needed to elucidate function of differential circRNAs during the initiation and progression of T1DM, eventually aiding with disease early prediction as well as conceptualization of therapeutic interventions.

ACKNOWLEDGMENTS

This research was funded by The National Key R&D Program of China (Grant No. 2016YFC1305000, 2016YFC1305005), National Natural Science Foundation of China (Grants No. 81570704 and 81770778 to Yu Liu, Grant No. 81800735 to Yangyang Li), the Key Medical Talents Project of Jiangsu Province (Grant No. ZDRCA2016088) and the Science and Technology Plan of Jiangsu Province-Social Development (Grant No. BE2017738). We thank Dr. Zejian Liu from the Department of Biochemistry and Molecular Biology at Nanjing Medical University for his assistance with manuscript editing.

Notes

CONFLICTS OF INTEREST

No potential conflict of interest relevant to this article was reported.

AUTHOR CONTRIBUTIONS

Conception or design: Y.L., C.Q.C., Y.L.

Acquisition, analysis, or interpretation of data: Y.L., Y.Z., M.Z., J.Z., Y.Z., X.Y.

Drafting the work or revising: Y.L., Q.L., H.C., C.Q.C., Y.L.

Final approval of the manuscript: Y.L., Y.Z., M.Z., J.Z., Y.Z., X.Y., Q.L., H.C., C.Q.C., Y.L.

References

1. Chiang JL, Kirkman MS, Laffel LM, Peters AL. Type 1 Diabetes Sourcebook Authors. Type 1 diabetes through the life span: a position statement of the American Diabetes Association. Diabetes Care 2014;37:2034–2054.
2. Tao B, Pietropaolo M, Atkinson M, Schatz D, Taylor D. Estimating the cost of type 1 diabetes in the U.S.: a propensity score matching method. PLoS One 2010;5e11501.
3. Bonifacio E. Predicting type 1 diabetes using biomarkers. Diabetes Care 2015;38:989–996.
4. Jin Y, She JX. Novel biomarkers in type 1 diabetes. Rev Diabet Stud 2012;9:224–235.
5. Chen LL. The biogenesis and emerging roles of circular RNAs. Nat Rev Mol Cell Biol 2016;17:205–211.
6. Salzman J. Circular RNA expression: its potential regulation and function. Trends Genet 2016;32:309–316.
7. You X, Vlatkovic I, Babic A, Will T, Epstein I, Tushev G, Akbalik G, Wang M, Glock C, Quedenau C, Wang X, Hou J, Liu H, Sun W, Sambandan S, Chen T, Schuman EM, Chen W. Neural circular RNAs are derived from synaptic genes and regulated by development and plasticity. Nat Neurosci 2015;18:603–610.
8. Lukiw WJ. Circular RNA (circRNA) in Alzheimer's disease (AD). Front Genet 2013;4307.
9. Conn SJ, Pillman KA, Toubia J, Conn VM, Salmanidis M, Phillips CA, Roslan S, Schreiber AW, Gregory PA, Goodall GJ. The RNA binding protein quaking regulates formation of circRNAs. Cell 2015;160:1125–1134.
10. Burd CE, Jeck WR, Liu Y, Sanoff HK, Wang Z, Sharpless NE. Expression of linear and novel circular forms of an INK4/ARF-associated non-coding RNA correlates with atherosclerosis risk. PLoS Genet 2010;6e1001233.
11. Fang Y, Wang X, Li W, Han J, Jin J, Su F, Zhang J, Huang W, Xiao F, Pan Q, Zou L. Screening of circular RNAs and validation of circANKRD36 associated with inflammation in patients with type 2 diabetes mellitus. Int J Mol Med 2018;42:1865–1874.
12. Zhao Z, Li X, Jian D, Hao P, Rao L, Li M. Hsa_circ_0054633 in peripheral blood can be used as a diagnostic biomarker of pre-diabetes and type 2 diabetes mellitus. Acta Diabetol 2017;54:237–245.
13. American Diabetes Association. Diagnosis and classification of diabetes mellitus. Diabetes Care 2013;36 Suppl 1:S67–S74.
14. Alberti KG, Zimmet PZ. Definition, diagnosis and classification of diabetes mellitus and its complications. Part 1: diagnosis and classification of diabetes mellitus provisional report of a WHO consultation. Diabet Med 1998;15:539–553.
15. Rabinovitch A, Suarez-Pinzon WL. Role of cytokines in the pathogenesis of autoimmune diabetes mellitus. Rev Endocr Metab Disord 2003;4:291–299.
16. Halling ML, Kjeldsen J, Knudsen T, Nielsen J, Hansen LK. Patients with inflammatory bowel disease have increased risk of autoimmune and inflammatory diseases. World J Gastroenterol 2017;23:6137–6146.
17. Liao KP, Gunnarsson M, Kallberg H, Ding B, Plenge RM, Padyukov L, Karlson EW, Klareskog L, Askling J, Alfredsson L. Specific association of type 1 diabetes mellitus with anti-cyclic citrullinated peptide-positive rheumatoid arthritis. Arthritis Rheum 2009;60:653–660.
18. Ludvigsson JF, Ludvigsson J, Ekbom A, Montgomery SM. Celiac disease and risk of subsequent type 1 diabetes: a general population cohort study of children and adolescents. Diabetes Care 2006;29:2483–2488.
19. Marrosu MG, Cocco E, Lai M, Spinicci G, Pischedda MP, Contu P. Patients with multiple sclerosis and risk of type 1 diabetes mellitus in Sardinia, Italy: a cohort study. Lancet 2002;359:1461–1465.
20. Fu Y, Yi Z, Li J, Li R. Deregulated microRNAs in CD4+ T cells from individuals with latent tuberculosis versus active tuberculosis. J Cell Mol Med 2014;18:503–513.
21. Li R, Shen Q, Wu N, He M, Liu N, Huang J, Lu B, Yao Q, Yang Y, Hu R. MiR-145 improves macrophage-mediated inflammation through targeting Arf6. Endocrine 2018;60:73–82.
22. Balasa B, Van Gunst K, Sarvetnick N. The microbial product lipopolysaccharide confers diabetogenic potential on the T cell repertoire of BDC2.5/NOD mice: implications for the etiology of autoimmune diabetes. Clin Immunol 2000;95:93–98.
23. Luo Q, Zhang L, Li X, Fu B, Deng Z, Qing C, Su R, Xu J, Guo Y, Huang Z, Li J. Identification of circular RNAs hsa_circ_0044235 in peripheral blood as novel biomarkers for rheumatoid arthritis. Clin Exp Immunol 2018;194:118–124.
24. Dolati S, Ahmadi M, Aghebti-Maleki L, Nikmaram A, Marofi F, Rikhtegar R, Ayromlou H, Yousefi M. Nanocurcumin is a potential novel therapy for multiple sclerosis by influencing inflammatory mediators. Pharmacol Rep 2018;70:1158–1167.
25. Kawasaki E. Type 1 diabetes and autoimmunity. Clin Pediatr Endocrinol 2014;23:99–105.
26. Mansouri L, Lundwall K, Moshfegh A, Jacobson SH, Lundahl J, Spaak J. Vitamin D receptor activation reduces inflammatory cytokines and plasma MicroRNAs in moderate chronic kidney disease: a randomized trial. BMC Nephrol 2017;18161.
27. Nakaoka H, Hirono K, Yamamoto S, Takasaki I, Takahashi K, Kinoshita K, Takasaki A, Nishida N, Okabe M, Ce W, Miyao N, Saito K, Ibuki K, Ozawa S, Adachi Y, Ichida F. MicroRNA-145-5p and microRNA-320a encapsulated in endothelial microparticles contribute to the progression of vasculitis in acute Kawasaki Disease. Sci Rep 2018;81016.
28. Eizirik DL, Darville MI. Beta-cell apoptosis and defense mechanisms: lessons from type 1 diabetes. Diabetes 2001;50 Suppl 1:S64–S69.
29. Jimenez JM, Boyall D, Brenchley G, Collier PN, Davis CJ, Fraysse D, Keily SB, Henderson J, Miller A, Pierard F, Settimo L, Twin HC, Bolton CM, Curnock AP, Chiu P, Tanner AJ, Young S. Design and optimization of selective protein kinase C θ (PKCθ) inhibitors for the treatment of autoimmune diseases. J Med Chem 2013;56:1799–1810.
30. Salzer E, Santos-Valente E, Keller B, Warnatz K, Boztug K. Protein kinase C δ: a gatekeeper of immune homeostasis. J Clin Immunol 2016;36:631–640.
31. Salzer E, Santos-Valente E, Klaver S, Ban SA, Emminger W, Prengemann NK, Garncarz W, Mullauer L, Kain R, Boztug H, Heitger A, Arbeiter K, Eitelberger F, Seidel MG, Holter W, Pollak A, Pickl WF, Forster-Waldl E, Boztug K. B-cell deficiency and severe autoimmunity caused by deficiency of protein kinase C δ. Blood 2013;121:3112–3116.
32. Hanes CM, D'Amico AE, Ueyama T, Wong AC, Zhang X, Hynes WF, Barroso MM, Cady NC, Trebak M, Saito N, Lennartz MR. Golgi-associated protein kinase C-ε is delivered to phagocytic cups: role of phosphatidylinositol 4-phosphate. J Immunol 2017;199:271–277.
33. Silva-Garcia O, Valdez-Alarcon JJ, Baizabal-Aguirre VM. The Wnt/β-catenin signaling pathway controls the inflammatory response in infections caused by pathogenic bacteria. Mediators Inflamm 2014;2014310183.
34. Suryawanshi A, Tadagavadi RK, Swafford D, Manicassamy S. Modulation of inflammatory responses by Wnt/β-catenin signaling in dendritic cells: a novel immunotherapy target for autoimmunity and cancer. Front Immunol 2016;7460.
35. Nilius B, Owsianik G, Voets T, Peters JA. Transient receptor potential cation channels in disease. Physiol Rev 2007;87:165–217.
36. Tsuji F, Aono H. Role of transient receptor potential vanilloid 1 in inflammation and autoimmune diseases. Pharmaceuticals (Basel) 2012;5:837–852.
37. Bluestone JA, Herold K, Eisenbarth G. Genetics, pathogenesis and clinical interventions in type 1 diabetes. Nature 2010;464:1293–1300.
38. Insel RA, Dunne JL, Atkinson MA, Chiang JL, Dabelea D, Gottlieb PA, Greenbaum CJ, Herold KC, Krischer JP, Lernmark Å, Ratner RE, Rewers MJ, Schatz DA, Skyler JS, Sosenko JM, Ziegler AG. Staging presymptomatic type 1 diabetes: a scientific statement of JDRF, the Endocrine Society, and the American Diabetes Association. Diabetes Care 2015;38:1964–1974.

SUPPLEMENTARY MATERIALS

Supplementary materials related to this article can be found online at https://doi.org/10.4093/dmj.2019.0151.

Supplementary Table 1

Sequences of quantitative real-time polymerase chain reaction primers for validation experiments

dmj-2019-0151-s001.pdf

Supplementary Table 2

List of differentially expressed plasma circRNAs identified in T1DM by microarray

dmj-2019-0151-s002.pdf

Supplementary Table 3

List of all the differentially expressed circRNAs and their interacting MREs

dmj-2019-0151-s003.pdf

Supplementary Table 4

The predicted miRNA-binding elements within validated circRNAs

dmj-2019-0151-s004.pdf

Supplementary Table 5

Role of predicted miRNAs in normal immune response and pathogenesis of other autoimmune diseases

dmj-2019-0151-s005.pdf

Supplementary Fig. 1

Scatter-plot for assessing variation of circular RNA (circRNA) expression. Points above the top green line and those below the bottom green line represent circRNAs that are up- and down-regulated in type 1 diabetes mellitus (T1DM) subjects, respectively.

dmj-2019-0151-s006.pdf

Supplementary Fig. 2

Gel electrophoresis analysis of quantitative real-time polymerase chain reaction (qRT-PCR) products. To demonstrate qRT-PCR specificity, we performed agarose gel electrophoresis using qRT-PCR products. All the circular RNA (circRNA) PCR products showed a single band, demonstrating high specificity and reliability of qRT-PCR results.

dmj-2019-0151-s007.pdf

Supplementary Fig. 3

Detailed annotation of circular RNA (circRNA)/microRNA (miRNA) interaction for hsa_circRNA_101062 and hsa_circRNA_100332. For each circRNA, five predicted miRNA response elements are included: (A) hsa_circRNA_101062, (B) hsa_circRNA_100332.

dmj-2019-0151-s008.pdf

Supplementary Fig. 4

Detailed annotation of circular RNA (circRNA)/microRNA (miRNA) interaction for hsa_circRNA_085129 and hsa_circRNA_103845. For each circRNA, five predicted miRNA response elements are included: (A) hsa_circRNA_085129, (B) hsa_circRNA_103845.

dmj-2019-0151-s009.pdf

Supplementary references

dmj-2019-0151-s010.pdf

Article information Continued

Fig. 1

Overview of differentially expressed plasma circular RNAs (circRNAs) identified in patients with newly diagnosed type 1 diabetes mellitus (T1DM) by microarray. (A) The box plot shows intensity distribution of expressed circRNA across all the samples after normalization. The central line within each box represents the median of the data, whereas the error bars represent the upper and lower quartiles. (B) Volcano plots show differentially expressed plasma circRNAs with fold-change greater than 2 and P value less than 0.05 between control and T1DM subjects. The upwards and downwards arrows indicate up- and down-regulated circRNA clusters, respectively. (C) Hierarchical cluster analysis (heat map) for visualizing differentially expressed circRNAs, wherein red and green colors denote high and low expression levels, respectively. C1–C3, healthy controls; D1–D3, T1DM patients.

Fig. 2

Classification of differentially expressed circular RNAs (circRNAs). (A) Genomic origins of differentially expressed circRNAs. (B) Chromosome distribution of differentially expressed circRNAs.

Fig. 3

Verification of microarray data by quantitative real-time polymerase chain reaction. Using independent human samples, our results confirmed that hsa_circular RNA (circRNA)_085129 (A), hsa_circRNA_100332 (B), hsa_circRNA_101062 (C) and hsa_circRNA_103845 (D) were all upregulated in type 1 diabetes mellitus (T1DM). However, expression of hsa_circRNA_005178 (E) remained unchanged. circRNA, circular RNA. Data are expressed as fold-change over healthy controls and represented as the mean±standard error (n=12). NS, not significant. aStatistically significant difference (P<0.05), bStatistically significant difference (P<0.01).

Fig. 4

The circular RNA (circRNA)-microRNA (miRNA)-messenger RNA (mRNA) network for the validated four circRNAs. Each circRNA interacts with its five miRNA response elements (MREs) and representative downstream mRNAs. For the convenience of visualization, only mRNAs with cumulative weighted context++ score no more than −0.7 were included. Specifically, hsa-miR-660-3p had the largest number of downstream target genes and hsa-miR-5189-5p exhibited the most interactions with other circRNA clusters. Red diamonds, circRNAs. Blue squares, miRNAs. Green ovals, mRNAs.

Fig. 5

Gene Ontology (GO) enrichment and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis using predicted target messenger RNAs (mRNAs) of validated four circular RNA (circRNA). (A) GO enrichment analysis comprises three categories: biological process (in blue), cellular component (in green) and molecular function (in orange). Top five terms of each category are displayed. (B) KEGG pathway analysis shows top 10 terms that may be involved in the regulatory network mediated by differentially expressed circRNAs in type 1 diabetes mellitus.

Table 1.

Characteristics and demographics of human participants

Characteristic Microarray screening
qRT-PCR validation
Control (n=3) T1DM (n=3) Control (n=12) T1DM (n=12)
Age, yr 16 (12–18) 15 (12–17) 12.6 (6–18) 12.8 (5–17)
Gender, female/male 0/3 2/1 6/6 7/5
Disease duration, mo - 3.7±2.1 - 3.0±1.9
HbA1c, % 5.6±0.1 9.7±0.8 5.4±0.4 9.9±1.2
Subjects positive for ≥1 islet Aba 0 3 0 12
Serum C-peptide during OGTT, ng/mLb
 0 min - 0.16±0.15 - 0.19±0.18
 60 min - 0.21±0.18 - 0.39±0.37
 120 min - 0.3±0.41 - 0.54±0.38

Values are presented as median (range) or mean±standard deviation.

qRT-PCR, quantitative real-time polymerase chain reaction; T1DM, type 1 diabetes mellitus; HbA1c, glycosylated hemoglobin; Ab, autoantibody; OGTT, oral glucose tolerance test.

a

Number of subjects who were positive for at least 1 of the following islet autoantibodies: glutamic acid decarboxylase autoantibodies (GADA), insulinoma associated-2 autoantibodies (IA-2A), and zinc transporter-8 autoantibodies (ZnT8A),

b

Serum C-peptide were measured 0 (fasting level), 60 and 120 minutes during OGTT.