MicroRNAs (miRNAs), a class of small noncoding RNAs, regulate the expression of other genes by targeting protein-coding transcripts in the post-transcriptional regulation (1). Approximately, 1,500 miRNAs were annotated in the miRNA database (miRBase) (2) and up to 60% of protein-coding genes were estimated to be regulated by miRNAs (3). miRNAs interact with 3'UTR of messenger RNAs (mRNAs) through base pairing and lead to their degradation, destabilization, or repression of translation through RNA-induced silencing complex (RISC), a complex of multiple proteins and miRNA-mRNA adduct (4,5). Occasionally, these miRNAs can also up-regulate gene expression (6).
It has been reported that miRNAs are essential for normal mammalian development and are involved in many biological processes, such as cell proliferation, differentiation, apoptosis, and metabolism, and their involvement in cancer has sparked increased interest in miRNA biology (7-9). Convincing evidence has indicated that miRNAs are key regulators of hematopoiesis (10,11). Leukemia is a group of hematopoietic neoplasms featuring impaired hematopoiesis and bone marrow failure caused by clonal expansion of undifferentiated myeloid precursors (12). The disease is frequently found to have recurrent chromosomal aberrations and gene mutations, and mostly carry driver mutations relevant to patient clinical outcomes (13). It is grouped by how quickly the disease develops (acute or chronic) as well as by the type of blood cell that is affected (lymphocytes or myelocytes). There are four main types of leukemia, including acute lymphocytic leukemia (ALL), chronic lymphocytic leukemia (CLL), acute myelocytic leukemia (AML), and chronic myelocytic leukemia (CML).
Recently, most researchers use next-generation sequencing (NGS) to profile miRNAs and then identify novel miRNAs associated with a disease of interest (14,15). However, most cohort studies were focused on the major populations such as Asia, America, Europe people etc. Little were known about the molecular profile changes (DNA mutations, RNA transcript and miRNAs, etc.) in the Chinese minorities. The Naxi people are one of the unique ethnic groups in Yunnan and live in a high concentration. In the present study, we analyzed differential expression profiles of miRNAs by quantitative small RNA NGS in six leukemia patients and one normal control. All these people are from Naxi ethnic minorities located in the northwest of Yunnan province, and considered to be an area of low incidence. This approach facilitated in identifying several miRNAs relevant to different leukemia pathogenesis. Additionally, we also verified six miRNAs by quantitative real time-polymerase chain reaction (qRT-PCR) in human leukemic cell lines.
Patient cohort and sample collection
With informed consent and ethics approval, six leukemia patients (case 2 and case 3 have AML; case 4, case 5, and case 7 have ALL; case 6 has CML) and one healthy person (case 1) were recruited from Tumor Hospital in Yunnan Province. All seven person are from Nakhi (Naxi) ethnic minorities, which considered to be an area of low incidence. The samples were collected in 1.5 mL microcentrifuge tubes, immediately snap-frozen in liquid nitrogen and stored at −80 °C until further processing.
RNA isolation and deep sequencing
Total RNA isolation was carried out from peripheral blood or cell lines using TRIzol® Reagent (Invitrogen) according to the manufacturer’s instructions. The RNA concentration was measured by Nanodrop ND-100 (Thermo Scientific, Waltham, USA), and RNA quality was assessed by the Agilent 2100 Bioanalyzer (Agilent, Santa Clara, USA). Small RNA sequencing libraries were constructed following IlluminaHTruSeqTM Small RNA Sample Preparation protocol. In brief, 3' and 5' adapters were ligated to sRNA population and ligated RNAs were purified by running on urea-PAGE. Modified sRNAs were reverse transcribed to generate cDNA libraries and PCR amplified to add unique index sequence to each library.
Sequencing data analysis
The sequences were converted into FASTQ format and de-multiplexing were performed using Illumina bcl2fastq2 software version 2.17 (http://support.illumina.com/downloads/). Adaptor trimming was performed using the FASTQ Toolkit App of Illumina BaseSpace (http://basespace.illumina.com/apps/). Quality of the sequenced reads before and after adapter trimming was evaluated using FastQC software (http//www.bioinformatics.babraham.ac.uk/projects/fastqc/). Cleaned sequences were aligned to the most recent mirBASE database release 21 (http://mirBase.org/) using the Small RNA App of Illumina BaseSpace. Quantification of miRNA expression was performed by counting aligned reads to miRNA genes.
Novel miRNA prediction
The novel miRNA prediction strategy is based on first removing all known RNAs including those derived from exonic regions and then identifying those that are derived from intronic and intergenic regions. Identifying novel miRNAs from sequencing data commonly apply a combination of (I) evaluating miRNA secondary structures; and (II) ranking miRNA candidates by utilizing existing annotation or evolutionary conservation of the mature miRNA sequence. Structure threshold parameters of most algorithms are often optimized based on animal miRNA precursor structures.
Prediction of miRNA targets
The most probable targets of the differentially regulated miRNAs were predicted by two criteria. (I) Prediction by established target prediction software, which facilitate in identifying all the target genes using intersection between miRanda and miTarget; (II) inverse correlation in expression pattern between miRNA and coding genes. This approach compares the expression of putative genes with the different regulated miRNAs, and only genes that show inverse correlation to the miRNA levels were considered as most putative targets of the select miRNA.
qRT-PCR was performed to determine expression of six miRNAs: hsa-miR-181b-5p, hsa-miR-181a-3p, hsa-miR-181a-5p, hsa-miR-342-3p, hsa-miR-450a-5p, hsa-miR-1255a in three different human leukemic cell lines (Jurkat, HL-60, K-562). cDNA synthesis was performed using M-MLV Reverse Transcriptase (Promega) according to the manufacturer’s protocol. All samples were run in three replicates. Ct values for miRNAs were normalized against U6 and the relative expression was calculated using 2−ΔCt (Gene-U6) method. All the RT-qPCR primers are provided by RIBOBIO (http://www.ribobio.com/), the item number was listed in Table S1.
All bioinformatics-associated statistical analyses were performed using the R package for statistical computing. One-way analysis of variance (ANOVA) was calculated and analyzed using GraphPad Prism. P<0.05 is considered as statistically significant.
Small RNAs sequencing and annotation
sRNAs from six leukemia patients (case 2 to case 7) and 1 normal control (case 1) were used and size selected by gel electrophoresis, then sequenced using a Solexa platform (Illumina, San Diego, USA). First, we removed the adaptor sequences from the sequence reads, and only those reads that were greater than 10 nucleotides were considered for further analysis. On average, ~13 million reads mapping to the human genome were obtained per sample (Figure 1A). The length of the detected sequences varied between 18 and 32 nucleotides in all these samples (Figure S1), similar to some other studies in animals (18–35 nucleotides) and plants (18–30 nucleotides). Then, we clustered these sequences on the basis of sequence similarity and performed similarity searches using specific databases (miRNAs, rRNA, tRNA, snRNA, snoRNAs). The frequency of reads mapping to miRNAs ranged from 59% to 71% and had a median of 65% (Figure 1B). In total, 1,392 known miRNAs (miRBase 21) were detected in at least one of seven samples sequenced (Table 1). To detect the proportion of mRNA derived by degradation, we checked if there are some sequences derived from intronic and exonic regions, and found that the median read frequency is 7.572% (Figure 1C). Additionally, a small fraction of reads mapped to rRNA, tRNA, snRNA, snoRNA, together comprising a frequency of ~0.8% (Figure 1C).
Identification of candidate novel miRNA genes
The representative hairpin structure of miRNA precursor is normally used to identify putative novel miRNAs. In this study, we integrated two miRNA prediction software, miRNAEVo (16) and mirdeep2 (17) for the new miRNA prediction. As shown in Table 2. Ninety-nine mature miRNAs, 26-star miRNAs, and 102 miRNA precursors were detected in at least one of the seven samples sequenced. Almost the same number of candidate novel miRNAs was obtained from each sample. During miRNA processing, Dicer enzyme plays an important role from precursor miRNAs to mature miRNAs. Due to the specificity of the enzyme cleavage sites, the first base of miRNAs might display bias. So, we then analyzed the first base distribution of miRNAs in different length, and the base distribution of all miRNAs. As observed, the nucleotide distribution of 18–30 nucleotides novel miRNAs at the first nucleotide position had a strong bias (Figures 2,S2). Additionally, the base distributions of all miRNAs display that there are bias at each position (Figures 2,S3).
miRNA expression pattern
In order to identify miRNAs relevant to different types of leukemia, we divided our samples into four groups, acute leukemia (jxzl), chronic leukemia (mxzl), myeloid leukemia (sxzl) and lymphocytic leukemia (lxzl). To compare the expression levels between different samples, it is necessary to normalize against the read count by calculating a normalization factor. Here, the absolute sequence reads were transformed into transcript abundance by normalizing the data in “transcripts per million (TPM)” for each library first. The expression levels ranged from less than 10 to more than 100,000 counts (Figure S4). Thus, the sequencing data revealed a wide range of expression levels spanning five orders of magnitude.
Using an arbitrary P value <0.05 and two-fold change cut-off, we found that in comparison to control, 21 and 125 miRNAs were up-regulated in jxzl and mxzl respectively, while 5 and 53 miRNAs were down-regulated in jxzl and mxzl respectively (Figure 3; Tables S2-S7). Venn diagram displayed that 18 miRNAs specifically regulated in jxzl and 170 miRNAs in mxzl (Figure S5). We then compared the differently expressed miRNAs between jxzl and mxzl, and found that 120 miRNAs were up-regulated and 203 miRNAs were down-regulated in jxzl. The up-regulated genes include six novel miRNAs (novel_581, novel_133, novel_692, novel_224, novel_401, novel_94) and several known miRNAs, such as hsa-miR-103a-2-5p, hsa-miR-320b, hsa-miR-320c, hsa-miR-1180-3p, hsa-let-7a-5p, hsa-let-7b-5p, hsa-let-7c-5p, and hsa-let-7d-5p. And the downregulated genes include novel_498, novel_377, novel_580, hsa-miR-582-5p, hsa-miR-582-3p, has-miR-450b-5p, and has-miR-450a-5p (Tables S2-S7). These differentially expressed miRNAs might play crucial role in regulating the rate of the disease development.
Next, we compared the differently expressed miRNAs between the sxzl or lxzl with control. As shown, 9 and 23 miRNAs were up-regulated in sxzl and lxzl, 2 and 6 miRNAs were down-regulated in sxzl and lxzl (Figure 3; Tables S2-S7). When compared the differently expressed miRNAs in sxzl and lxzl, the results showed that 29 known miRNAs were up-regulated and 10 known miRNAs were downregulated, including hsa-miR-181-5p, hsa-miR-181a-3p, hsa-miR-181b-5p, hsa-miR-342-3p, hsa-let-7f-1-3p, and hsa-miR-30a-5p (Tables S2-S7). These miRNAs might play an important role in the specific type of blood cell development.
Predicted targets of differentially regulated miRNAs
In order to understand the biological role of the differentially regulated miRNAs, it is necessary to figure out the putative targets of the detected miRNAs. In this study, two approaches were employed for target prediction. First, predict target using two popular software tools, miRanda, and miTarget. Then, compare the expression levels between miRNA and putative target genes. Those genes that coexist in two databases, and show inversion correlation to its miRNAs were considered to be the targets. KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway analysis was performed for predicted targets. We analyzed the pathways of differently expressed genes between jxzl and control, and found that the most significant gene-enrichment included “pathways in cancer”, “Rap1 signaling pathways”, “Ras signaling pathways”, and “endocytosis”. However, when compared between mxzl and control, the most significant gene-enrichment was metabolic pathways. Similar results were found in jxzl and mxzl, which suggests that the metabolic pathway is responsible for the development of leukemia. Then we compared the pathways of differently expressed genes between sxzl and control, between lxzl and control (Figure 4). The results display that cancer development is the most significant pathway, which indicates that the miRNAs play an important role in determining the type of cancer and its development.
Validation of candidate miRNAs by qRT-PCR
To validate these findings, we examined the expression of six miRNAs, hsa-miR-181b-5p, hsa-miR-181a-3p, hsa-miR-181a-5p, hsa-miR-342-3p, hsa-miR-450a-5p, and hsa-miR-1255a in three different human leukemic cell lines (Jurkat, HL-60, K-562) by RT-qPCR (Figure 5). As known, Jurkat is acute lymphoid leukemia, HL-60 is AML, and K-562 is CML and consistent with the small RNAseq results, we found that hsa-miR-181a-5p and hsa-miR-181b-5p displayed higher expression in Jurkat and HL-60 than that in K-562. However, the expression level of hsa-miR-342a-3p in Jurkat was higher than other two (HL-60 and K-562). And hsa-miR-181a-3p showed different expression in three cell lines (Jurkat had highest and K-562 had lowest). Additionally, both hsa-miR-450a-5p and hsa-miR-1255a displayed low expression in all three cell lines. These results suggest that hsa-miR-181b-5p, hsa-miR-181a-3p, hsa-miR-181a-5p, and hsa-miR-342-3p showed different expression patterns in different types of cancer cells, and hsa-miR-450a-5p, hsa-miR-1255a were dysregulated in all the leukemia cells.
The study of miRNA expression has a certain guiding role in clinical practice as the abnormal expression of miRNA may be used as a prognostic indicator for leukemia patients (18,19). However, this is based on a large number of races and there was very little research on ethnic minorities. To provide accurate treatment for minority people, the differences between Han people in genetic characteristics and molecular maps need to be analyzed. In this study, 6 patients from the Naxi ethnic minority were selected for miRNA analysis, the Naxi ethnic group lives in the western Yunnan Plateau in the northwest of Yunnan Province, with unique national culture and regional characteristics. This is the first report of the Naxi ethnic miRNA expression. In spite of the small number of examples, it has laid a certain foundation for the follow-up study.
In this study, six patients from the Naxi ethnic minority were selected for miRNA analysis, the Naxi ethnic group located in the western Yunnan plateau, the northwest of Yunnan province, with a unique national culture and regional characteristics. Our paper, for the first time, studied the miRNA expression in Naxi people. Although the number of samples is limited, it lays the foundation for future research.
Certain investigations indicated that some miRNAs can affect the directional differentiation of hematopoietic cells, and associated with the pathogenesis of leukemia and lymphoma. Understanding and clarifying the relationship between miRNAs and leukemia cells can help to reveal the molecular mechanism of leukemia, and miRNAs and specific target genes can also be used to achieve targeted therapeutics. Till now, several miRNAs expressed tags, as well as miRNAs, responsible for the AML and CLL have been reported (18,20,21). For example, it was found that frequent deletions and down-regulation of miRNA15 and miRNA16 were happened in ALL (22); Garzon et al. found that the expression of miR-191 and miR-199a was up-regulated in AML, and its high expression was associated with poor prognosis in AML patients (18). Dixon-McIver et al. also confirmed that t(15;17) chromosome heterotopia is related to the upregulation of miR-127, miR-154, miR-299, miR-323, miR-368 and miR-370 (19).
In the present study, deep sequencing technology was used to quantify miRNA expression in six leukemia patients from Naxi ethnic. From these samples, we observed several interesting findings, including the discovery of novel miRNAs and a valuable list of differentially expressed miRNAs in chronic and acute leukemia, as well as lymphocytic and myelocytic leukemia. Our major findings are a list of expressed miRNAs (1,392 known and 125 candidate novels) in leukocytes and the discovery that the differential expression of miRNAs is most abundant in the pathway of cancer. These 125 candidate miRNAs should be further studied to exclude the false positive results and find the specific miRNAs after more samples would be collected. Our qRT-PCR results show convincingly that such miRNAs are likely to be real and not a result of an artifact of sequencing.
The known miRNAs found to be enriched in acute leukemia were members of the let-7 family, miR185 and miR320, which were also reported in peripheral blood mononuclear cells using “Taqman miRNA assay” (23). Therefore, our results are in agreement with previous studies based on a different methodology. There are multiple mechanisms that are likely to regulate miRNA levels similar to that of mRNAs.
Additionally, we found for the first time that novel-94, 692, 581, 580, 498, 401, 377, 224 and 133 showed abnormal expression in leukemia patients. We also validate our RNAseq data by RT-qPCR in three leukemia cell lines, the expression of hsa-miR-181b-5p, hsa-miR-181a-3p, hsa-miR-181a-5p, hsa-miR-342-3p, hsa-miR-450a-5p, and hsa-miR-1255a were consistent with our sequence data. Most interestingly, WT1 gene was predicted to be the target of miR-1255a by TargetScan v7.1. As known, WT1 was highly expressed in peripheral blood or narrow bone of most AML patients, and supposed to be the marker of detecting minimal residual disease (MRD). Our results indicated that the expression of miR-1255a was significantly down-regulated in AML patients, which suggested that the WT1 gene might be regulated by miR-1255a. This new finding provides a new direction for the further research.
Funding: This research was supported by grants from National Natural Science Foundation of China (81360089) and Scientific and Technological Commission of Yunnan Province (2015FB072).
Conflicts of Interest: The authors have no conflicts of interest to declare.
Ethical Statement: The study was approved by the Tumor Hospital in Yunnan Province (2015FB072).
- Kusenda B, Mraz M, Mayer J, et al. MicroRNA biogenesis, functionality and cancer relevance. Biomed Pap Med Fac Univ Palacky Olomouc Czech Repub 2006;150:205-15. [Crossref] [PubMed]
- Griffiths-Jones S, Grocock RJ, van Dongen S, et al. miRBase: microRNA sequences, targets and gene nomenclature. Nucleic Acids Res 2006;34:D140-4. [Crossref] [PubMed]
- Friedman RC, Farh KK, Burge CB, et al. Most mammalian mRNAs are conserved targets of microRNAs. Genome Res 2009;19:92-105. [Crossref] [PubMed]
- Rana TM. Illuminating the silence: understanding the structure and function of small RNAs. Nat Rev Mol Cell Biol 2007;8:23-36. [Crossref] [PubMed]
- Hammell CM. The microRNA-argonaute complex: a platform for mRNA modulation. RNA Biol 2008;5:123-7. [Crossref] [PubMed]
- Ørom UA, Nielsen FC, Lund AH. MicroRNA-10a binds the 5'UTR of ribosomal protein mRNAs and enhances their translation. Mol Cell 2008;30:460-71. [Crossref] [PubMed]
- Tsuchiya S, Okuno Y, Tsujimoto G. MicroRNA: biogenetic and functional mechanisms and involvements in cell differentiation and cancer. J Pharmacol Sci 2006;101:267-70. [Crossref] [PubMed]
- Tang JT, Fang JY. MicroRNA regulatory network in human colorectal cancer. Mini Rev Med Chem 2009;9:921-6. [Crossref] [PubMed]
- Sassen S, Miska EA, Caldas C. MicroRNA: implications for cancer. Virchows Arch 2008;452:1-10. [Crossref] [PubMed]
- Ling H, Fabbri M, Calin GA. MicroRNAs and other non-coding RNAs as targets for anticancer drug development. Nat Rev Drug Discov 2013;12:847-65. [Crossref] [PubMed]
- Larsen MT, Häger M, Glenthøj A, et al. miRNA-130a regulates C/EBP-ε expression during granulopoiesis. Blood 2014;123:1079-89. [Crossref] [PubMed]
- Lee YG, Kim I, Oh S, et al. Small RNA sequencing profiles of mir-181 and mir-221, the most relevant microRNAs in acute myeloid leukemia. Korean J Intern Med 2019;34:178-83. [Crossref] [PubMed]
- Papaemmanuil E, Gerstung M, Bullinger L, et al. Genomic Classification and Prognosis in Acute Myeloid Leukemia. N Engl J Med 2016;374:2209-21. [Crossref] [PubMed]
- Schee K, Lorenz S, Worren MM, et al. Deep Sequencing the MicroRNA Transcriptome in Colorectal Cancer. PLoS One 2013;8:e66165. [Crossref] [PubMed]
- Vaz C, Ahmad HM, Sharma P, et al. Analysis of microRNA transcriptome by deep sequencing of small RNA libraries of peripheral blood. BMC Genomics 2010;11:288. [Crossref] [PubMed]
- Wen M, Shen Y, Shi S, et al. miREvo: an integrative microRNA evolutionary analysis platform for next-generation sequencing experiments. BMC Bioinformatics. 2012;13:140. [Crossref] [PubMed]
- Friedländer MR, Mackowiak SD, Li N, et al. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res 2012;40:37-52. [Crossref] [PubMed]
- Garzon R, Volinia S, Liu CG, et al. MicroRNA signatures associated with cytogenetics and prognosis in acute myeloid leukemia. Blood 2008;111:3183-9. [Crossref] [PubMed]
- Dixon-McIver A, East P, Mein CA, et al. Distinctive patterns of microRNA expression associated with karyotype in acute myeloid leukaemia. PLoS One 2008;3:e2141. [Crossref] [PubMed]
- Luo HC, Zhang ZZ, Zhang X, et al. MicroRNA expression signature in gastric cancer. Chin J Cancer Res 2009;21:74-80. [Crossref]
- Navarro A, Gaya A, Martinez A, et al. MicroRNA expression profiling in classic Hodgkin lymphoma. Blood 2008;111:2825-32. [Crossref] [PubMed]
- Calin GA, Dumitru CD, Shimizu M, et al. Frequent deletions and down-regulation of micro- RNA genes miR15 and miR16 at 13q14 in chronic lymphocytic leukemia. Proc Natl Acad Sci U S A 2002;99:15524-9. [Crossref] [PubMed]
- Liang Y, Ridzon D, Wong L, et al. Characterization of microRNA expression profiles in normal human tissues. BMC Genomics 2007;8:166. [Crossref] [PubMed]