Comprehensive circular RNA expression profiles in a mouse model of nonalcoholic steatohepatitis.

Recent studies have revealed that circular RNAs (circRNAs) are involved in the development of various liver diseases. However, the regulatory role of circRNAs in nonalcoholic steatohepatitis (NASH) has not been fully elucidated. In the present study, the circRNA profiles in a NASH mouse model were investigated, and their functions in NASH were predicted using bioinformatics analysis, with the aim of providing novel clues for delineating the mechanisms of action. A NASH mouse model was established by feeding mice with a methionine and choline‑deficient diet. The liver circRNA profile was screened using a circRNA microarray, and the differentially expressed circRNAs were verified by reverse transcription‑quantitative polymerase chain reaction. Subsequently, circRNA‑microRNA (miRNA) interactions were predicted. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were performed to annotate the biological functions of host linear transcripts of circRNA. A total of 450 circRNAs were revealed to be dysregulated, with 298 circRNAs upregulated and 152 circRNAs downregulated in the NASH model mice. circRNA_29981 was identified as a significantly differentially expressed circRNA. The results from the circRNA‑miRNA pathway interaction analysis revealed that circRNA_29981 was a potential regulator of hepatic stellate cell activation. The host linear transcripts were also analyzed, and the top 10 enriched GO entries and KEGG pathways were annotated. These findings suggested that circRNAs may be important regulators of NASH. Taken together, the results of the present study demonstrated that the circRNA profile in NASH may provide potential candidates for future studies aimed at elucidating the pathogenic mechanism(s) involved in the disease.

Introduction

Non-alcoholic fatty liver disease (NAFLD) is the most commonly occurring type of chronic liver disease in developed and developing countries (1,2). The prevalence of this disease has considerably increased in recent years. NAFLD has been reported as the manifestation of metabolic syndrome in the liver (3,4). Nonalcoholic steatohepatitis (NASH), which is a pivotal stage in NAFLD, is more likely to progress to fibrosis, cirrhosis and hepatocellular carcinoma (HCC) in certain cases (5,6). Currently, NASH is the third most common indicator of liver transplantation in the USA, although it is predicted to become the most common one between 2020 and 2025. Concerning the mechanism of NASH, the so-called ‘two hit theory’ is the current hypothesis. The first hit is insulin resistance, which causes accumulation of lipids in hepatocytes. The second hit involves increased oxidative stress, lipid peroxidation and release of proinflammatory cytokines, which lead to the occurrence and development of NASH (7,8). However, the exact pathogenic mechanism of NASH has not been fully elucidated.

Circular RNA (circRNA) is a recently discovered type of endogenous noncoding RNA molecule with a closed loop structure formed through covalent bonding, without the 5′-end cap and the 3′-end poly (A) tail, and is attracting a lot of attention in the RNA family as the latest research ‘hot spot’, following microRNAs (miRNAs) and long noncoding RNAs (lncRNAs) (9). CircRNAs are highly stable and conservative across species (10). Accumulating evidence has demonstrated that circRNAs may function as ‘miRNAs sponges’, which bind with and inhibit the function of miRNAs. CircRNAs are able to regulate gene expression at the transcriptional and post-transcriptional levels (11–13). An increasing number of studies have reported that circRNAs are closely associated with disease development, and this class of molecule is expected to become a novel target for disease diagnosis and treatment.

Even though the functions of the majority of circRNAs in various types of disease remain unclear (11), increasing evidence suggests that the expression profile of circRNAs is altered in numerous diseases, including neurological diseases, diabetes and cardiovascular diseases (14–17). Recently, a growing number of circRNAs have been reported to be dysregulated in liver diseases, and their clinical significance has been investigated (18–23). However, to date, very few studies have focused on the expression profile of circRNAs in NASH. Jin et al (24) investigated the expression profile of circRNAs in a methionine and choline-deficient (MCD)-induced NASH mouse model, and identified four circRNA-miRNA-mRNA pathways; furthermore, it was demonstrated that the circRNA profile in NASH is able to provide potential candidates for future mechanistic studies.

In the present study, changes in the circRNA expression profile in a well-developed NASH mouse model were investigated, with the aim of investigating the potential pathogenic mechanisms of circRNAs during NASH development. Combining circRNA microarray screening and reverse transcription-quantitative polymerase chain reaction (RT-qPCR) validation with bioinformatics methods, a novel circRNA-miRNA interaction pathway that is closely associated with the disease processes of NASH was identified, enabling us to predict the potential roles of circRNAs in NASH development.

Materials and methods

Animals

A total of 12 male C57BL/6 mice, aged 6 weeks, weighing 20–25 g, were purchased from the Shanghai Xipuer-Bikai Laboratory Animal Co., Ltd (Shanghai, China). All mice were housed in cages under standard conditions with an ambient temperature of 24°C under a 12 h light/dark cycle. Mice had free access to food and water (as detailed in the subsequent section). All experimental protocols were approved by the Animal Ethics Committee of the Eighth People's Hospital of Shanghai (Permit no. 2016-02-2).

NASH animal model construction

Animals were randomly divided into two groups: Control group (n=6) and NASH group (n=6). All mice received food and water ad libitum. Mice in the NASH group were provided with an MCD diet. Mice in the control group were fed the same diet, but with sufficient quantities of DL-methionine (3 g/kg) and choline bitartrate (2 g/kg). At the end of 8 weeks of administration of the diet, all the mice were anesthetized and sacrificed, and their livers were removed. Serum and liver samples were frozen in liquid nitrogen prior to biochemical analysis. A small piece of the liver sample was fixed in 4% paraformaldehyde (96 h at 4°C) for subsequent histological analysis.

The levels of serum triglycerides (TG) and total cholesterol (TC), and the activities of the liver-associated enzymes, alanine aminotransferase (ALT) and aspartate aminotransferase (AST), were estimated with an automatic biochemical analyzer (Roche/Hitachi P-800 Modular Analytics Diagnostic System; Roche Diagnostics, Indianapolis, IN, USA). Hepatic homogenates were used for determination of the TG and TC content using a kit provided by Nanjing Jianchen Bioengineering Institute (Nanning, China) according to the manufacturer's protocol. Tissue lipids were extracted with methanol/chloroform (1:2).

The liver samples removed from each mouse were fixed in 4% paraformaldehyde for 96 h at 4°C. Following dehydration in graded alcohol (once in 50, 70 and 95%, and twice in 100% for 2–3 h each) and embedding in paraffin wax, the sections were cut to a thickness of 4 µm, and stained with hematoxylin and eosin (H&E) (25). The liver histology was scored for ballooning (0–2), steatosis (0–3), and inflammation (0–3), and the sum of these scores was used to create the NAFLD Activity Score. Picrosirius red staining (0.1% Sirius Red in saturated picric acid stained for 75 min at room temperature) was used to assess collagen fibers in the liver tissues. Each section was assessed under 10×20 light microscopic fields. For lipid staining, 6 µm frozen liver sections were air-dried for 30 min, followed by fixation in 4% formaldehyde. Oil red O staining was performed with staining solution (0.5% Oil Red O in anhydrous isopropanol, sections were stained for 30 min at 37°C) (Sigma-Aldrich; Merck KGaA, Darmstadt, Germany) according to the manufacturer's protocol. Then, six slides from the control and NASH groups were used for quantitative analysis. The stained area was quantified in six randomly selected ×200 microscopic fields per mouse section using ImageJ software (ImageJ2×; National Institutes of Health, Bethesda, MD, USA).

Isolation and quality control of RNA

Total RNA was isolated from the different groups using TRIzol® reagent (Invitrogen; Thermo Fisher Scientific, Inc., Waltham, MA, USA) and purified with an RNeasy Mini kit (Qiagen GmbH, Hilden, Germany) according to the manufacturer's protocol. The purity and concentration of RNA samples were determined with the NanoDrop® ND-1000 spectrophotometer (Thermo Fisher Scientific, Inc.). The integrity of the RNA was assessed by electrophoresis on a denaturing gel.

CircRNA microarray analysis

Three liver RNA samples from each group were selected for microarray studies. Sample labeling and array hybridization were performed according to the manufacturer's protocol (Arraystar, Inc., Rockville, MD, USA). In brief, the extracted total RNA was treated with RNase R (Epicentre®; Illumina, Inc., Madison, WI, USA) to remove linear RNAs and enrich the population of circRNAs. The enriched circRNAs were amplified and transcribed into fluorescent complementary RNA (cRNA) using random primers, according to the instructions provided by the Arraystar Super RNA Labeling kit (Arraystar, Inc.). Following purification of the labeled cRNAs with the RNeasy Mini kit (Qiagen GmbH), the concentration and specific activity of the labeled cRNAs were measured using the NanoDrop® ND-1000 spectrophotometer. The labeled cRNA was hybridized on to Arraystar mouse circRNA array v1.0. Following washing of the slides, the arrays were scanned using the Agilent G2505C scanner (Agilent Technologies, Inc., Santa Clara, CA, USA).

RT-qPCR validations

RT-qPCR was used to confirm the circRNA expression profiles obtained from the microarray data. Total RNA was isolated from liver tissues using TRIzol® reagent (Invitrogen; Thermo Fisher Scientific, Inc.) and reverse transcribed into cDNA using SuperScript™ III Reverse Transcriptase (Invitrogen; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. The relative gene expression was determined using the Agilent Mx3000P qPCR system. Amplifications were performed using a SYBR Green PCR kit (Takara Biotechnology Co., Ltd., Dalian, China). The qPCR program comprised an initial denaturation step of 3 min at 95°C, followed by 40 cycles of 15 sec at 95°C, 30 sec at 60°C, and 30 sec at 72°C. All samples were normalized to the signal generated from GAPDH. Results were obtained from three independent wells and shown as the fold change, according to the 2−ΔΔCq method (26). Primers were designed to amplify the circRNA-specific back splice junctions (Table I).

Table I.

Primer sequences used for reverse transcription-quantitative polymerase chain reaction validation.

circRNA ID	Sequence
mGAPDH-Q-F	TTCCTACCCCCAATGTGTCC
mGAPDH-Q-R	GGTCCTCAGTGTAGCCCAAG
mmu_circRNA_29981-Q-F	ACCAAAACCTGCATTGGCAC
mmu_circRNA_29981-Q-R	ACAGAACATGGCGATCTGGG
mmu_circRNA_004372-Q-F	GCTGAGCTACACAACTGCATC
mmu_circRNA_004372-Q-R	TGGTCCATGAATCGCTCTGG
mmu_circRNA_36096-Q-F	CCCTGCATTGGAAGTTTGCT
mmu_circRNA_36096-Q-R	TGAAGGTCCCACCATTTCCTG
mmu_circRNA_43021-Q-F	AGGAACTACCGTGGGGATGT
mmu_circRNA_43021-Q-R	CACCATGGGCCAAGATAGGT
mmu_circRNA_28028-Q-F	AGCGAGATCCCTACTGGTGT
mmu_circRNA_28028-Q-R	GTGATCTTCTCGTGCAACGC
mmu_circRNA_22646-Q-F	GAAGTGCTACGGGCAACAGA
mmu_circRNA_22646-Q-R	CATCGAAACACGGAAACGCC
mmu_circRNA_008234-Q-F	CGTGAGATGTGTTGCTGAGG
mmu_circRNA_008234-Q-R	CAGCAGGTGCTAAAGCATGAC

F, forward; R, reverse; circRNA, circular RNA.

Microarray data and bioinformatics analysis

Scanned images were imported into the Agilent Feature Extraction software (version 11.0.1.1) for raw data extraction. Quantile normalization of raw data and subsequent data processing were performed using the R software package (R version 3.1.2). The circRNAs that exhibited fold changes (FC) ≥1.5 (P≤0.05) between the two groups were categorized as significantly differentially expressed. The selected differentially expressed circRNAs were identified through volcano plotting and FC filtering. Hierarchical clustering analysis was performed to show the distinguishable circRNA expression pattern among samples.

The circRNA-miRNA interaction analysis was performed using Arraystar's miRNA target prediction software, which includes TargetScan (www.targetscan.org/) and miRanda (www.microrna.org/). Potential miRNA response elements (MREs) for circRNA were predicted using the miRNA support vector regression (mirSVR) which obtained from miRanda as previously described (27). Gene Ontology (GO) analysis was performed to explore the functional roles of host linear transcripts in the terms biological processes (BPs), cellular components (CCs) and molecular functions (MFs). Kyoto Encyclopedia of Genes and Genomes (KEGG) was used to determine the target genes in different biological pathways.

Statistical analysis

Data are presented as the mean ± standard deviation for triplicate measurements. Statistically significant differences between groups were estimated using the Student's t-test with SPSS software, version 13.0 (SPSS, Inc., Chicago, IL, USA). P<0.05 was considered to indicate a statistically significant difference.

Results

Successful construction of the NASH mouse model

At the end of 8 weeks of feeding the mice on the prescribed MCD diet, the body weight of the mice in the model group was significantly decreased compared with that in the control group (P<0.01; Table II). The H&E-stained sections exhibited typical histopathological features of NASH. As presented in Fig. 1, sections from the NASH livers exhibited disordered hepatic lobules, full fat vacuoles in lobule cells, infiltration of inflammatory cells, and cell swelling (Fig. 1A). The livers of MCD diet-fed mice exhibited marked lipid accumulation and liver fibrosis (Fig. 1B and C). Serum AST and ALT levels indirectly reflect the failure of liver function. As indicated in Table II, the serum ALT (406.7±132.1 U/l) and AST (573.3±169.5 U/l) activities were significantly increased following the administration of the MCD diet, as compared with the normal group (25.0±7.5 and 175.3±29.5, respectively; P<0.05). The TG and TC levels were lower in the MCD diet-fed mice compared with that in control mice (P<0.05 and P<0.01, respectively). The hepatic TG content (705±69.3 µmol/g protein) was higher in the MCD diet-fed mice compared with the control mice (179±17.6 µmol/g protein), whereas the hepatic TC content did not exhibit any difference between the two groups (Table II). These changes in serum biomarkers and liver tissue indicated the successful establishment of the NASH model in mice.

Table II.

Biochemical parameters in NASH mice.

Variables	Control	NASH
Initial body weight (g)	22.9±0.1	22.3±0.4
Final body weight (g)	30.1±2.3	13.8±0.1[b]
Serum alanine transferase (U/l)	25.0±7.5	406.7±132.1[a]
Serum aspartate transferase (U/l)	175.3±29.5	573.3±169.5[a]
Serum TG (mmol/l)	1.0±0.2	0.7±0.2[a]
Serum TC (mmol/l)	4.1±0.1	1.2±0.1[b]
Hepatic TG (µmol/g protein)	179±17.6	705±69.3[b]
Hepatic TC (µmol/g protein)	50.6±7.3	58.7±13.3

Control mice were MCS-fed; NASH mice were MCD-fed. Data are expressed as the mean ± standard deviation (n=3 in each group).

P<0.05 vs. control group

P<0.01 vs. control group. NASH, nonalcoholic steatohepatitis; TC, total cholesterol; TG, triglyceride.

Figure 1.

Liver morphology in NASH mice. Representative microscopic photographs of liver sections stained with H&E, Oil red O and picrosirius red staining. (A) H&E staining (magnification, ×200) and NAFLD activity scores in different groups. (B) Oil red O and (C) picrosirius red staining of liver sections from different groups. Bar graph shows quantification of positive areas. Data are expressed as the mean ± standard deviation (n=6 in each group). **P<0.01. NAFLD, non-alcoholic fatty liver disease; NASH, nonalcoholic steatohepatitis; H&E, hematoxylin and eosin.

Expression profiles of circRNA in NASH mice

To screen for circRNAs that were differentially expressed between the NASH and control mice livers, circRNA microarrays were performed in three samples from each group. A box plot showed that the distributions of circRNA for the tested samples were essentially identical following normalization (Fig. 2A). Hierarchical clustering, scatter plot visualization and volcano plot analyses revealed the differences in circRNA expression levels in NASH and control mice liver (Fig. 2B-D; FC >2.0; P<0.05). Consequently, 450 circRNAs were identified as dysregulated in the livers of NASH mice, with FC >1.5 and P<0.05. Among these, 298 circRNAs were upregulated, and 152 circRNAs were downregulated, in the NASH group compared with those in the control group (Fig. 3A). Additionally, 21 circRNAs were upregulated, and 12 circRNAs were downregulated, in the livers of NASH mice with FC >3.0 and P<0.05 (Fig. 3A). CircRNAs on mouse chromosomes are illustrated in the pie chart format in Fig. 3B. Types and counts of differentially regulated circRNAs are presented in Fig. 3C. The top 10 upregulated and downregulated circRNAs sorted by their FC values are summarized in Table III, along with their potential miRNA targets. Taken together, these data indicated that the circRNA expression patterns in NASH livers were different from those in the control group.

Figure 2.

Bioinformatics analysis of differentially expressed circRNAs in NASH mice. (A) Box plot, revealing the distributions of circRNAs in a more direct way. Following normalization, the distributions of log2 ratios among the samples were nearly the same. (B) Scatter plot showing the distributions of circRNAs, which revealed the differences in circRNA expression comparing between livers in the NASH and control groups. The values of the x- and y-axes in the scatter plot were the averaged normalized signal values of the group (on the log2 scale). (C) Heat maps revealed different circRNA expression profiles between livers in the NASH and control groups. The red color indicates a higher FC value; green indicates a smaller FC value. (D) Volcano plots comparing significantly dysregulated circRNAs in the NASH and control groups [the circRNA expression log2-transformed FC values (on the x-axis) were plotted against the-log10 P-values (on the y-axis)]. The red dots represent the circRNAs having an FC value >1.5 and P<0.05 comparing between the NASH and control groups. FC, fold change; circRNA, circular RNA; NASH, nonalcoholic steatohepatitis.

Figure 3.

Identification of differentially expressed circRNAs in NASH mice livers. (A) Counts of differentially expressed circRNAs between the NASH and control groups (FC ≥1.5 and FC ≥3, P<0.05). (B) The pie chart indicates where the circRNAs are located on the respective human chromosomes. (C) Types and counts of differentially regulated circRNAs (FC ≥1.5, P<0.05). circRNA, circular RNA; NASH, nonalcoholic steatohepatitis.

Table III.

Annotation for differentially expressed circRNAs/miR interaction.

circRNAs	Regulation	Fold-change	Type	MRE1	MRE2	MRE3	MRE4
mmu_circRNA_28028	Up	6.8	Exonic	mmu-miR-9769-3p	mmu-miR-8097	mmu-miR-194-2-3p	mmu-miR-693-3p
mmu_circRNA_29981	Up	6.1	Exonic	mmu-miR-7661-3p	mmu-miR-181b-5p	mmu-miR-6952-3p	mmu-miR-3077-3p
mmu_circRNA_39832	Up	5.9	Exonic	mmu-miR-6344	mmu-miR-130a-5p	mmu-miR-7091-3p	mmu-miR-7214-5p
mmu_circRNA_29980	Up	4.4	Exonic	mmu-miR-1249-5p	mmu-miR-6923-3p	mmu-miR-7058-3p	mmu-miR-6976-3p
mmu_circRNA_21454	Up	4.2	Exonic	mmu-miR-1903	mmu-miR-7028-3p	mmu-miR-7032-3p	mmu-miR-3059-5p
mmu_circRNA_34829	Up	4.0	Exonic	mmu-miR-6943-3p	mmu-miR-7064-5p	mmu-miR-691	mmu-miR-7032-5p
mmu_circRNA_40141	Up	3.9	Exonic	mmu-miR-1903	mmu-miR-223-3p	mmu-miR-881-3p	mmu-miR-7646-3p
mmu_circRNA_30398	Up	3.8	Exonic	mmu-miR-7061-3p	mmu-miR-103-1-5p	mmu-miR-103-2-5p	mmu-miR-136-5p
mmu_circRNA_004372	Up	3.7	Exonic	mmu-miR-103-1-5p	mmu-miR-103-2-5p	mmu-miR-107-5p	mmu-miR-143-5p
mmu_circRNA_21453	Up	3.7	Exonic	mmu-miR-1903	mmu-miR-7028-3p	mmu-miR-6946-3p	mmu-miR-6919-3p
mmu_circRNA_28143	Down	6.0	Intronic	mmu-miR-7118-5p	mmu-miR-1894-3p	mmu-miR-7665-5p	mmu-miR-1249-5p
mmu_circRNA_28144	Down	5.8	Intronic	mmu-miR-7118-5p	mmu-miR-1894-3p	mmu-miR-7665-5p	mmu-miR-1249-5p
mmu_circRNA_22646	Down	4.7	Exonic	mmu-miR-196a-5p	mmu-miR-6906-5p	mmu-miR-6909-3p	mmu-miR-196b-5p
mmu_circRNA_31317	Down	4.3	Exonic	mmu-miR-3099-5p	mmu-miR-6976-5p	mmu-miR-5627-3p	mmu-miR-497b
mmu_circRNA_008234	Down	4.2	Antisense	mmu-miR-361-3p	mmu-miR-504-5p	mmu-miR-1291	mmu-miR-3090-3p
mmu_circRNA_21487	Down	3.8	Exonic	mmu-miR-5621-5p	mmu-miR-7665-5p	mmu-miR-3069-5p	mmu-miR-6907-5p
mmu_circRNA_35127	Down	3.4	Exonic	mmu-miR-3618-5p	mmu-miR-6973b-5p	mmu-miR-295-5p	mmu-miR-677-5p
mmu_circRNA_40509	Down	3.3	Exonic	mmu-miR-6932-3p	mmu-miR-26a-2-3p	mmu-miR-7074-5p	mmu-miR-3075-5p
mmu_circRNA_29948	Down	3.1	Intronic	mmu-miR-6946-3p	mmu-miR-7651-3p	mmu-miR-6399	mmu-miR-7116-3p
mmu_circRNA_39830	Down	3.1	Exonic	mmu-miR-6908-3p	mmu-miR-6975-3p	mmu-miR-149-5p	mmu-miR-500-5p

CircRNA, circular RNA; miR, microRNA; MRE, miR response elements.

Annotations for the circRNA-miRNA interactions

To further validate the results of the present study, Arraystar's miRNA target prediction software, based on TargetScan and miRanda, was used to predict the circRNA-miRNA interactions. A total of 20 differentially expressed exonic circRNAs with the highest FC values were selected (FC >3; P<0.05) to predict their miRNA response elements (MREs), including 10 upregulated circRNAs and 10 downregulated circRNAs. Four MREs that were identified to have ‘good’ mirSVR scores for each circRNA are presented in Table III.

Validation of selected circRNAs by RT-qPCR

To validate the microarray results, seven differentially expressed exonic circRNAs were randomly selected (FC >3.5; P<0.05), including five upregulated circRNAs and two downregulated circRNAs. Primers for circRNAs were designed (shown in Table I), and the appearance of a single peak in the melting curve confirmed the specificity of the amplified products for each circRNA (Fig. 4). The results demonstrated that five circRNAs (circRNA-28028, circRNA-29981, circRNA-004372, circRNA-0043021 and circRNA-36096) were upregulated, and these findings were broadly in accord with the data obtained from the microarray analysis. However, only circRNA-29981 reached the required level of statistical significance to become a candidate for further study (Fig. 5). By contrast, circRNA-008234 was observed to have the opposite expression trend in microarray and RT-qPCR analysis, although the difference was revealed not to be statistically significant. Additionally, the expression levels of circRNA-22646 were not significantly different between the NASH and control groups (Fig. 5). Furthermore, the detailed annotations for circRNA-29981 interaction with various miRNAs (miR-7661-3p, miR-181b-5p, miR-6952-3p, miR-3077-3p, and miR-1912-3p) are presented in Fig. 6. The effect of the circRNA-29981 regulation pathways in NASH is worth further exploration.

Figure 4.

Melting curves derived from the reverse transcription-quantitative polymerase chain reaction experiments for each circRNA. circRNA, circular RNA.

Figure 5.

Validation of selected circRNAs by RT-qPCR. Each RT-qPCR assay was performed at least three times. *P<0.05 in RT-qPCR verification when comparing data between the NASH and control groups. RT-qPCR, reverse transcription-quantitative polymerase chain reaction; circRNA, circular RNA.

Figure 6.

mmu_circRNA-29981 and its top five predicted miRNA targets. The analysis included MRE sequences, miR seed types, AU-richness surrounding the seed sequences, and relative positions of MREs in the linearized sequences of circRNAs. A/U (marked with red bars) near the MRE facilitates accessibility, whereas G/C (marked with black bars) means low accessibility. The heights of the bars represent the extent of accessibility. Binding sites of conserved miRs with good mirSVR scores are represented. miR, microRNA; circRNA, circular RNA; MRE, miR response element.

Bioinformatics analysis predicting biological function of host linear transcripts

CircRNAs are generally formed by alternative splicing of pre-mRNA species, which are usually spliced to generate a linear mRNA (11). A more recent study has demonstrated that a class of circRNAs, termed exon-intron circRNAs, enhances the expression of their pre-mRNAs (28). CircRNAs may also be potentially in competition with their linear transcripts (24). To provide a comprehensive landscape of the origination of circRNAs, GO and KEGG analyses were performed to annotate the biological functions of host linear transcripts. The top 10 generally altered GO terms between NASH and control groups were classified according to the BP, CC and MF (Fig. 7). For increased mRNAs, the most enriched and meaningful BP terms were ‘cellular process’, ‘single-organism process’ and ‘single-organism cellular process’. The most enriched CC terms were ‘cell’, ‘cell part’, ‘intracellular’, and ‘intracellular part’. As for MF, the most enriched terms were ‘binding’, ‘protein binding’, and ‘organic cyclic compound binding’. In addition, the data from the GO analysis on the downregulated transcripts revealed that the top two GO processes were identical with those in the increased transcripts (Fig. 7). Furthermore, KEGG analysis revealed that certain pathways, including regulation of the actin cytoskeleton, 2-oxocarboxylic acid metabolism, and the endocytosis signaling pathway, were associated with the upregulated transcripts, whereas pathways such as the pyruvate metabolism, fatty acid biosynthesis and fatty acid metabolism signaling pathways were associated with the downregulated transcripts (Fig. 8).

Figure 7.

GO analysis of increased (left panel) and decreased (right panel) host linear transcripts. Prediction terms with P<0.05 were selected. The top ten generally changed GO terms were classified by BP, CC, and MF. GO, Gene Ontology; BP, biological process; CC, cellular component; MF, molecular function.

Figure 8.

Kyoto Encyclopedia of Genes and Genomes pathways of up (left panel) and down (right panel) regulated host linear transcripts. The -log10 (P-value) yields an enrichment score indicating the significance of pathway correlations. DE, differentially expressed; ABC, ATP-binding cassette; TCA, tricarboxylic acid.

Discussion

The mechanism of NASH has yet to be fully elucidated. Because of a lack of effective diagnosis and therapeutic approaches, NASH has become the predominant cause of cryptogenic cirrhosis. An in-depth analysis of the pathogenic mechanism of NASH is essential for proper diagnosis and the implementation of effective therapies. CircRNAs have been associated with numerous diseases, cell processes and gene expression studies (10). The characteristics of circRNAs, such as tissue- and development-specific expression, serving as a ‘miRNA sponge’ to inhibit the activity of miRNAs, highlight the role of circRNAs in disease development. Accumulating evidence has suggested that circRNAs are associated with liver diseases, including HCC, liver regeneration and hepatic steatosis (18–21,29). Certain circRNAs have been revealed to have important roles in regulating HCC growth, migration, and invasion (19,22). Recently, Chen et al (30) demonstrated that inhibiting the expression of has_circ_0071410 affected irradiation-induced hepatic stellate cell (HSC) activation by regulating the expression of miR-9-5p. Guo et al (29) demonstrated that the dysregulation of circRNAs is associated with hepatic steatosis. These studies suggested that circRNAs are closely associated with liver diseases, and further in-depth analyses of circRNAs may provide novel clues to obtain a deeper understanding of the pathogenesis and processes of NASH.

Since it is difficult to obtain liver biopsy samples from patients with NASH, either a NASH mouse model or a hepatic steatosis cell model has been used to investigate the circRNA profile (24,29). These studies demonstrated that dysregulated circRNAs were associated with hepatic steatosis and NASH development. However, currently there is no ideal NASH animal model (31); therefore, a classic NASH mouse model that is induced by an MCD diet was used in the current study (32). To confirm the success of the mouse model, H&E, Oil red O and picrosirius red staining were performed. Additionally, the serum and liver biochemistry also demonstrated the pathological features of NASH. The aforementioned staining and biochemistry analysis confirmed a successful simulation of the pathological characteristics of NASH. In the present study, 450 significantly dysregulated circRNAs (FC ≥1.5, P<0.05) were identified, including 298 upregulated and 152 downregulated circRNAs in the NASH group. Furthermore, the most significantly differentially expressed circRNAs comparing between the NASH and control groups were annotated in detail, and the circRNA-miRNA interaction was predicted. Following validation by RT-qPCR, a novel circRNA, circRNA_29981, was confirmed to be differentially expressed in NASH model mice, and the circRNA-miRNA interaction information was predicted. With respect to the target miRNAs of circRNA_29981, miR-181b has been validated as a key regulator in steatosis via the targeting of sirtuin 1 (SIRT1) (33). SIRT1 has been reported to have an essential role in guiding the transition of the HSC phenotype (34). miR-181b was also reported to promote HSC proliferation, and has an important role in liver fibrosis (23,35,36). Another miRNA, mmu-miR-130a, has been demonstrated to enhance the activation of HSCs by suppressing the expression of peroxisome proliferator-activated receptor γ (37). Therefore, one may surmise that circRNA_29981 could be involved in regulating the activation of HSCs. Further research on the role of circRNA_29981 is required, including experiments designed to validate the target miRNAs of circRNA_29981, and investigate the molecular changes of a hepatic cell line following circRNA knockdown or overexpression. Further exploration of the circRNA-miRNA interaction networks will be the subject of future investigations.

Jin et al (24) previously used BALB/c mice to create a NASH model by providing an MCD diet for 2 weeks. In the present study, C57BL/6 mice received an MCD diet for 8 weeks. No common circRNAs were identified in the top 10 up- or downregulated circRNAs when comparing between the two studies. This may be due to the difference in the time for which the MCD diet was supplied, or that different mouse strains were used. However, following the GO pathway analysis, the results of the most enriched BP, CC, and MF terms from the host linear transcripts in the present study were identical with those of the study by Jin et al (24). This may potentially suggest that there are certain associations between circRNAs and linear mRNAs in NASH processes. The data obtained from the KEGG pathway analysis revealed that host linear transcripts were associated with pathways of energy metabolism, including ‘amino acid biosynthesis’ and ‘fatty acid biosynthesis and metabolism’, which have been credited as being the important activities in NASH processes. In addition, glycerolipid metabolism and inflammatory mediator regulation of transient receptor potential channels were revealed to be important in terms of NASH progression in the study by Jin et al (24). Taken together, all these findings reflect the complicated and critical mechanism underlying circRNA-associated regulation.

Although circRNA_29981 in the mouse shares 89% sequence homology with the sequence of circRNA_29981 in humans (Fig. 9), the expression profile of its homologous sequences in humans has yet to be determined. Therefore, whether this result can be verified in human NASH requires further confirmation. In addition, due to the low validation rate of circRNA microarray data compared with the RT-qPCR experiments, it is not possible to exclude that there may be other dysregulated circRNAs that also participate in NASH development. Further studies on these issues are required. As circRNAs are stably expressed in blood, identifying the circulating circRNAs in patients with NASH may provide a promising avenue for research. This will also be an integral part of future studies. Taken together, the present study has produced circRNA expression profiles during the NASH process in mice, and a novel candidate circRNA, circRNA_29981, that may be a potential regulator of HSC activation was identified. The results of the present study have corroborated the regulatory role of circRNAs in NASH, and provided novel clues for understanding the mechanism of NASH pathogenesis.

Figure 9.

The results of the homologous sequence analysis of circRNA_29981 comparing between mouse (upper) and human (lower) are shown. Note that circRNA_29981 has two homologous sequences in humans. The mouse circRNA_29981 was shown to have 89% sequence homology compared with its sequence in humans. CircRNA, circular RNA.