Transcriptome analysis identifies key regulators and networks in Acute myeloid leukemia.

ABSTRACT Objectives: Acute myeloid leukemia (AML) is a heterogeneous and highly recurrent hematological malignancy. Studies have shown an association between microRNAs and drive genes in AMLs. However, the regulatory roles of miRNAs in AML and how they act on downstream targets and the signaling pathway has been little studied. Methods: As to understand the mechanism of mRNA-miRNA interaction in the blood malignancy from a large scale of transcriptomic sequencing studies, we applied a comprehensive miRNA-mRNA association, co-expression gene network and ingenuity pathway analysis using TCGA AML datasets. Results: Our results showed that his-mir-335 was a critical regulatory of homeobox A gene family. PBX3, KAT6A, MEIS1, and COMMD3-BMI1 were predicted as top transcription regulators in the regulatory network of the HOXA family. The most significantly enriched functions were cell growth, proliferation, and survival in the mRNA-miRNA network. Conclusion: Our work revealed that regulation of the HOXA gene family and its regulation played an important role in the development of AML.


Introduction
Acute myeloid leukemia (AML), is the most common hematological malignancy characterized by the abnormal proliferation and differentiation of myeloid progenitor cells with high incidence and recurrence rates [1]. AML is shown to be a type of disorder in the clonal form but is also relevant heterogeneous concerning the genotypic, phenotypic, and clinical characteristics.
With the development of sequencing techniques, targeted therapy facilitated by molecular biology and genetics is a promising area for therapy of AML [2]. High throughput technologies such as RNA sequencing have enabled the revelation of global transcriptomic changes in AML, especially key molecules in the pathogenesis. Currently, studies indicated that during the pathogenesis of AML regulation of a variety of signaling pathways were overexpressed which was caused by mutations contributing to hematopoietic transformation in terms of transcriptional targets during pathogenesis, such as TET2, DNMT3A, ASXL1, IDH1 and IDH2 [3]. The different signaling pathways are also involved in the pathogenesis of AML including PI3K/ AKT pathway, MLL B-ALL gene expression signatures, and HOX gene family MicroRNAs (miRNAs) are a group of conserved short non-coding RNAs (19-25nt) [4] that regulate mRNA through targeting their seed sequences [5], playing crucial roles in biological functions during hematopoiesis and apoptosis in blood cells [6]. Although previous studies showed that miR-126, miR-196b, miR-10b, miR-21 played important roles in mutated FLT3-ITD, RUNX1, RUNX1-RUNX1T1, and CBFB-MYH11 associated AMLs, understanding of the distinctive pattern of miRNAs is still not enough to thoroughly reveal the mechanism of AML. A combined microRNA and mRNA analysis to deeper understand the pathogenesis of AML is needed.
The Cancer Genome Atlas (TCGA) project is aimed to profile genomic changes in more than 33 different cancer types, including AML, through a better understanding of the genetic basis of the disease. In this study, we investigated the TCGA transcriptomic dataset of AML to reveal the regulatory events between miRNA and mRNA, which were related to AML. Herein, as to comprehensively understand the mechanism of mRNA-miRNA interaction in the blood malignancy, we applied the TCGA AML dataset for miRNA-mRNA association analysis, co-expression gene network analysis, and pathway analysis.

Data description and preprocess
The sequencing-based mRNA and miRNA expression data of 179 AML patients were obtained from the TCGA project. The gene expression FPKM (fragment per kilobase per million mapped reads) values were calculated. Gene with a zero FPKM value in all samples was discarded. After normalization 17,300 mRNAs and 612 miRNAs remained. To avoid infinite values, a value of 0.001 was added to the FPKM value before log2 transformed for each gene.

miRNA-mRNA association analysis
The pair-wise Pearson correlation coefficient between mRNA and miRNA expression was calculated based on 17,300 mRNAs and 612 miRNAs, yielding a correlation coefficient data matrix with 17,300 mRNAs in a row and 612 miRNAs in columns. A Pearson Correlation Coefficient less than −0.6 was used as a cutoff to obtain the most probable biologically relevant miRNA-mRNA regulations. We only select the strongest negative correlations (−0.6 or less) to draw in miRNA-mRNA regulations network graphs with Cytoscape. Data preprocessing, the Pearson correlation coefficient was performed in the R (www.r-project.org/) environment with its 'base' function and 'stat' packages.

Integrity pathway analysis of target mRNA
Target mRNA identified by miRNA-mRNA association analysis (cor < −0.6) was used for functional annotation and pathway analysis to determine the biological function of these genes. The enriched pathway was determined by both significant fisher exact test (p-value < 0.05), and at least three genes were involved in the pathway.

Correlation analysis identified the miRNA-mRNA interaction network in AML
The microRNAs are regulatory, play an important role in many cancer types which may be involved in cancer initiation and development. To identify miRNA associated with AML, we download transcriptome sequencing data of AML, including both mRNA and miRNA gene expression profiles of 179 tumor samples. Given the regulatory relationship between miRNA and mRNA, we assumed that the existence of a negative correlation between miRNA and the expression of its target genes. By performing correlation analysis, we obtained a group of miRNA-mRNA pairs. We applied a cutoff of negative regulation with a correlation coefficient <−0.6 to select the most significant miRNA association pairs, which might be strongly associated with the initiation and progression of AML. The network in Figure 1 shows the miRNA-mRNA interaction network in AML. We observed that several miRNA functions as hub regulators in the network, such as has-mir-335, has-mir-146a, and has-mir-21. These miRNA have been reported to be up-regulated in AML, but the clear mechanism of the regulatory in cancer remains unknown. Notably, the expression of has-mir-135 in AML was associated with five family members of homeobox genes. As shown in Figure 2, we observed that expression of has-mir-335 was significantly associated with HOXA with correlation <−0.6, which might indicate has-mir-335 targeted the HOXA family members.

Function annotation for miRNA regulated mRNA in AML
To further characterize the role of miRNA and their target mRNA in AML, we perform a function annotation for all the mRNA in the regulation network. As shown in Table 1, most of the mRNA in the network of Figure 1 is involved in cell growth, proliferation, and survival.
We also perform an IPA pathway analysis by using these target mRNA. The enriched pathways were listed in Table 2. Salvage Pathways of Pyrimidine Deoxyribonucleotides was the top pathway with a −log (pvalue) 2.63, ATP and TK2 gene in the regulation network that was involved in this pathway. Previous studies [7,8] have demonstrated the important role of deoxyribonucleotides metabolism in leukemia and other cancer types. Besides, Glutamyl Cycle and Colanic Acid Building Blocks Biosynthesis were also significantly enriched in leukemia. Notably, we observed that ATP was involved in most enriched pathways. It was indicated ATP associated with energy metabolism might play an important role in AML.

Gene-expression regulator prediction for Targeted mRNAs in the regulation network
We also use IPA to predict the potential regulator of these targeted genes. As shown in Table 3, PBX3 was predicted as a top transcription regulator of six target genes in the regulatory network HOXA10, HOXA3, HOXA5, HOXA6, HOXA7, and HOXA9. HOXA and PBX3 are independent predictors of poor survival in patients with acute myeloid leukemia [9,10] and also associated with the chemotherapy of AML. Except for PBX3, there are also 12 genes that were predicted as a transcription regulator of several family members of homeobox A, including the KAT6A gene which was reported as a commonly rearranged gene in AML [11]. Our results suggested regulation of HOXA gene expression might be associated with the development of AML.

Discussion
In this study, we investigated a transcriptome dataset of AML to reveal the regulatory events between miRNA and mRNA, which were related to AML. Our results showed that his-mir-335 was a key regulatory of homeobox A gene family. Regulation of the HOXA gene family may play an important role in the development of AML. Studies have demonstrated that genetic alterations accumulation might contribute to the initiation of AML. Despite an increased understanding of AML, our knowledge about transcriptome regulation in AML still needed to be improved. Firstly, our study predicted a series of miRNA regulators in AML by performing correlation analysis between miRNA and mRNA, including has-mir-130a, has-mir-146a, has-mir-335, has-has-mir-21, has-mir-126, and has-mir-181a/b/c/d. Previous studies had revealed these miRNAs played important roles in the development of AML. Has-mir-130a was involved in the megakaryocytic differentiation and erythropoiesis, and its expression was also related to the differentiation patterns in AML [12]. Has-mir-146a was a regulator of apoptosis [13], and its expression was associated with patient prognosis of AML. The increase of has-mir-146a can enhance the sensitivity of leukemic blast cells to cytotoxic drugs [14]. Besides, has-mir-335, has-mir-21, and has-mir-126 were significantly higher expressed. Expression of has-mir-335 was associated with patient clinical outcomes [15]. Has-mir-181 family members were reported to be associated with the myeloid differentiation, and has-mir-181a could enhance cell proliferation [16] and correlated with the morphological sub-type of AML [17]. However, research also demonstrated that has-mir-181a acts as a tumor suppressor in the pathogenesis of AML [18]. Our study provided further evidence for the regulatory role of these miRNAs in AML.
Next, we observed that the target genes of hub miRNAs were related to cell growth, proliferation, death, and survival. Pathway analysis demonstrated that the ATP enzyme was significantly affected by miRNA regulation. The ATP was involved in biosynthesis pathways which were significantly enriched. Salvage Pathways of Pyrimidine Deoxyribonucleotides was a pathway related to the biosynthesis of deoxyribonucleotides. ATP played a very important role in this biological procedure [19][20][21][22][23]. Deoxyribonucleotides were essential to DNA synthesis/repair. The previous study demonstrated that the biosynthesis of deoxyribonucleotides was an important target for anticancer therapy [24,25].
Finally, we reveal that expression six homeobox gene family members were related to has-mir-335. The dysregulation of HOX genes is associated with AML [26]. Recent work has shown that their overexpression is important in the formation of malignancy [27,28]. Our regulatory prediction results showed that PBX3 was a key regulatory of HOXA genes. Homeobox family, including four members (A, B, C, and D) along with 39 homeotic or Hox genes, are conserved and Table 3. IPA gene-expression regulator prediction for Targeted mRNAs in the regulation network.

Data availability
The sequencing-based mRNA and miRNA expression data of 179 AML patients were obtained from the TCGA project (http://cancergenome.nih.gov/).

Disclosure statement
No potential conflict of interest was reported by the authors.

Funding
This work was supported by the Science Technology Bureau of Shaoxing (2018C30022) and the Shaoxing Shangyu Municipal Funding (Funding-ID-2015-11-07).