Alterations in cirRNA Expression Profile during the Development of Acute-on-Chronic Liver Failure Induced by Hepatitis B

Background . Hundreds of millions of people worldwide su ﬀ er from chronic hepatitis B (CHB), making it one of the most serious public health problems. CHB can lead to serious acute-on-chronic liver failure (ACLF). Since ACLF has a high mortality rate, predicting the risk of ACLF is a critical issue in clinical practice. Methods . To investigate the occurrence of ACLF in CHB, we used high-throughput RNA sequencing to detect the expression of a novel endogenous noncoding RNA (circRNA) in healthy controls (HC) as well as CHB and ACLF patients. Di ﬀ erentially expressed circRNAs were selected and analyzed by Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) biological pathway analyses, circRNA-miRNA coexpression networks, and statistical analyses. Results . A total of 21,101 circRNAs were identi ﬁ ed in HC, CHB, and ACLF subjects. Compared with the HCs, 4 circRNAs were upregulated and 14 circRNAs were downregulated in ACLF subjects, with a consistent trend in all three groups. GO analysis revealed that regulation of lymphocyte activation and macrophage tolerance induction were the most important biological processes in ACLF progression. KEGG pathway analysis revealed that primary immunode ﬁ ciency and NOD-like receptor signaling pathways were the most enriched terms. The circRNA-miRNA coexpression network and statistical analyses showed that circRNA_07734, circRNA_08533, and circRNA_16083 were the most important circRNAs in the process of ACLF. Conclusions . Immune dysfunction and de ﬁ ciency may play a key role in the development of ACLF, especially circRNA_07734, circRNA_08533, and circRNA_16083.


Introduction
Chronic hepatitis B (CHB) is one of the most serious public health problems worldwide. More than 2 billion people were infected with hepatitis B virus (HBV), 360 million of whom were chronically infected [1]. Although CHB is initially asymptomatic, it can lead to serious acute-on-chronic liver failure (ACLF, HBV-ACLF) after several years, which has a poor prognosis and a high mortality rate (>70%) [2]. ACLF is serious in CHB patients, who constitute approximately 70% of all ACLF cases in most Eastern countries [2,3]. In China, due to the high prevalence of chronic HBV infection, HBV-associated ACLF patients account for approximately 80% of all ACLF cases [4]. Since ACLF has a very high mortality rate, prediction of ACLF risk is a major issue in the medical field.
At present, the model for end-stage liver disease (MELD) is recognized as a reliable method to measure the mortality risk of patients with liver failure. Shen et al. assessed six widely used short-and long-term prognostic models for patients with ACLF. They concluded that the integrated MELD model may be optimal [5]. Liu et al. suggested that incorporating the γ-glutamyl transpeptidase to platelet ratio into the MELD model may allow more accurate prediction of survival in ACLF patients [6]. Ma et al. established a scoring system offering superior predictive performance, in terms of both specificity and sensitivity, for ACLF patients compared with the MELD [7]. Wei et al. reported that serum Golgi protein 73 (GP73) may be useful in the diagnosis of ACLF in populations with chronic HBV infections [8]. Cui et al. evaluated the prognostic utility of serum Hcy levels, which may serve as a biomarker for 3-month mortality rate in ACLF patients [9]. Although many studies predicting the risk of ACLF have been performed, no method had been proven effective and the pathogenesis of ACLF remains unclear.
Circular RNA (circRNA) is a novel endogenous noncoding RNA characterized as a widespread, abundant, stable, conserved, and tissue-specific molecule in mammalian cells [10]. Different from linear RNA, circRNA is not degraded by RNase R; it is composed of a class of RNA developing covalently closed loop structures without 5 ′ -3 ′ polarities [10,11]. circRNAs are important competing endogenous RNAs. They act as sponges for microRNA (miRNA) and regulate the expression of miRNA-targeted transcripts. Moreover, circRNAs have the potential to serve as biomarkers for the diagnosis and prognosis of several cancers [12,13]. For example, circRNAs associated with the occurrence and development of hepatocellular carcinoma (HCC) [14], such as circRNA_100338 and circRNA_0005075, serve as biomarkers for HCC and as targets for HCC therapeutics associated with HBV [15][16][17].
The numbers of circular RNAs are different between CHB and ACLF and increase with symptom severity [18]. So far, there is no report on circRNAs in CHB and ACLF, and whether aberrant expression of circRNAs plays a role in ACLF is still unknown. Our study is aimed at investigating the numbers of circRNAs in CHB and ACLF and at identifying the important circRNAs in ACLF development. Blood samples (5 mL) were collected into heparinized tubes and diluted 2 : 1 with phosphate-buffered saline. Mononuclear cells (MNCs) were obtained with equal amounts (1.077 g/cm 3 ) of Ficoll-Paque (Ficoll-Paque™ PLUS; GE Healthcare, Chicago, IL, USA) and density gradient centrifugation. Samples with an RNA integrity number ðRINÞ ≥ 7 were subsequently subjected to analysis.

circRNA Library Construction and Identification of
Differentially Expressed circRNAs. For the circRNA library construction, we used the circRNA analysis software CIRI. First, raw reads were filtered to remove connector sequences, sequences containing more N, and low-quality reads, and clean reads were obtained. Then, by comparing with the ribosomal database, the ribosomal RNA sequences were removed and effective reads were obtained. On the one hand, effective reads were compared with the reference genome by HISAT2 (mapping) and the mapping reads were assembled into mRNA by string to generate the mRNA expression matrix. Effective reads, on the other hand, predicts the back-splicing junction of circRNA through CIRI, The predicted circRNA was obtained, and then, the predicted circRNA was rematched back to the reads of the reverse shearing site to obtain the expression matrix of cir-cRNA. Meanwhile, the circRNA shear ratio could be obtained by the combined analysis of the two aspects.
We use circSeq data to compare and analyze the differential expression of the same circRNA in two samples, and two criteria were selected: First is fold change, that is, the change multiple of the expression level of the same circRNA in two samples. The second is P value or FDR (adjusted P value). The P value of each circRNA is calculated by FDR error control method, and then, the P value is corrected by multiple hypothesis testing. The screening condition of differentially expressed circRNAs is P < 0:05 and fold change > 2.
2.4. Identification of circRNAs and mRNA. The concentration and purity of total RNA extracted from MNCs were measured on a NanoDrop ND 1000 spectrophotometer (Thermo Fisher Scientific, Inc.) at 260 and 280 nm wavelengths. Before the synthesis of firstand second-strand cDNA, Ribo-Zero (ribosomal RNA) was depleted and RNA-fragmented. cDNA was purified followed by processing of adenylate 3 ′ ends, ligation of adapters, and enrichment of DNA fragments. Finally, 1 μL was loaded onto a 2100 Bioanalyzer (Agilent Technologies) for validation. High-throughput RNA sequencing was performed according to the method of Memczak et al. [19].
Raw data were extracted using Agilent Feature Extraction software and predicted by circBase. circRNAs and mRNA were quantitatively analyzed by Shanghai OE Biotech Ltd., Co. (Shanghai, China). The circRNAs and mRNA of the HC, CHB, and ACLF subjects were verified with CIRI software (Shanghai OE Biotech), and the identified sequence was selected for further analysis. The expression levels of cir-cRNAs among the three groups were measured by mapped back-splice junction reads per million mapped reads (RPM), and the mRNA_expression was measured by FPKM (Fragments Per Kilobase per Million). The location of the chromosome where the circRNA sequence overlapped was annotated. The criteria were chosen according to the fold change and P value of the circRNA and mRNA in the same sample between two groups.

Pathway Analyses and miRNA Target Prediction.
Gene functions were classified into three subgroups: biological process (BP), cellular component (CC), and molecular function (MF). The most enriched GO terms among the CHB and ACLF groups (ranked by enrichment score) were selected, and the enriched GO terms between HC and ACLF were ranked by the number of differentially expressed linear transcripts. The different biological pathways of linear transcripts were analyzed by the KEGG pathway.
circRNA can be used as a target molecule of miRNA and regulated by miRNA. Analysis of the circRNA-miRNA interaction could help elucidate the function and mechanism of circRNAs in contact with miRNAs. The circRNA-miRNA interactions were predicted by Arraystar. The target interaction network diagram of the circRNA-miRNA/gene interaction network was built by Cytoscape (available on http:// www.cytoscape.org/). Spearman's correlation was performed to evaluate the relation between circRNA with mRNA and liver function indexes. All statistical analyses were conducted with SPSS software (ver. 22.0; IBM Corp., Armonk, NY, USA), and P < 0:05 was considered statistically significant.

Overview of circRNA Profiles.
To investigate the expression profile of circRNAs involved in HC, CHB, and ACLF, tissues of these three groups were analyzed by circular RNA high-throughput sequencing. Results showed that approximately 86.61% of circRNAs had a predicted spliced length of less than 2,000 nt; 48.44% had a length less than 500 nt, and 26.53% had a length between 500 nt and 1,000 nt (Figure 1(a)). Moreover, Venn diagram of circRNA expression showed that a total of 21,101 circRNAs in the HC, CHB, and ACLF were detected. Among the 21,101 cir-cRNAs, 3581 circRNAs were differentially expressed in CHB tissues compared with HC and ACLF; 3214 circRNAs were differentially expressed in ACLF tissues compared with HC and ACLF (Figure 1(b)). There are no outliers of circRNA expression in the HC, CHB, or ACLF group due to instrument or operational anomalies. The distribution of cir-cRNAs on the genome is shown in Supplemental Figure 1.

Identification of Differentially Expressed circRNAs and
mRNA. According to the circRNA criteria (fold change > 2 and P < 0:05), 550 differentially expressed circRNAs in the CHB and ACLF groups were identified compared with HC. Among them, 176 differentially expressed circRNAs were detected in the CHB group (89 upregulated and 87 downregulated) and 374 differentially expressed circRNAs were detected in the ACLF group (181 upregulated and 193 downregulated) compared with HC. The differentially expressed circRNAs between ACLF and HC groups are shown in a volcano plot (Figure 2(a)). The unsupervised hierarchical cluster analysis results are shown in Figure 2(b).
The differentially expressed candidate circRNAs among the three groups showing progressive up-or downregulation were further analyzed. Compared with HC, there were 4 upregulated and 14 downregulated circRNAs (Table 1, Figure 3) in ACLF. mRNAs were analyzed by the same method. The results showed that there were 7 upregulated and 9 downregulated mRNA in ACLF compared with HC. The expression information and statistic differences of cir-cRNA and mRNA are shown in the Supplemental Table 2, 3. 3.3. GO and KEGG Analyses. The functions of the target genes of differentially expressed circRNAs in the ACLF group were further analyzed by GO enrichment and KEGG pathway analysis. Significant terms in the GO and KEGG pathway analyses (P < 0:05) were selected and ranked by P value.
After obtaining the differentially expressed circRNA, we conducted GO enrichment analysis for the differentially expressed circRNA using the information of circRNA source genes. The number of differential transcripts included in each GO entry was counted, and the significance of differential circRNA enrichment in each GO entry was calculated using hypergeometric distribution test. The results of calculation will return a P value of enrichment significance, small P value indicates enrichment of circRNA in this GO item, and its calculation formula is as follows: Enrichment score is calculated as follows: N is the number of circRNAs with GO annotation in all circRNAs; N is the number of circRNAs with GO annotation in source genes of differentially expressed circRNAs in N; M is the number of circRNAs whose source genes are annotated as a specific GO term among all circRNAs; M is the number of differentially expressed circRNAs annotated by the source gene for a specific GO term. The biological The top-10 downregulated GO terms, classified by BP, CC, and MF, are shown in Supplemental Figure 2A. The most enriched BP, CC, and MF terms were regulation of lymphocyte activation, integrin alphaL-beta2 complex, and RNA polymerase I core binding, respectively. Among the upregulated GO terms, the most enriched BP, CC, and MF terms were positive regulation of macrophage tolerance induction, Golgi lumen, and extracellular matrix structural constituent, respectively (Supplemental Figure 2B).
Identified circRNAs were annotated using the KEGG public database. Pathway analysis was carried out on differentially circRNA-derived genes (combined with KEGG annotation results) using the KEGG database, and the significance of differentially circRNA enrichment in each pathway item was calculated using the hypergeometric     Figure 3) revealed that the most enriched term in downregulated circRNAs was primary immunodeficiency, while in upregulated circRNAs, NOD-like receptor signaling pathway was the most enriched term. As shown in Table 2, the KEGG pathway hsa05340 was the most important.

Correlation
Analysis. Spearman's correlation analysis indicated that the RPM values ofcircRNA_07734, cir-cRNA_08533, and circRNA_16083 were related to ALT, AST, TBIL, and IBIL (Table 3). circRNA_07734 and cir-cRNA_08533 showed strong negative correlation with ALT, AST, TBIL, and IBIL, which indicated if ALT, AST, TBIL, and IBIL increased higher, the expression of cir-cRNA_07734 and circRNA_08533 would be lower. cir-cRNA_16083 is different with circRNA_07734 and circRNA_08533, which had positive correlation with ALT, AST, TBIL, and IBIL. In clinical practice, ALT, AST, TBIL, and IBIL are used to evaluate the liver function. Commonly, the higher the ALT, AST, TBIL, and IBIL levels were, the more severe the liver injury will be. Therefore, circRNA_ 07734, circRNA_08533, and circRNA_16083 may be useful as the sensitive indicators of liver injury.

Discussion
circRNAs have received increasing attention due to their highly stable and specific spatiotemporal expression patterns. circRNAs are a transcriptional product in various tissue and cell types in humans, mice, drosophila, and nematodes, among other species. circRNAs can be classified into three types: exonic (circRNA arising from the exons of the linear transcript), intergenic (circRNA located outside of a known gene locus), and sense-overlapping (circRNA whose gene locus overlaps with linear RNA but is transcribed from the opposite strand). The differentially expressed hosting genes of the circRNAs and the disordered splicing machinery will lead to circRNAs differentially expressed. So far, numerous circRNAs have been reported in different animals. The majority of circRNAs showing specific expression patterns are related to tissue development and disease [19]. Therefore, circRNA has potential as a biomarker for some diseases [11]. Thus far, many studies have reported abnormal expression levels of circRNAs in a number of cancers, such as gastric [20,21], colorectal [22,23], breast [24,25], and lung cancers [26,27]. However, few studies have focused on the development of ACLF.
In this study, we detected several circRNAs in HC, CHB, and ACLF subjects. The circRNAs were initially identified by reference to a database or website. circBase, (http://circbase .org/) which includes merged genomic data on circRNAs from humans, mice, nematode worms, pike, and coelacanth, provides evidence of their expression [28]. Therefore, we identified the circRNAs detected herein using circBase. circRNA-seq showed that there were 4 upregulated and 14 downregulated circRNAs in ACLF compared with HC. Dysregulation of circRNAs in patients with chronic hepatitis B prompted us to investigate the functional role of these circRNAs.
CHB is mainly caused by aberrant cellular immunity, due, for example, to natural killer (NK) cells, cytotoxic T cells, and antibody-dependent lymphocytes. When infected with hepatitis B, three types of lymphocytes target antigens of the hepatocyte membrane: HBsAg, HBcAg, and other immunoregulatory cells. Auxiliary and inhibitory T cells take part in the immune response, which can lead to  [29] reported that tetramer + CD 4+ T cell numbers were reduced during CHB. Zhang et al.
reported that HBV-specific T cell responses to recombinant HBV core protein were reduced in CHB patients and positively correlated with HBV viral load in coinfected, chronic HBV patients [30]. Therefore, the progression of ACLF from CHB may be related to the BP of immunodeficiency or dysfunction, which can in turn be caused by HTLV-I infection, Staphylococcus aureus infection, or other diseases. The process of immunity dysfunction may be related to the NODlike receptor signaling pathway, the NF-kappa B signaling pathway, and/or the interleukin-1-mediated signaling pathway.
Since circRNAs contain multiple miRNA-binding sites, the function of individual circRNAs can be predicted by their miRNA target gene and annotated according to the function of their miRNA target gene. The circRNA-miRNA coexpres-sion network identified two upregulated circRNAs and nine downregulated circRNAs in this study. After combining the results from the GO and KEGG analyses, circRNA_00389, circRNA_07734, circRNA_08533, circRNA_ 15699, and cir-cRNA_16083 were identified as the most important cir-cRNAs during the progression of ACLF. Among them, circRNA_00389 had the highest enrichment score in the KEGG analysis, consistent with the GO analysis. The top term in BP analysis was regulation of lymphocyte activation, which corresponded to downregulated circRNA_00389. However, the expression of circRNA_00389 had no statistical significance between the ACLF and HC group. Therefore, combined with the results of Spearman's correlation analysis, circRNA_07734, circRNA_08533, and circRNA_16083, which are related to immune dysfunction and correlate with ALT, AST, TBIL, and IBIL, are important circRNAs in ACLF.
This study was conducted to initially explore the role of circRNAs in CHB and acute liver failure due to CHB. The etiology of ALCF is CHB. Compared with healthy subjects, many circRNAs changed in CHB and ACLF. We speculated that circRNAs was an intermediate regulatory molecule that alters the activity of downstream miRNA and thus regulates mRNA levels of some genes. Therefore, we did not find significant changes in circRNAs between the CHB and ACLF groups.
The limitation of this study is that there were only three samples in each group, but the selected circRNAs all had the same expression tendency in CHB and ACLF (Supplemental Table 1). The results of U test analysis between groups were consistent with the expression tendency of selected circRNAs. Besides, the Spearman correlation analysis of circRNA_07734, circRNA_08533, and circRNA_16083 showed that there was statistical significance with differentially expressed mRNA and ALT, AST, TBIL, and IBIL (Table 3). Therefore, our results were credible and can be used as the basis for further studies on the regulatory role of circRNA downstream target genes in CHB and ACLF. In further work, we should perform qRT-PCR to determine the expression levels of these circRNAs and carry out research on mechanisms.

Conclusion
In this work, 21,101 circRNAs were detected in CHB, HC, and ACLF groups. In total, 4 upregulated and 14 downregu-lated circRNAs were identified in the ACLF group compared with the HC group, which showed a consistent trend of change in the three groups. GO and KEGG biological pathway analyses and circRNA-miRNA coexpression network and statistical analyses revealed that circRNA_07734, cir-cRNA_08533, and circRNA_16083 are the most important circRNAs in ACLF. Their RPM values correlated with ALT, AST, TBIL, and IBIL; thus, immune dysfunction or immunodeficiency may be critical for ACLF development.

Data Availability
The datasets generated and/or analyzed during the current study are not publicly available but are available from the corresponding author on reasonable request.

Disclosure
An earlier version, named as "Circular RNA expression profiles in the development of acute-on-chronic liver failure in hepatitis B virus-infected patients," has been presented as preprint in http://www.researchsquare.com, in the following link: https://www.researchsquare.com/article/rs-38521/v3. The funder had no role in the study design, data collection, analysis, decision to publish, or manuscript preparation. Table 3: Spearman's correlation analysis of circRNA_07734, circRNA_08533, and circRNA_16083 with differentially expressed mRNA and alanine transaminase (ALT), aspartate aminotransferase (AST), total bilirubin (TBIL), and indirect bilirubin (IBIL) (n = 9).