Identification and Characterization of Liver MicroRNAs of the Chinese Tree Shrew via Deep Sequencing

1. Background

The tree shrew (Tupaia belangeri chinensis) is a kind of bred small animal widely distributed in regions of Southeast Asia and Southwest China. Due to its high-affinity to primates, tree shrew has been proposed as animal model, alternative to primates, in human diseases and drug screening researches, for several decades (1). Currently, there are many attempts to employ tree shrew to create animal models for the research of viral disease. Previous studies have demonstrated that tree shrew could be infected by several hepatitis viruses, especially hepatitis C virus (HCV) (2). However, a great number of obstacles in the establishment of tree shrew hepatitis C animal model still remain, especially the low and unstable infection rates. Therefore, exploring the molecular mechanism on tree shrew’s susceptibility of HCV infection, in the genomic, transcriptional and post-transcriptional levels, is of great importance.

MicroRNAs (miRNAs) are endogenous non-coding RNA molecules found in plants and animals. They functions in transcriptional and post-transcriptional regulation levels, by targeting a diverse set of cellular genes, and regulate multiple biological processes, such as cell proliferation, cell cycle, apoptosis, immune response and cellular response to stress (3). There is increasing evidence that miRNAs play important roles in virus infection, through regulating gene expression in hosts. In vitro studies have found that miRNAs, including miR-122, miR-155, miR-27a, let-7b etc., were able to impact HCV replication, host’s post-infection response and related diseases progression (4). Although the genome and transcriptome data of tree shrew have been made available (5), the sRNA data for tree shrew has not been fully investigated. The data obtained from small RNA analysis may provide valuable information for further understanding of miRNAs in the regulation of HCV infection and will be a useful resource for researches taking tree shrew as hepatitis C animal model.

2. Objectives

To present the liver microRNA data of three second filial generation tree shrews.

3. Materials and Methods

3.1. Tissue Sample Collection and RNA Isolation

Three second filial generation tree shrews were used in this study, which were provided by the Institute of Medical Biology, in Yunnan, China. Those tree shrews selected were all six-month-old and in good physical health. The fresh liver tissues isolated from the three animals were placed immediately into liquid nitrogen and then stored at -80°C, until used. All procedures involving animals were approved by the Animal Care and Use Committee of Kunming University of Science and Technology, Kunming, China. Total RNA was extracted from the liver of each tree shrew using Trizol reagent (Invitrogen, Carlsbad, CA, USA), according to the manufacturer’s instructions. The purity and concentration of isolated RNA were measured using a Thermo Scientific NanoDrop 2000c spectrophotometer (ThermoFisher Scientific, Waltham, MA, USA) and then the quality was confirmed by electrophoresis on a 1% agarose gel.

3.2. Small RNA Library Construction and Sequencing

All the high-quality RNA samples obtained were quantified and equalized, to guarantee equal amounts of RNA from each animal. The fragments of 18 - 30 nt in length were isolated and purified from the mixed RNA samples, using 15% denaturing polyacrylamide gel electrophoresis (PAGE). Subsequently, a 5’ RNA adaptor and 3’ RNA adaptor were ligated to the purified sRNA, using RNA ligase and reverse-transcribed to cDNA, using the adaptor primers. The cDNA was further amplified and fragments of around 75 bp (small RNA + adaptors) were collected from agarose gels and purified. The DNAs were confirmed by sequencing with Solexa sequencing technology (BGI, Shenzhen, China).

3.3. Data Analysis

All the raw datasets generated by the Solexa deep sequencing were analyzed, as follows. First, the low quality reads, contaminated reads involving 5’ primer contaminants, reads without 3’ primer, reads without the insert tag, reads with poly (A), and reads shorter than 18 nt were removed. The length distribution of the clean reads obtained was then summarized. Next, the clean reads were aligned against the tree shrew genome, using the SOAP (Short Oligonucleotide Alignment Program) (Beijing Genomics Institute, Beijing, China) to analyze their expression and distribution on the genome. Then, the clean reads were compared against the non-coding RNAs (ncRNAs), like ribosomal RNA (rRNA), transfer RNA (tRNA), small nuclear RNA (snRNA), and small nucleolar RNA (snoRNA), deposited in the NCBI GenBank and the Rfam databases, using Basic Local Alignment Search Tool (BLAST) to annotate the small RNA sequences. Conserved miRNAs in the tree shrew were identified using BLAST at the miRBase (http://www.mirbase.org/). Finally, the unannotated sequences were used to predict potential novel miRNA candidates by Mireap (http://sourceforge.net/projects/mireap/). The potential novel miRNA candidates generated should meet the following criteria, as reported previously (6): minimal miRNA sequence length of 18 nt, maximal miRNA sequence length of 26 nt, minimal miRNA reference sequence length of 20 nt, maximal miRNA reference sequence length of 24 nt, minimal depth of Drosha/Dicer cutting site of three nt, maximal copy number of miRNAs on reference of 20 nt, maximal free energy allowed for a miRNA precursor of -18 kcal/mol, maximal space between miRNA and miRNA* of 35 nt, minimal base pairs of miRNA and miRNA* of 14 nt, maximal bulge of miRNA and miRNA* of four nt, maximal asymmetry of miRNA/miRNA* duplex of five nt and the flank sequence length of miRNA precursor of 10 nt. The secondary structure prediction of miRNA precursors was carried out using RNA Structure 4.6 software (Mathews Lab Software, University of Rochester Medical Center, Rochester, NY, USA).

3.4. Prediction and Functional Analysis of MicroRNA Targets

The targets of the differentially expressed miRNAs were predicted using the database miRecords (http://mirecords.umn.edu/miRecords). To infer the potential functions of the differentially expressed miRNAs, we performed the functional analysis of their target genes using the database for annotation, visualization, and integrated discovery (DAVID). Functional categories were clustered using the functional annotation clustering tool, and representative gene ontology (GO) categories and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways, from each clustered set, with a P < 0.05, were selected and taken into consideration for further analysis.

3.5. Novel MicroRNAs Validation

Stem-loop real time-polymerase chain reaction (RT-PCR) was used to validate novel miRNAs. Total RNA was extracted from 10 livers of tree shrews. The cDNA was synthesized from total RNA with specific stem-loop primers. All stem-loop reverse transcriptase-polymerase chain reaction (RT-PCR) primers and miRNA-specific PCR primers are listed in Table 1. Superscript reverse transcriptase was used to process reverse transcriptase reactions that contained 1.5 mg of total RNA and 10 mM stem-loop RT-primer. The 20 μL reactions were incubated in an Applied Biosystems 2720 Thermal Cycler (Applied Biosystems, Waltham, MA, USA) at 25°C for 5 minutes, 42°C for 1 hour, and 70°C for 15 minutes. The cDNA were diluted 1:10 to perform PCR for validation. The PCR mixtures included 2 μL cDNA, 0.8 mM forward and reverse primers. The 20 μL reactions were incubated in an Applied Biosystems 2720 Thermal Cycler (Applied Biosystems, Waltham, MA, USA) in a 96-well plate, at 95°C for 5 minutes, followed by 35 cycles of 95 for 15 seconds, 62°C for 30 swconds, and 30 seconds at 72°C. The PCR products were then gel-purified and T-A cloned, and at least three independent clones were sequenced to determine the potential novel miRNAs.

4. Results

4.1. Expression Profile of MicroRNAs in Tree Shrew Liver

High-throughput sequencing was performed using small RNAs from the liver of the second filial generation tree shrew. After eliminating 3’ adapter null, 5’ adapter contaminants, insert null, smaller than 18 nucleotides (nts) and poly A sequences, we obtained 17503248 high-quality clean reads, accounting for 94.37% of the total raw reads, among which 14207005 reads were mapped to the Tupaia belangeri genome. Further analysis showed the majority of the small RNAs were 21 - 23 nt in size, accounting for 76.29% of the total clean reads (Figure 1A). MicroRNAs comprised 79.8% of the total small RNAs (Figure 1B), which accounted for only 12.8% of the total number of unique small RNAs (Figure 1C).

Profile of the Solexa Sequencing Data

A, length distribution of sequenced small RNAs. Small RNA reads with a length of 22 nt were the most abundant, followed by 21 and 23 nt reads. B, total sRNA annotations. C, unique sRNA annotations.

Figure 1.

Profile of the Solexa Sequencing Data

4.2. The Similarity of HCV Related MicroRNAs From Tree Shrew and Human

To identify conserved miRNAs in our data, the sorted clean sRNA reads were aligned with the mature miRNAs in miRBase 20.0 (http://www.mirbase.org/). As a result, 2060 unique small RNAs were identified to be matched with the known miRNAs. Expression levels of miRNAs that occupied the top 20 of all identified miRNAs are shown in Figure 2A. The most abundant conserved miRNA identified in the library was miR-122, with more than 9000000 reads, which accounted for 69.1% of the total miRNAs (Figure 2B). Four types of unique small RNA, belonging to the miR-122 family, were identified and nominated as tsh-miR-122, tsh-miR-122-5p, tsh-miR-122b and tsh-miR-122-3p, respectively. Expression levels of these four miRNAs were shown in Figure 2C. The abundance of tsh-miR-122 and tsh-miR-122-5p was much higher than that of tsh-miR-122b and tsh-miR-122-3p. The precursor of tsh-miR-122 was subjected to further analysis, by comparing with miR-122 precursors of 18 animal species, derived from the miRBase. It was shown that the sequence of the precursor of miR-122 of tree shrew was most identical with that of human (Figure 3). In addition, the secondary structure and the energy of primary miRNA of tree shrew for miR-122 were also similar with that of human (Figure 2D).

Characteristics of the Conserved MicroRNA and the MiR-122 Family in Tree Shrew

A, expression of miRNAs of which the expression level occupied the top 20 of all conserved miRNAs. B, the pie chart of the percentage of miR-122 and Let-7 in the tree shrew liver miRNAs. C, expression of four members of miR-122 family in tree shrew. D, the secondary structures of the four precursors of miR-122 in human, tree shrew, mice and rat.

Figure 2.

Characteristics of the Conserved MicroRNA and the MiR-122 Family in Tree Shrew

Sequences Alignment of the MiR-122 Precursors in 18 Tree Shrew Related Species

The sequences variation compared with human miR-122 precursors was shown in capital letters.

Figure 3.

Sequences Alignment of the MiR-122 Precursors in 18 Tree Shrew Related Species

The Let-7 family was ranked the second most abundant miRNA in tree shrew liver, with a presence of 16.9%. Among the Let-7 family, the Let-7f accounted for the largest proportion (57.32%), as shown in Figure 2B, while the Let-7a, Let-7c, Let-7b, Let-7e accounted for 19.91%, 7.61%, 5.02% and 3.60%, respectively. As reported in previous in vitro studies, Let-7b, miR-196, miR-199a and miR-130 were involved in inhibiting HCV infection, whereas miR-122, miR-141, miR-130a and miR-21were found to participate in promoting HCV infection (7). These HCV related miRNAs mentioned above were all found in miRNA data of the tree shrew liver, and the expression level is shown in Figure 4.

Tree Shrew MicroRNAs That Promote Hepacivirus Infection and Inhibit Hepacivirus Infection

The detailed read number of miRNA sequences is shown in brackets.

Figure 4.

Tree Shrew MicroRNAs That Promote Hepacivirus Infection and Inhibit Hepacivirus Infection

4.3. Identification of Novel MicroRNAs in Tree Shrew Liver

Furthermore, the sRNA sequences were aligned to the genome of Tupaia belangeri to identify potential novel miRNAs. Using the miRNA prediction software Mireap (http://sourceforge.net/projects/mireap/), sequences that met the criteria for miRNA were considered to be miRNA precursors. Ultimately, 80 novel miRNAs were identified. However, the expression level of these predicted novel miRNAs was significantly lower than that of the conserved miRNAs. The novel miRNAs, whose read number was > 10, was listed in Table 2 The read number that occupied the top 3 of all the new miRNAs was 233,133,120, respectively (Table 2). The secondary structures of these three microRNAs were shown in Figure 5A. The stem-loop RT-PCR and sequencing methods were used to further validate the novel miRNAs, using liver samples of 10 adult tree shrews. The results showed that the new three miRNA could be amplified by PCR. Moreover, the sequences of the PCR products were consistent with the microRNA data (Figure 5B).

The Secondary Structures of the Precursors of Three new MicroRNAs

A, the secondary structure of the precursors of the new microRNAs was analyzed using the software RNA structure. B, the polymerase chain reaction and sequencing results of the top 3 abundant novel miRNAs in other tree shrew samples.

Figure 5.

The Secondary Structures of the Precursors of Three new MicroRNAs

4.4. Functional Analysis of the MicroRNAs

To identify the functions and pathways that the miRNAs and their targets may affect, GO and KEGG pathway analyses were performed. GO analysis showed that 32% of the genes were annotated to biological processes, 34.4% of the genes mapped to molecular function and 33.6% of the genes were composed of cellular component. These genes ranged from elevation of cytosolic calcium ion concentration to regulation of signaling pathway, genes regulating sodium channel activity to genes regulating phosphatase binding and genes regulating microvillus to genes regulating plasma membrane (Figure 6). The findings revealed that predicted targets of miRNAs were associated with a wide variety of biological processes. Furthermore, for KEGG pathway analysis, some classic pathways reported to involve in immunoregulation were found in this study, including Toll-like receptor signaling pathway, B cell receptor signaling pathway, T cell receptor signaling pathway, Natural killer cell mediated cytotoxicity and so on (Figure 7). Intriguingly, Hepatitis C was also enriched from the predicted targets, suggesting that these targets play important roles in the process of infection and inflammation. Taken together, GO and KEGG analyses supply a better understanding of the cellular components, molecular functions and biological processes of target genes, and provided a reference for future study.

The Gene Ontology Enrichment Analysis of Predicted Targets of Differentially Expressed MicroRNAs With P < 0.05

A total of 36007 putative target genes were assigned to the GO annotations of Biological Process (11519), Cellular Component (12110), and Molecular Function (12378). The GO terms were sorted by the enrichment P value calculated by the calculating formula as previously described (8). A, biological process, B, molecular function, C, cellular component.

Figure 6.

The Gene Ontology Enrichment Analysis of Predicted Targets of Differentially Expressed MicroRNAs With P < 0.05

The Kyoto Encyclopedia of Genes and Genomes Pathway Analysis of Predicted Targets Differentially Expressed MicroRNAs with P < 0.05

The enrichment P value calculated by the calculating formula as previously described (8).

Figure 7.

The Kyoto Encyclopedia of Genes and Genomes Pathway Analysis of Predicted Targets Differentially Expressed MicroRNAs with P < 0.05

5. Discussion

It was estimated that HCV infection causes the liver-related death in more than 300000 people, annually (9). Once exposure to HCV occurs, only a minority of patients clear the acute infection by HCV-specific interferon-γ responses (10), whereas about 80% persist with life-long chronic viremia, if not successfully treated (11). Up to date, there is no vaccine for prevention of HCV infection. Interferon/ribavirin therapy is currently the mainstay in the management of HCV infection. However, the drug has been plagued with excessive toxicity, therapeutic limitations, and restricted patient availability (12). Recently, an oral antiviral agent sofosbuvir has been approved by the US food and drug administration (FDA) and European Union, for the treatment of chronic HCV infection. Although the new anti-HCV drug looks promising, the development of new therapy or preventative vaccines that definitely reduce the global HCV burden remains a formidable challenge.

Chimpanzee is the only animal model for studies of HCV infection and related innate and adaptive host immune responses (9). However, their large adult body size and long reproductive cycle create challenges for in-depth study on mechanism of HCV infection and pathogenicity. Most efforts have been directed towards the study of development of small animal models of HCV infection. Although T-cell and B-cell deficient mice grafted with human hepatocytes could support HCV infection robustly, it could not be used in studies of adaptive immunity (13). Recent studies in developing genetically humanized mice are in process, although these animal models only permit studies of specific steps in the HCV life cycle and have limited or no viral replication (14).

Genome sequencing of the Chinese tree shrew and comparison with 14 other species, including six primate species, showed that the tree shrew was clustered with primate species (5). Our previous studies had proved that CD81, scavenger receptor class B1, claudin-1 and occludin, the considered HCV receptors of tree shrew, have the mediating function to support HCV infection (15). Therefore, multiple attempts are made to employ tree shrew to create animal models for studying HCV infection (16, 17). However, although the tree shrew has the ability to be infected, the low infection rate and instable infection effect are still complicated issues, urgently needing to be addressed. Here, three second filial generation tree shrews were used to generate the liver miRNA data, which could help us in a future study to investigate the mechanisms of HCV infection in tree shrew.

MicroRNAs are important players in the establishment of HCV infection replication and its propagation in infected hepatocytes. Deregulated expression of miRNAs has been linked to the pathogenesis associated with HCV infection, by controlling the process of cell proliferation, apoptosis and migration (4). Therefore, the miRNA data will provide us more information for the utilization of tree shrew, as an animal for HCV research. The current study was designed to profile the general expression of miRNAs in the liver of tree shrew. A total of 2060 conserved miRNAs and 80 novel microRNAs were identified, and miR-122 was the most abundant miRNA in the liver of tree shrew. The miR-122 was also identified as a liver-specific, well-conserved and the most abundant microRNA in liver, accounting for about 70% in the adult liver. Many studies have previously demonstrated that miR-122 is bind to two different sites in the 5’ untranslated region (UTR) of the HCV genome and promotes viral RNA accumulation. The miR-122 plays an important positive role in the regulation of HCV replication and is a crucial target for anti-HCV therapy (18). Therefore, several miR-122 inhibitors become potential drug candidates for the treatment of HCV infection. However, miR-122 also has an essential role in maintaining liver homeostasis and differentiation, while a downregulated expression of miR-122 has been associated with liver disease. Recent studies demonstrated that decreased miR-122 have been associated with poor prognosis and metastasis of liver cancer and can even promote hepatocarcinogenesis (19). The findings mentioned above indicate that the role of miR-122 is controversial. To clarify the role of miR-122 in liver disease, it is imperative to find an applicable animal-model of HCV infection, to further research the molecular mechanisms of miR-122 in liver disease. In the current study, we demonstrated that the miRNA data generated in tree shrew was consistent with previous studies in humans. Further comparison of the secondary structures of miR-122 precursors, among tree shrew, human, mouse and rat, revealed the close relationship between tree shrew and human. Therefore, based on the evidences of the high expression level of miR-122 in the liver of tree shrew, the similarity of the secondary structure of miR-122 between tree shrew and human, and the link between miR-122 and HCV processing, it further supported tree shrew as an potential animal model for HCV infection. In addition, our findings also suggest that tree shrew may be a potential animal model to better understand the molecular mechanisms of miR-122 in liver disease caused by HCV infection in vivo.

Another interesting finding in our study was that half of the top 20 small RNAs, sorting by expression level, belonged to let-7 family. Previous studies found that the let-7 family was also strongly associated with HCV infection. Several members of let-7 family (has-let-7a, has-let-7b, has-let-7c, has-let-7d, has-let-7e, has-let-7f, has-let-7g, has-let-7i) had the potential to repress the expression of the HCV core protein in vitro (20). Another research found that Let-7b repressed HCV replication by interaction with the coding sequences of NS5B and 5’-UTR of HCV genome that were conserved among various HCV genotypes (21). Therefore, further studies on the roles of the let-7 family in HCV replication, in tree shrew, are helpful for the research on HCV infection.

Taken together, tree shrew has been proposed as a valid experimental animal to replace primates for studying HCV infection. However, there are limited usages of this animal in the field, because of the unclear mechanism on HCV infection. The miRNA data we generated and the recently annotated genome sequence of the Chinese tree shrew offer an opportunity to decipher the genetic basis of the tree shrews’ suitability, as an animal model in HCV research.

The second filial generation tree shrew, used in our study, underwent elaborate biological evolution. The indexes of microbiology, virology and parasitology can be precisely controlled. In addition, they have the advantages of relatively clear genetic background and strong disease resistance. However, there are several limitations in our current study. Due to the poor breeding potential, only three livers of tree shrews, which were 6-month-old, were used to generate the relatively limited miRNA data. To accurately characterize the expression profiles of miRNA in tree shrew, further studies, involving more liver tissue samples of different growth and development periods of tree shrew, would be necessary.