Regulatory network of two circRNAs and an miRNA with their targeted genes under astilbin treatment in pulmonary fibrosis

Abstract Circular RNAs (circRNAs) are becoming new therapeutic drug targets. However, their profiles under astilbin treatment have not been reported yet. In this study, we analysed the global reprogramming of circRNA transcriptome and a regulatory network of circRNAs with their targeted genes under astilbin treatment in pulmonary fibrosis. A total of 145 circRNAs were differentially expressed in the astilbin‐treated group compared with the bleomycin‐treated group using RNA sequencing. In the bleomycin‐ and astilbin‐treated groups, 29 coexpressed circRNAs were found. The maximum number of circRNAs was distributed on chromosome two, and their length varieties were mainly within 1000 bp. Four differentially expressed circRNAs (circRNA‐662, 949, 394 and 986) were tested to validate the RNA sequencing data, and their targeted microRNAs and genes were analysed by qRT‐PCR, Western blot, Pearson correlation coefficient, a dual‐luciferase reporter system and anti‐AGO2 RNA immunoprecipitation. The results showed that circRNA‐662 and 949 can act as “miR‐29b sponges” targeting Gli2 and STAT3 to exert their functions. Our work suggests that the transcriptome complexity at the circRNA level under astilbin treatment. These circRNAs may be potential molecular targets for drug action.

Astilbin has attracted considerable attention because of its antiinflammatory and anti-oxidative properties. 6,7 It can alleviate psoriasis-like skin lesions by inhibiting Th17-related inflammation 8 and ameliorate cisplatin-induced nephrotoxicity by reducing ROS accumulation. 9 In addition, astilbin can promote angiotensin-converting enzyme inhibition in vitro 10 and protect diabetic rat heart against ischaemia-reperfusion injury. 11 In a previous study, we reported that astilbin ameliorates fibrosis via blockade of the Hedgehog signalling pathway, in which we mainly focused on the role of coding RNAs under astilbin treatment. 12 However, high-throughput sequencing technology on RNA transcriptomes revealed that most of these transcripts are not translated to protein which are named as non-coding RNAs (ncRNAs), 13 such as microRNA (miRNA), long ncRNA (lncRNA) and circular RNA (circRNA). Moreover, increasing evidence revealed that these ncRNAs play diverse roles in the development and progression of human diseases and are becoming new therapeutic drug targets. [14][15][16] However, the changes of ncRNAs under astilbin treatment have not been reported yet.
In our previous study, we found that ncRNAs play key roles in fibrosis and can be used as therapeutic targets for pulmonary fibrosis using microarray technology. [17][18][19] However, microarray technology lacks sensitivity for low-abundance RNA and may miss RNAs related to disease risk and progression. Moreover, compared with miRNA and lncRNA, circRNA is a novel class of RNA transcript and a new clinical diagnostic and prognostic biomarker. 20,21 Therefore, in this study, we first combined RNA sequencing and other experiments to analyse the global reprogramming of circRNA transcriptome under astilbin treatment in pulmonary fibrosis. We found that different types of RNAs can collaborate to establish a sophisticated regulatory network to exert their function. We hope by characterizing the interactions among mRNAs, miRNAs and circRNAs can enhance our understanding of disease pathogenesis and drug action, and create new theories for the pathogenesis and treatment of pulmonary fibrosis.

| Animal model
Mice with mean weight of 25 ± 2 g were purchased from the Model Animal Research Center of Nanjing University (Nanjing, China). Experiments on the animals studied were approved by the Committee on the Ethics of Animal Experiments of Binzhou Medical University. All animals were bred and maintained in a 12 hours light/dark cycle and allowed free access to food and water, and randomly divided into four groups: sham group, 40 mg/ kg astilbin alone group (AST group), bleomycin-treated group (BLM group) and BLM + 40 mg/kg AST-treated group. BLM animal model was administered with 5 mg/kg BLM dissolved in saline via single intratracheal instillation under anaesthesia as previously described. 12 At day 14, astilbin was orally administered every day. At day 30, all mice were killed, and lung tissue sections were collected and immediately frozen in liquid nitrogen for further studies.

| Cell model and treatment
The mouse lung fibroblast cell line L929 was purchased from the Cell Bank of the Chinese Academy of Sciences. The cells were maintained in advanced minimum essential medium and supplemented with 10% newborn calf serum, 100 U/mL penicillin and 100 µg/mL streptomycin at 37°C under a humidified atmosphere of 5% CO 2 and 95% air. The cells were first administered with 5 ng/mL transforming growth factor beta 1 (TGF-β1) for 6 hours and then co-treated with 40 µg/mL AST for 48 hours as previously described. 12
Transverse sections of 4 μm thickness were stained with H&E or Masson's trichrome following the manufacturer's standard protocol. Three sections collected from each lung were analysed in the experiment.

| RNA sequencing
Total RNA was obtained from sham group, BLM-treated group and astilbin-treated group, respectively. mRNA was purified from total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in NEBNext First Strand Synthesis Reaction Buffer (5×). These samples were used for sequencing. Sequencing libraries were generated using NEBNext ® Ultra ™ RNA Library Prep Kit for Illumina ® (#E7530L, NEB) following the manufacturer's recommendations and index codes were added to attribute sequences to each sample. Following quantile normalization of the raw data, circRNA with at least two out of two present or marginal flags were selected for further data analysis.
Differentially expressed circRNAs were identified through fold change filtering.

| CircRNA identification pipeline
CircRNAs were predicted by CIRCexplorer. 22 CIRCexplorer utilized the unmapped reads exported by TopHat and aligned these reads using TopHat-Fusion. 23 The alignment results were imported into CIRCexplorer to predict circRNAs. All circRNA candidates with junction reads <2 were discarded. The remaining circRNA candidates were considered bonafide circRNAs.

| Gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) analysis
The main functions of the target circRNAs were determined according to GO analysis. Pathway analysis was performed to determine the significant pathway of the differential genes according to KEGG.
Fisher's exact test and chi-squared were used to classify the GO category and select the significant pathway. The significance threshold was defined by P-value and false discovery rate.

| Quantitative real-time reverse transcriptase polymerase chain reaction (qRT-PCR)
Approximately 2 mg RNA was extracted from lung tissue samples for the synthesis of first-strand cDNA. qRT-PCR analysis was per- Primers used in qRT-PCR are shown in Table 1. 20  were then washed with TBST, incubated with ECL reagent and exposed. Then, membranes were subsequently stripped and reprobed with glyceraldehyde 3-phosphate dehydrogenase antibody (GAPDH), which served as loading control.

| Statistical analysis
Data are expressed as the mean ± standard deviation (SD), and statistical significance was assessed by one-way ANOVA with Student-Newman-Keuls post-test for experiments comparing three or more groups. Probability values of <0.05 were considered significant.
Statistical analyses were performed using SPSS version 19.0 software.

| AST ameliorated pulmonary fibrosis in vivo and in vitro
To assess circRNA expression under astilbin treatment, H&E and Masson's trichrome staining were first used to evaluate the fibro-

| Profiles of differentially expressed circRNAs
A computational pipeline for the systematic identification of circR-NAs is shown in Figure 2. The results of heat map analysis showed that the circRNAs are differentially expressed. A total of 145 circRNAs were differentially expressed in the AST-treated group compared with the BLM-treated group. Of these circRNAs, 58 were up-regulated, and 87 were down-regulated ( Figure 3A). Then, we further analysed the distribution of these circRNAs on chromosomes. The circRNAs were distributed on chromosomes 1-19. The maximum number of circRNAs was mainly distributed on chromosome two ( Figure 3B). Their length varieties widely varied, and most of them were within 1000 bp ( Figure 3C).

F I G U R E 4
Analysis of circRNAs' host gene function. A, Biological process, molecular function, and cell component of these host genes were analysed using the GO database. B, Biological pathways of these host genes were analysed using the KEGG database

| Analysis of the differentially expressed circRNAs' host genes
One of circRNA mechanisms is that circRNA exert its regulatory function by specifically interacting with its host genes. Therefore, we analysed the host genes of differentially expressed circRNAs by using GO and KEGG analysis. According to the GO annotation tool, these host genes were mainly enriched for GO terms related to biological adhesion, developmental process, cell junction, extracellular region, signal transducer activity and structural molecule activity ( Figure 4A). Biological pathways associated with the host genes of differentially expressed circRNAs were analysed using the KEGG database. These pathways were associated with cGMP-PKG, MAPK and vascular smooth muscle contraction ( Figure 4B).

| Analysis of differentially expressed circRNAs' targeted genes
The other circRNA mechanism is that circRNA act as "miRNA sponges" to regulate the targeted genes, namely ceRNA. Regulatory network of two circRNAs and an miRNA with their targeted genes under astilbin treatment was explored in detail.
Firstly, we analysed the differentially coexpressed circRNAs between BLM vs sham group and AST vs BLM group to validate the RNA-sequencing data. In total, 29 differentially coexpressed circRNAs were obtained in the BLM and AST groups ( Figure 5A).
Therefore, the miR-29b-2-5p targeted genes were predicted based on TargetScan, MiRanda data and miRBase ( Figure 7A). Among the miR-29b-2-5p targeted genes, we focused on STAT3 and Gli2 because they are the main regulatory factors of lung fibrosis in our previous study. 12,18 STAT3 and Gli2 decreased under astilbin treatment compared with BLM treatment ( Figure 7B). Their expressions were inversely correlated with miR-29b-2-5p ( Figure 7C). A dual-luciferase reporter system was conducted to confirm the identity of their direct target. The results showed that the results showed that STAT3 is the direct target of miR-29b-2-5p ( Figure 7D).

| D ISCUSS I ON
In this study, we described 145 differentially expressed circRNAs and analysed the relationships between these circRNAs and their protein-coding genes under astilbin treatment in pulmonary fibrosis.
We found that (a) circRNAs have important functions in regulating gene expression in pulmonary fibrosis and affect drug activity and that (b) various types of RNAs and their protein-coding genes can collaborate to establish a sophisticated regulatory network.
F I G U R E 6 Expression of the targeted miR-29b-2-5p. A, Binding sites of circRNA-394, circRNA-662, circRNA-949 and circRNA-986 for miRNAs. B, Expression of miR-29b-2-5p by using RNA sequencing and qRT-PCR in vivo and in vitro. C, Inverse relationship of circRNA-949 and circRNA-662 with miR-29b-2-5p was analysed by Pearson correlation coefficient. D, CircRNA-949 and circRNA-662 are the direct targets of miR-29b-2-5p by using a dual-luciferase reporter system experiment. E, RIP data confirmed the result of the dual-luciferase experiment. Each bar represents the mean ± SD, n = 6, *P < .05, **P < .01 With the development of high-throughput sequencing, emerging evidence has indicated that ncRNAs play key roles in fibrosis.
To date, most studies on the role of ncRNAs during pulmonary fibrogenesis have focused on miRNAs and lncRNAs, which regulate mRNA translation via the RNA interference pathway. 25 Pandit et al reported that let-7d inhibition is a key regulatory event in preventing lung fibrosis. 26 Jiang et al summarized the role of lncRNAs in the process of lung fibrosis. 27 However, the occurrence of cir-cRNAs in fibrosis remains largely unknown. CircRNAs are derived from transcripts that are back-spliced and joined head to tail at the splice sites. 28,29 Their covalently closed ring structure make them highly stable, even exceeding the stability of the housekeeping gene GAPDH. 30 Bachmayr-Heyda et al discovered a negative correlation of circRNA abundance and proliferation in lung fibrosis. 31 27 Based on our previous studies, competing endogenous RNAs are present not only in lncRNAs but also in F I G U R E 7 Expression of the targeted genes STAT3 and Gli2. A, MiR-29b-2-5p targeted genes were predicted. B, Western blot analysis showed that STAT3 and Gli2 protein levels decreased after astilbin treatment. C, Inverse relationships of miR-29b-2-5p with STAT3 and Gli2 were analysed by Pearson correlation coefficient. D, The result of the dual-luciferase experiment showed that STAT3 was the direct target gene of miR-29b-2-5p. Each bar represents the mean ± SD, n = 6, *P < .05, **P < .01 circRNAs. 19,43 In this study, we showed the existence of competing  47 Lee et al analysed the transcriptome complexity through RNA sequencing in normal and failing murine hearts. 48 However, few studies have identified the profiles of differentially expressed circRNAs under drug action using RNA sequencing in progressive lung fibrosis. Our findings regarding the crosstalk and molecular consequences of these circRNAs provide candidate therapeutic targets for drug treatment. From a clinical perspective, the targeting of circRNAs as a novel therapeutic approach will require a deeper understanding of their function and mechanism of action.
However, in the short term, changes in circRNA expression are likely to be used as biomarkers for disease stratification and/or assessment of drug action.

CO N FLI C T O F I NTE R E S T
The authors declare no conflicts of interest.

AUTH O R CO NTR I B UTI O N S
XD Song, ME Li, P Xu, JJ Zhang and CJ Lv participated in the conception, hypothesis and design of the study. GP Lu and WB Liu performed the experiments. GH Cao, XY Liu and XL Zhang carried out the statistical analyses. All authors contributed to interpretation of the data. XD Song wrote the manuscript and all authors made critical revisions. All authors read and approved the final manuscript.

DATA ACCE SS I B I LIT Y
The analysed data sets generated during the study are available from the corresponding author on reasonable request.