{'fraction_non_alphanumeric': 0.035786423114077454, 'fraction_numerical': 0.030385527355756566, 'mean_word_length': 5.21418826739427, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 2, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
CircRNA expression profiles and circRNA-miRNA-mRNA crosstalk in allergic rhinitis☆

Background Circular RNAs (circRNAs) are involved in inflammation; however, their role in allergic rhinitis (AR) remains unclear. In this study, we analyzed circRNA expression and identified a circRNA-miRNA-mRNA network through which circRNAs regulate AR pathogenesis. Methods We analyzed circRNA, miRNA, and mRNA expression profiles in the nasal mucosa by high-throughput sequencing (HTS), using a fold-change >1.5 and p-value < 0.05 to pinpoint significantly differentially expressed (DE) circRNAs, miRNAs, and mRNAs in AR. A DEcircRNA-DEmiRNA-DEmRNA crosstalk network was then constructed using bioinformatics and statistical analysis. Gene ontology and Kyoto encyclopedia of genes and genomes pathway analyses were performed to identify the biological terms enriched in the network; whereas RT-PCR and Sanger sequencing were used to confirm the circRNAs. Results A total of 264 DEcircRNAs were identified by HTS, including 120 upregulated and 144 downregulated in AR compared to controls. A DEcircRNA-DEmiRNA-DEmRNA crosstalk network was constructed with 17 miRNAs, 11 circRNAs, 29 mRNAs, and 64 interaction pairs. These genes were involved in the Wnt signaling pathway, TNF biosynthesis, inflammatory responses, the PI3K-Akt signaling pathway, and Toll-like receptors. Of the 11 DEcircRNAs, hsa_circ_0008668 and circTRIQK were upregulated, whereas hsa_circ_0029853 and circRNA_01002 were downregulated in AR tissues. Sanger sequencing confirmed the back-splicing junctions of these circRNAs. Conclusions We constructed a novel DEcircRNA-DEmiRNA-DEmRNA network for AR that provides a basis for future studies to investigate its underlying molecular mechanisms.


INTRODUCTION
Allergic rhinitis (AR) is a chronic multi-factorial inflammatory disease, with a prevalence rate of 17.6% in adults that is expected to increase yearly in China. 1,2 Despite the discovery of some genetic factors responsible for AR progression, the underlying molecular mechanisms have not yet been fully elucidated.
Previous studies have demonstrated that microRNAs (miRNAs) and long non-coding RNAs (lncRNAs) play critical roles in AR 3,4 and circular RNAs (circRNAs), a newly-discovered family of endogenous non-coding RNAs, have been found to exert a variety of functions. 5 For instance, largescale studies have confirmed that circRNAs are involved in neuronal function, tumor development, and innate immune responses, 6 while others are associated circRNAs with inflammatory diseases. 7, 8 Fang et al 9 have reported that circANKRD36 expression positively correlates with the level of inflammatory factor IL-6 and may be involved in type 2 diabetes mellitus and inflammation-related pathways. Moreover, Liu et al 10 have found circRNA Ddx17 can exert protective effects by inhibiting miR-17-5p expression in ovalbumin-induced AR mice. Similarly, Zhu et al 11 have found circHIPK3 is highly expressed in AR mucosa, and that circHIPK3 can upregulate Th2 cell-specific transcript factor GATA-3 via modulating miR-495. Wang et al 12 have found circARRDC3 promotes the development of AR by regulating the miR-375/KLF4 axis.
Competing endogenous RNA (ceRNA) networks are able to reveal previously unknown functions of the human genome in diseases. 13 For instance, circRNAs have been shown to regulate mRNAs by sponging miRNA, ie, competitively binding to miRNAs via common miRNA response elements (MREs). 14 Zhang et al 8 have constructed a network of differentially expressed (DE) circRNAs, miRNAs, and mRNAs, with gene enrichment analysis revealing that they regulate cell adhesion, cell activation, and immune responses in atherosclerosis. CircRNA-related ceRNA networks have also been established for illustrating the pathogenesis of various diseases, including liver cancer, 15 glioma, 16 and osteoarthritis. 17 However, no such networks have been elucidated for AR. Therefore, in this study, we described circRNA, miRNA, and mRNA expression profiles in the nasal mucosa from patients, with and without AR, using high-throughput sequencing (HTS), and constructed a DEcircRNA-DEmiRNA-DEmRNA crosstalk network to illustrate the potential molecular mechanisms of AR progression.

Clinical data and nasal mucosa tissues
The nasal mucosa tissues were collected from a cohort of 20 adult subjects (10 with perennial AR and 10 non-allergic controls) aged 18-53 years who received endoscopic septoplasty from August 2019 to January 2020. Tissue samples were immediately frozen and stored in liquid nitrogen, and then subjected to transcriptomic profiling analysis. AR was diagnosed according to the clinical guidelines. 1 Serum total IgE and allergenspecific IgE were measured by the Automated Immuno-Strip Analysis System (LG Life Sciences Ltd., Korea) in all the subjects. Patients with perennial AR enrolled were sensitized mainly by house dust mites, including Dermatophagoides pteronyssinus and/or Dermatophagoides farinae (Table S1). The controls had negative IgE test results in sera. All subjects had ceased treatment with glucocorticoids, antihistamines, antileukotrienes, and other immunomodulating drugs for at least 4 weeks prior to surgery. Subjects who had an immune deficiency, upper respiratory tract infection, chronic rhinosinusitis, nasal polyps, sinonasal neoplasms, asthma and other allergic diseases, and a history of nasal surgery or immunotherapy were excluded.
Total RNA isolation, library construction, and sequencing Total RNA was extracted from nasal mucosa samples using TRIzol reagent (Life Technologies, CA, USA) and its quantity and quality were measured using a NanoDrop ND-1000 (Thermo Scientific, DE, USA). RNA integrity was assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, CA, USA) and RNAs with an integrity value of !7 were retained.
Six randomly selected samples (three AR and three controls) from the 20 subjects were used to construct circRNA, miRNA, and mRNA libraries, and the 6 samples and the remaining 14 samples were sent to reverse transcription-polymerase chain reaction (RT-PCR). Approximately 1 mg of RNA was extracted from each sample, and ribosomal and linear RNAs were depleted from the samples using Ribo-zero (Illumina, CA, USA) and RNase R (Epicenter, CA, USA), respectively, according to the manufacturers' instructions. CircRNA libraries were constructed using a TruSeq Stranded Total RNA HT/LT Sample Prep Kit (Illumina) according to the manufacturer's instructions. Library quality was assessed using an Agilent 2100 system (Agilent Technologies) and the libraries were sequenced using an Illumina sequencing platform (HiSeqTM 2500).
RNA (1 mg) from each of the 6 samples was used as input material for small RNA library generated using TreSeq Small RNA Sample Prep Kits (Illumina), according to the manufacturer's instructions. Library quality was assessed using an Agilent 2100 system and libraries were sequenced using an Illumina sequencing platform (HiSeqTM 2500). All sequencing procedures and analyses were performed by OEbiotech (Shanghai, China).

Identification of miRNA-target interactions
MiRNA-mRNA and miRNA-circRNA pairs were selected using miRanda, which predicts miRNA targets based on sequence complementarity, conserved target sites, and free formation energy. Next, we filtered out miRNAs with crosstalk with circRNAs and mRNAs and used the 3 types of RNAs to produce a circRNA-miRNA-mRNA network. To identify circRNA-mRNA pairs, we used a ceRNA score that evaluates the significance (p-value) of miRNAs shared by each circRNA and mRNA as the ratio of common miRNA MRE to circRNA-miRNA binding MRE, as follows: where M1 represents the quantity of all miRNAs, M2 represents the quantity of miRNAs regulating mRNAs, M3 represents the quantity of miRNAs with regulatory effects on the circRNA, and M4 represents the quantity of common miRNAs. pvalues < 0.05 were considered statistically significant.

Bioinformatics and data analysis
Gene ontology (GO) analysis was used to annotate gene function. Enrichment scores were calculated as the negative logarithm of the p-value and used to determine the statistical significance of GO term clusters targeted by DE genes. Kyoto encyclopedia of genes and genomes (KEGG) analysis was used to determine the biological pathways involving genes in the network, as quantified by the enrichment score.

Sanger sequencing and RT-PCR
Sanger sequencing was performed by Tsingke Biotechnology (Nanjing, China) to confirm the back-splicing junctions of each circRNA. Sanger sequencing data that matched the HTS results and displayed a unimodal dissolution curve were retained. To ensure that each PCR amplification fragment contained junction sites, the designed circRNA primers were divergent, instead of the frequently used convergent primers (Table S2). The experiment was repeated 3 times for each sample from the 20 subjects (10 with perennial AR and 10 controls). The RT-PCR dissolution curve was unimodal, proving that each pair of specific reverse primers had good repeatability and a single product.

Statistical analysis
All statistical analyses were carried out using SPSS 20.0 (IBM SPSS, Chicago, IL, USA). Betweengroup differences were analyzed using Student's t-tests. CircRNA expression was displayed as a fold-change calculated using the 2 -OOct method on the RT-PCR data. Data from at least 3 independent experiments were presented as the mean AE standard deviation (SD). p-values < 0.05 were considered statistically significant.

CircRNA expression and characteristics in the nasal mucosa
Three pairs of nasal mucosa tissues from 10 patients with AR and 10 controls were randomly selected for HTS analysis. The clinical characteristics of the subjects are shown in Table S1. All  Some circRNAs were annotated using circBase. Novel circRNAs were described in terms of their genome location. Green and red nodes represent circRNAs and miRNAs, respectively. circRNA, circular RNA; miRNA, microRNA; DE, differentially expressed normalized sample datasets were found to display similar circRNA distribution densities (Fig. 1A). CIRI software 18 was used to identify and annotate 30,936 circRNAs from the 2 groups, including 12 376 known and 18 560 novel (Fig. 1B), and widely distributed circRNA transcripts were found on all chromosomes (Fig. 1C). CircRNA transcript expression was quantified using the spliced reads per million (RPM) distribution density, revealing that almost all circRNA transcripts were expressed at low levels (Fig. 1D).
Among the DEcircRNAs, 1119 were antisense and 29 335 were sense, including 28 694 exonic, 641 intronic, and 1119 intergenic circRNAs. In particular, 93% of the circRNAs were sense exons and most circRNA transcripts carried 1-7 exons, with a length of 200-400bp or >2000bp. Hierarchical clustering revealed DE circRNAs between samples ( Fig. 2A) and HTS analysis identified 264 DEcircRNAs, including 120 upregulated and 144 downregulated in AR. Based on the number of possible miRNAbinding sites, structure, energy score and pvalue, the top 300 relevant circRNA-miRNA pairs were screened and their network illustrated using Cytoscape software (Fig. 3).

MiRNA and mRNA expression profiles in the nasal mucosa
Next, we used HTS to analyze the miRNA and mRNA expression profiles in the nasal mucosa tissues from the AR patients and controls. Of the 1784 miRNAs available, 13 were upregulated and 15 were downregulated. On the other hand, of the 2030 mRNAs, 93 were upregulated and 79 were downregulated. Hierarchical clustering revealed that these DEmiRNAs (Fig. 2B)
To determine the biological roles of the DEc-ircRNAs, we performed GO and KEGG analyses on the DEmRNAs in the DEcircRNA-DEmiRNA-DEmRNA crosstalk network. These genes converged in 159 terms, including 172 biological processes, 47 cell components, and 46 molecular functional components. The GO terms mainly involved the Wnt signaling pathway, lung epithelial cell development, and TNF biosynthesis (biological processes), endoplasmic reticulum lumen, flotillin complex, and intercellular canaliculus (cell components), fibronectin binding, collagen binding, and calcium: sodium antiporter activity (molecular functions; Fig. 6A). In addition, these genes were enriched in the inflammatory response and various signaling pathways (such as cell surface receptor signaling pathways), indicating that circRNAs may participate in multiple pathways related to AR pathogenesis. KEGG pathways analysis revealed 24 significantly involved pathways, of which local adhesion, ECM receptor interaction, protein degradation and absorption, and PI3K-Akt signaling pathways were the most enriched (Fig. 6B). However, other pathways were related to AR, such as the Toll-like receptor (TLR) signaling pathway, suggesting that circRNAs may Fig. 5 Depiction of the DEcircRNA-DEmiRNA-DEmRNA crosstalk network. The network includes 17 miRNAs (red nodes), 11 circRNAs (green nodes), and 29 mRNAs (yellow nodes). CircRNAs and miRNAs were annotated using circBase and miRbase, respectively. Novel circRNAs and miRNAs were described in terms of their genome location. circRNA, circular RNA; miRNA, microRNA; DE, differentially expressed act as ceRNAs in AR to regulate mRNAs by competitively binding to miRNAs.

Sanger sequencing verification and RT-PCR quantification of circRNAs
The DEcircRNA-DEmiRNA-DEmRNA crosstalk network included 11 circRNAs (Table 1), of which two were known (hsa_circ_0008668 and hsa_circ_0029853) and nine were novel. To verify these circRNAs, hsa_circ_0008668, hsa_circ_0029853, circTRIQK (a senseoverlapping circRNA), and circRNA_01002 (an intergenic circRNA) were subjected to RT-PCR amplification and Sanger sequencing to confirm their back-splicing sites (Fig. 7A-D). The RT-PCR dissolution curve had a single peak, suggesting that the primers had good repeatability and yielded a single product type. In addition, RT-PCR revealed that hsa_circ_0008668 (Fig. 7E) and circTRIQK (Fig. 7G) were upregulated in AR, whereas hsa_circ_0029853 (Fig. 7F) and circRNA_01002 (Fig. 7H) were downregulated. These findings were consistent with the HTS data and suggested that the HTS results are highly reliable.

DISCUSSION
Recent epigenetic research has yielded considerable advances in our understanding of AR; however, these studies focused on protein-coding genes, miRNAs, and lncRNAs. Research on circR-NAs has shown that they can act as new diagnostic markers for cancers. 19 Indeed, circRNAs can function as ceRNAs to sponge miRNAs 14 and thereby manipulate transcription and alternative splicing, 20 cooperate with RNA-binding proteins, 21 and affect protein translation. 22 CircRNAs can also control the cellular levels of specific miRNAs, 23 the majority of which harbor highly conserved sequences and act in multiple cellular processes via post-transcriptional regulation. Thus, disturbed miRNA activity can be partly explained by circRNA involvement. 24 With the help of MREs, circRNAs, miRNAs and mRNAs participate in large-scale ceRNA crosstalk regulating post-transcriptional genes in various manners. Previous studies have found that Ddx17, circHIPK3 and circARRDC3 can function as ceRNAs to sponge miRNAs in AR. [10][11][12] However, each of these studies has just proved one circRNA related to immunity regulation. To the best of our knowledge, our study is the first to investigate circRNA expression profiles and ceRNA network in AR.
We predicted DEmiRNA-DEcircRNA interactions using miRanda and screened 116 pairs with negative expression relationships, including 17 miRNAs (such as hsa-miR-98-5p and hsa-miR-455-5p), while HTS revealed that hsa-miR-98-5p was significantly upregulated. Takuse et al 25 have found that hsa-miR-98-5p negatively correlates with peripheral IL-10 þ cell count and IL-10. Moreover, an increasing number of studies have demonstrated imbalance between CD4 þ T cells and B cells in peripheral blood plays a role in AR. 26 Moreover, studies have reported reduced numbers of regulatory T cells (Treg) 27 and regulatory B cells (Breg) 28,29 in AR, which are major sources of IL-10, an important molecule in allergen tolerance. 30,31 Hence, hsa-miR-98-5p may affect allergen tolerance in AR by negatively regulating Tregs, Bregs, and IL-10 secretion. Our DEmiRNA-DEcircRNA network revealed that hsa-miR-98-5p competed with several circRNAs, such as hsa_circ_0029853, suggesting that circRNAs may interact with hsa-miR-98-5p during AR immune tolerance. Back-splice junction sequences for hsa_circ_0008668, hsa_circ_0029853, circTRIQK, and circRNA_01002 were confirmed by Sanger sequencing, respectively. Arrows represent divergent primers binding to the genomic region of circRNAs. (E-H) RT-PCR results for circRNAs in the nasal mucosa of ten patients with AR and ten control subjects. Hsa_circ_0008668 (E) and circTRIQK (G) were upregulated in AR. whereas hsa_circ_0029853 (F) and circRNA_01002 (H) were downregulated, consistent with the HTS data. Data represent the mean AE SD of at least three independent experiments.*p < 0.05, **p < 0.01, ****p < 0.0001. circRNA, circular RNA Similarly, we identified 379 DEmiRNA-DEmRNA pairs and many DE genes (such as DOCK8 and MUC19) that were significantly associated with AR. In particular, DOCK8 expression was decreased in the nasal mucosa of patients with AR. DOCK8 deficiency plays a key role in allergic inflammation, susceptibility to infection, autoimmunity, and malignant tumors 32,33 as it exerts negative effects on the migration, function, and survival of immune cells as well as innate and adaptive immunity. 32 For instance, DOCK8 participates in the migration of interstitial dendritic cells 34 and maintains protective immunity by repressing the function of RORgt þ ILCs cells. 35 DOCK8 has also been shown to regulate adaptive immune cells, while Janssen et al 33 have found that DOCK8 is related to Treg cell homeostasis and function and enhances peripheral B cell tolerance. In addition, DOCK8 modulates MyD88-dependent TLR signal transduction and activates the STAT3 transcription factor. 32 By examining DEmiRNA-DEmRNA relationships, we found that DOCK8 could bind to hsa-miR-548aq-3p, hsa-miR-582-5p and hsa-miR-96-5p, and that its expression correlated negatively with that of these miRNAs. We also found that MUC19, SLC8A3 36 and SPP1 37 in the ceRNA network were associated with AR, while PTPN22 and FCRL5 expression was downregulated in the nasal mucosa. Since single nucleotide polymorphisms in PTPN22 38 and FCRL5 39 is associated with asthma and AR in Chinese children, future studies should focus on the mechanisms via which DE genes regulate AR.
Finally, we constructed a novel DEcircRNA-DEmiRNA-DEmRNA network for AR consisting of 64 ceRNA pairs, including 11 circRNAs bound to 17 miRNAs, such as hsa-miR-98-5p, to regulate the expression of 29 genes. Previous studies have found that AR is regulated by many signaling pathways, including the PI3K-Akt 40 and TLR 41 signaling pathways. Consistent with these findings, we found that the 29 AR-related genes were involved in various important biological functions and pathways, including the PI3K-Akt and TLR signaling pathways, suggesting that these pathways play key roles in AR development.
The DEcircRNA-DEmiRNA-DEmRNA network also revealed that hsa_circ_0029853 could competitively bind to hsa-miR-98-5p, hsa-miR-140-5p, and hsa-miR-3065-5p to regulate ERRFI1 and NRCAM (Fig. 5), whose expressions were significantly downregulated in the nasal mucosa of patients with AR. ERRFI1 is also downregulated in human breast, skin, pancreatic, and ovarian cancers. 42 Further, Izumchenko et al 43 have found that ERFFI1 can significantly inhibit LPSinduced TNF-a, IL-1b, and IL-6 mRNA and protein expression. Moreover, our previous studies revealed that these pro-inflammatory cytokines were upregulated in the nasal mucosa of patients with AR, 41 suggesting that their overexpression may be related to ERRFI1 downregulation in these patients. In addition, ERRFI1 has been implicated in the epidermal growth factor receptor (EGFR) signaling pathway 43 and studies have demonstrated an interplay between EGFR and eosinophilic airway inflammation. 44,45 Moreover, intranasal administration of the EGFR inhibitor AG1478 has been found to significantly suppress GM-CSF, eotaxin-1, and RANTES secretion by nasal mucosal cells. 45 Since ERRFI1 is a negative regulator of EGFR, 43 it may participate in AR pathogenesis by regulating the EGFR signaling pathway in the nasal mucosa. Additionally, ERFFI1 is regulated by multiple miRNAs, including miR-126, 46 miR-205 47 and miR-2355-5p 43 ; therefore, hsa_circ_0029853 may sponge hsa-miR-98-5p, hsa-miR-140-5p, and hsa-miR-3065-5p to inhibit ERFFI1 expression and thereby promote AR progression.
There are some limitations in the present study. First, we did not perform a functional analysis for the miRNAs (such as hsa-miR-98-5p, hsa-miR-455-5p) and genes (such as DOCK8, ERRFI1), especially their ceRNA network, to determine their roles in AR pathogenesis. The second is that we did not explore the mechanism of circRNAs in sponging miRNAs to regulating mRNAs, or other functions of circRNAs, such as acting as coding protein.
In conclusion, circRNA expression profiles and the ceRNA network for AR could provide a basis for future studies on the underlying molecular mechanisms of AR.
Ethics approval and consent to participate This study was approved by the Ethics Committee of the First Affiliated Hospital of Nanjing Medical University (Approval No. 2019-SR-023), and all the participants have signed written informed consent.

Availability of data and materials
The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

Consent for publication
All authors have approved the manuscript for submission.
Author contributions C-Y Qiu and X-Y Cui performed data analysis and wrote the draft manuscript. M-P Lu, M Yin, W-Y Xu, and X-J Zhu collected samples and performed specimen quality screening. Q Yang performed RT-PCR. L Cheng designed the study and finalized the manuscript. All authors approved the final version of the manuscript.

Funding
This study was supported by a grant (No. 81970862) from the National Natural Science Foundation of China.

Potential competing interests
The authors report no competing interests.

Confirmation of unpublished work
A part of data and materials in this manuscript has been posted to the AUTHOREA preprint server at https:// authorea.com/doi/full/10.22541/au.160193459.95947821/ v1. The authors confirm that the final version of the manuscript is original, has not been published before, and is not currently being considered for publication elsewhere.