2.2. Differentially-Expressed miRNAs and Target Genes’ Enrichment Analysis

In this work, 112 differentially-expressed (DE) miRNAs were identified between before and after FZHY-treated HBC patients, including 63 upregulated miRNAs and 49 downregulated miRNAs (Figure 2A,B). Using the prediction programs, the miRNA target genes were predicted. After redundancy analysis, the final database of 112 DE miRNAs with 5161 target genes was built.

Open in a separate window

Figure 2

The bioinformatics analysis of miRNAs and target genes between before and after FZHY treatment for hepatitis B-caused cirrhosis (HBC) patients. (A) Heat-map of the differential expression of miRNAs between before and after FZHY treatment. Red color represents that the miRNA is upregulated, and green color represents downregulated. The relationship among the samples was divided by binary tree classification and is shown in the upper portion. The hierarchical cluster of miRNAs is displayed at the bottom; (B) The distribution of differential expressed (DE) miRNA expression levels; the right side represents the up-expressed miRNAs between before and after FZHY treatment; the left side represents the down-expressed miRNAs between before and after FZHY treatment; (C–E) MiRNA target genes related Gene ontology (GO) terms, Kyoto Encyclopedia of Genes and Genomes (KEGG) terms and disease terms.

To understand these miRNAs holistically, we conducted functional enrichment analysis for the target genes using DAVID (The Database for Annotation, Visualization and Integrated Discovery) online analysis [20,21]. GO (Gene ontology) term (10% top terms) analysis reveals that the miRNA-targets are mainly associated with metal ion binding, transition metal ion binding, regulation of transcription, regulation of the RNA metabolic process, DNA binding, transcription regulator activity, regulation of transcription from the RNA polymerase II promoter, positive regulation of the RNA metabolic process, positive regulation of transcription, phosphate and the phosphorus metabolic process (Figure 2C). KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways show an impressive functional association with various signal-related pathways, such as the neurotrophin signaling pathway (p = 5.45 × 10⁻¹¹), the MAPK signaling pathway (p = 3.02 × 10⁻⁹), the ErbB signaling pathway (p = 1.37 × 10⁻⁸), the insulin signaling pathway (p = 1.44 × 10⁻⁷), the mTOR signaling pathway (p = 5.98 × 10⁻⁶), the TGF-β signaling pathway (p = 9.16 × 10⁻⁶), the T cell receptor signaling pathway (p = 9.77 × 10⁻⁶), the Wnt signaling pathway (p = 5.26 × 10⁻⁵), as well as the B cell receptor signaling pathway (p = 5.34 × 10⁻⁴) and the Fc epsilon RI signaling pathway (p = 9.20 × 10⁻⁴) (Figure 2D). Furthermore, axon guidance, ubiquitin-mediated proteolysis, renal cell carcinoma, acute myeloid leukemia, non-small cell lung cancer, focal adhesion, long-term potentiation, neuroactive ligand-receptor interaction and endocytosis also have a high correlation with these target genes.

Furthermore, the disease terms show that these target genes are mainly associated with developmental, neurological and psychological diseases, such as schizophrenia, attention-deficit hyperactivity disorder (ADHD), alcoholism, autism, osteoporosis, brain aneurysm, spinocerebellar disease, normal variation and bipolar disorder (Figure 2E). Although none of the disease terms appear relevant for hepatitis B or fibrosis, however, these terms in fact are related to the complications and often occur with liver cirrhosis. This result suggests that FZHY can obviously improve the quality of life of patients.

2.3. miRNA-Target Network Building and Potential miRNA Marker Screening

Using the DE miRNAs and target genes, the miRNA-target network was constructed. As shown in Figure 3A, the topological profile of the network is more likely similar to the “Medusa” model [22], which consists of a regulatory core structure by kernel nodes and represents most prominently the functions in the network. This suggests that the kernel nodes are determinants in the realized network profiles, while the periphery nodes only should be regulated [23]. For this reason, 12 clusters with 13 kernel miRNAs were isolated from the miRNA-target network based on the ClusterONE algorithm [24] (Figure 3B, Table 1). These kernel miRNAs might play important roles in HBC treatment progression and also might act as potential biomarkers for the efficacy evaluation of FZHY.

Open in a separate window

Figure 3

The miRNA-target network was constructed based on DE miRNAs and target genes, and the related clusters were identified using the ClusterONE algorithm. (A) The global profiles of the FZHY treatment-related miRNA-target network were built. The target genes of miRNAs were predicted using the TarBase (v7.0), miRecords and miRTarBase databases and the miRanda, miRDB, miRWalk and RNAhybrid prediction programs; (B) The kernel miRNA-related clusters isolated from the network using the ClusterONE algorithm, which was defined as node > 5, density > 0.05, quality > 0.60 and p < 0.0001. There are 12 clusters, containing 13 kernel miRNAs.

2.4. Establishing and Validating the Potential miRNA Markers

To explore whether the potential miRNA markers play important roles in FZHY treatment, we performed qPCR to measure the expression levels of 13 miRNAs from FZHY and placebo-treated HBC patients. At the transcript expression level, a statistical significance was noted for nine miRNAs before compared to after FZHY treatment, such as hsa-miR-1182, -18a-5p, -18b-5p, -193b-5p, -23a-5p, -326, -378a-5p, -564 and -760 (p < 0.05, Figure 4A). ROC analysis shows that these miRNAs are suitable to diagnose the clinical efficacy of FZHY (The area under the ROC-curve >0.60; Figure 4B). These results suggest that the miRNAs’ expressed levels are a benefit to the efficacy evaluation of FZHY formula in the HBC treatment process.

Open in a separate window

Figure 4

The potential kernel miRNAs were validated using the qPCR method. (A) Expression levels of 13 kernel miRNAs in the FZHY group. Illustrated p-values are based on pair-wise comparisons by the Mann–Whitney U-test. The result show that hsa-miR-1182, -18a-5p, -18b-5p, -193b-5p, -23a-5p, -326, -378a-5p, -564 and -760 (p < 0.05) have a statistically significant difference, while has-miR-1225-3p, -18b-3p, -1915-3p and -324-3p are insignificant between before and after trial; (B) Expression levels of 13 kernel miRNAs in the placebo group; hsa-miR-18b-5p and -193b-5p (p < 0.05) have significant differences between before and after trial. Illustrated p-values are based on pair-wise comparisons by the Mann–Whitney U-test; (C) Receiver operating characteristic (ROC) curves for hsa-miR-1182, -18a-5p, -18b-5p, -193b-5p, -23a-5p, -326, -378a-5p, -564 and -760 in the FZHY group; (D) ROC curves for hsa-miR-18b-5p and -193b-5p in the placebo group; (E) Hierarchical cluster and heat map of nine kernel miRNAs. The pathway analysis was performed using the miRPath program, and the significant pathways were determined when the p-value <0.001.

Furthermore, the miRPath (v2.0) analysis [25] shows that these nine miRNAs are mainly associated with alcoholism, systemic lupus erythematosus, lysine degradation, transcriptional misregulation in cancer, the metabolism of xenobiotics by cytochrome P450, human T-cell lymphotropic virus I (HTLV-I) infection, the Wnt signaling pathway, the PI3K-Akt signaling pathway and the MAPK signaling pathway (Figure 4E). This indicates that the FZHY can increase the expression levels of relative kernel miRNAs in HBC treatment progression and subsequently affect several biological processes.

As the placebo-control, we notice that the miRNAs hsa-miR-18b-5p and -193b-5p have significant differences in the placebo group (Figure 4C); however, their lower AUC (AUC_{(hsa-miR-193b-5p)} = 0.540 and AUC_{(hsa-miR-18b-5p)} = 0.557; Figure 4D) suggested that the placebo might have a temporary and limited effect for HBC treatment.

4.2. Clinical Specimens

In the trial, participants have been randomly given FZHY and placebo (2% FZHY), and then, 35 medical indexes (MIs) were detected from 90 FZHY- and 90 placebo-treated HBC patients, respectively. After 6 months of treatment, the same works were performed in these 180 cases. The Child–Pugh total score was calculated from ascites (Asc), bilirubin (BIL), albumin (ALB), prothrombin time (PT) and hepatic encephalopathy (HE). The clinical data of subjects are presented in Table 2.

4.3. Clinical Data Analysis

Using the R package, a random variance model was designed to test the p-value of medical indexes. In the model, i represent the distribution of each medical index expression level between before and after treatment.

where the numerator is the difference in the means of the log value i in the two classes (FZHY and placebo), n₁ and n₂ are the number of samples in the two classes, ${\tilde{\mathsf{\sigma }}}_{\mathrm{i}}^2$ is a weighted average of the usual medical index-specific variance and the average variance for all indexes and ${\widehat{\mathsf{\sigma }}}_{\mathrm{i}}$ denotes the square root of the usual within-class variance for index i. The 1/ab is an estimate of the expected variance for the inverse γ model, and 2a represents the average variance for all indexes.

4.4. miRNA Profiles Detection

The Child–Pugh score was used to identify the selected samples from FZHY-treated HBC patients. If the Child–Pugh score has a statistically-significant decrease after the trials, this patient was considered to respond to FZHY and to act as a candidate subject. Based on this, 6 pairs of serum samples (6 samples before trial and 6 samples after trial) were selected for miRNA profile detection. The miRNA profiles were generated using Agilent Human miRNA microarray V3 (Agilent Technologies Inc., Santa Clara, CA, USA). All raw data were transformed to log2, and each expression was normalized by zero mean and unit sample variance. Using the SAM (significance analysis of microarrays) method of the R package, the differentially-expressed (DE) miRNAs were identified before compare to after FZHY treatment, where the fold-change >1.5 and p < 0.001 were considered to be significant.

4.5. miRNA Target Genes Prediction

The miRNA target genes were predicted using three databases, TarBase (v7.0) [42], miRecords [43] and miRTarBase [44], which contained the largest collection of manually-curated experimental data. Furthermore, miRanda, miRDB, miRWalk and RNAhybrid programs were used to predict the non-experimental targets, where the p < 0.001 was considered to be significant. The predicted target genes were analyzed using DAVID online [20,21]; significance analysis was defined as the p-value adjusted by the false discovery rate (FDR), and gene sets containing less than 5 genes overlapping were removed. In this study, Gene Ontology (GO), pathway and the disease term with an FDR-adjusted p-value of less than 0.05 were retained.

4.6. miRNA-Target Network Construction

The differentially-expressed (DE) miRNAs and predicted targets were combined to construct the miRNA-target network. miRNAs were weighted by the fold change (|log2|); target genes were weighted based on the distribution of degrees. Subsequently, all nodes were ranked according to their weights, and the similarity was tested; then, the obtained deregulated nodes were used to remap the network. In the network, the node represents the miRNA or target; the edge represents the connection strength. The topological profiles of the network were calculated by the ClusterONE algorithm [24]; consecutive clusters were generated, which were defined as a p-value <0.0001, node size >5, quality >0.60 and density >0.05. The determined miRNA clusters were considered as kernel miRNAs and play important roles in HBC treatment progression.

4.7. Quantification of Potential Marker miRNA

The qRT-PCR was used to identify the kernel miRNAs in 360 serum samples (FZHY before n = 90 and after n = 90; placebo before n = 90 and after n = 90). The quantification of miRNA was performed with SYBR Green PCR Master Mixture (TOYOBO Biotech Co., Ltd., Osaka, Japan) according to the manufacturer’s instructions using a Rotor-Gene 6000 Real-time PCR machine (Corbett Life Science, Sydney, Australia). The specificity of each PCR product was validated by the melting curve at the end of PCR cycles. All miRNAs were validated in triplicate; the C_t was considered as the number of cycles required and for the fluorescent signal to reach the threshold. The levels of miRNAs were calculated using 2ΔC_t, where ΔC_t = C_t of the internal reference −C_t of the target miRNA. The differences in miRNAs’ expression levels between groups were compared using the Student’s t-test, and a p-value <0.05 was considered to be a statistically-significant difference.

This work was supported by the Program of the National Natural Science Foundation of China (81272295 and 81473443), the Natural Science Foundation of Shanghai Science and Technology Commission (12ZR1402200 and 14ZR1438800), the National Natural Science Foundation of China for Young Scholar (81100589), the Scheming Project of Shanghai Municipal Education Commission (2012JW25) and the Leading Project of Integrated Traditional and Western Medicine of Shanghai University of Traditional Chinese Medicine (2013). The clinical trial registration number was NCT00543426, which can be found on http://www.clinicaltrials.gov/.