Expression Profiles of tRNA-Derived Fragments and Their Potential Roles in Multiple Myeloma.

BACKGROUND: Multiple myeloma is an incurable hematologic malignancy. The discovery of mechanisms may help to find new therapeutic targets and prolong survival. tRNA-related fragments (tRF) and tRNA halves (tiRNA) are small RNA derived from tRNAs and implicated in a wide variety of pathological processes, including cancer initiation and progression. However, there are almost no research reporting the role of these tRNA derived fragments in myeloma as far as we know.
METHODS: In this study, we proposed the expression profiles of tRFs/tiRNAs in myeloma through RNA-sequencing in new diagnosed myeloma and healthy donors. The expression of selected tRFs/tiRNAs was further validated in clinical specimens by qPCR. Bioinformatic analysis was performed to predict their potential roles in multiple myeloma.
RESULTS: A total of 33 upregulated tRFs/tiRNAs and 22 downregulated tRFs/tiRNAs were identified. GO enrichment and KEGG pathway analysis were performed to analyze the functions of one significantly up-regulated and one significantly down-regulated tRNA-derived fragments. tRFs/tiRNAs may be involved in MM initiation and progression.
CONCLUSION: We concluded that tRFs/tiRNAs were dysregulated and could be potential biomarkers for multiple myeloma.

Introduction

Multiple myeloma (MM) is the second most common hematological malignancy and is characterized by the accumulation of malignant plasma cells in the bone marrow.1 With the advent of new drugs, such as proteasome inhibitors and immunomodulatory agents, the outcomes of patients with MM have steadily improved. Despite improvements in response and survival rates, MM remains a costly and incurable disease.2,3 Elucidating the mechanisms governing MM development is pivotal to optimizing the clinical efficacy and prolonging survival. Myeloma cells are dependent on multiple factors at all stages of tumor development and progression. Among these factors, noncoding RNAs (ncRNAs) play critical roles.4,5

The largest portion of the mammalian genome is transcribed into RNA species with little or no coding potential, which are known as ncRNAs. An explosion of interest in RNA biology has identified a variety of regulatory functions for ncRNAs.6,7 With the rapid development of high-throughput sequencing in recent years, a new class of small ncRNAs derived from tRNAs has attracted increasing attention from researchers.

According to the cleavage site, these tRNA derived fragments can be broadly classified into two main groups: tRNA related fragments (tRFs) generated from mature or precursor tRNA and tRNA halves (tiRNAs) generated by specific cleavage in the anticodon loops of mature tRNA.8 On the basis of mapped positions on the precursor or mature tRNA transcript, tRFs/tiRNAs are subdivided into five types: tRF-5, tRF-3, tRF-1, tRF-2 and tiRNA. tRFs/tiRNAs are involved in various molecular processes through translation regulation, epigenetic regulation and RNA silencing.9,10

Accumulating evidence supports the significance of tRFs/tiRNAs in multiple human diseases, including cancer. For example, some sex hormone-dependent tRNA-derived RNAs were found to enhance cell proliferation in breast cancer cells and to act as oncogenes in breast cancer.11 tRF-1001, which is present in the cytoplasm and produced by the tRNA 3ʹ-endonuclease ELAC2, is essential for the proliferation of prostate cancer cells and colorectal cancer cells.12 Upregulated or downregulated tRFs/tiRNAs are both involved in the occurrence and development of solid tumors.13,14 In hematological malignancies, limited data have demonstrated the potential roles played by tRFs/tiRNAs in leukemia and lymphoma.15,16 To the best of our knowledge, few studies have described the roles played by tRFs/tiRNAs in the mechanism underlying MM. One study showed that smoldering multiple myeloma patients with lower expression levels of 3′-tRF-LeuAAG/TAG, i-tRF-GluCTC, and i-tRF-ProTGG were at great risk of suffering from MM.17 However, the expression profile of tRFs/tiRNAs in myeloma patients compared to healthy people has not been determined.

In this study, we investigated the expression profiles of tRFs/tiRNAs in paired newly diagnosed MM and healthy donor samples through RNA-sequencing and quantitative real-time PCR (qPCR) verification. Next, using bioinformatics analysis, we characterized the functions of these RNAs in the initiation and progression of MM. These findings may provide potential biomarkers and therapeutic targets for MM.

Materials and Methods

Clinical Samples

This study was approved by ethics committee of the third Xiangya hospital. Twenty new diagnosed MM(NDMM) patients and 18 healthy donors were enrolled, respectively. The diagnosis of NDMM was based on the diagnostic criteria of symptomatic multiple myeloma defined by National Comprehensive Cancer Network (NCCN). Bone marrow was obtained from NDMM patients and healthy donors after informed consent.

RNA Extraction and Quantitative RT-PCR Validation

Firstly, plasma cells were enriched from bone marrow samples using the anti-human CD138 magnetic-activated cell sorting (MACS) beads (Miltenyi, Germany). The purity of separated plasma cells was above 90% (identified by flow cytometry). Total RNA was then extracted from enriched plasma cells using TRIzol (Invitrogen, USA) following the instruction manual. Prepared RNA was preserved at −80°C. RNA samples were qualified and quantified by agarose gel electrophoresis and NanoDrop ND-1000 (NanoDrop, USA). qRT-PCR was implemented using ViiA 7 Real-time PCR System (Applied Biosystems) and 2 × PCR Master Mix. U6 small nuclear RNA (snRNA) was used as internal reference. The fold changes were calculated and normalized by U6. The primers used are listed in Table 1. RNA quality is shown in . RNA integrity number (RIN) was evaluated by Agilent BioAnalyzer 2100 and was all above 7.

Table 1

Primers for qRT-PCR Validation of Candidate tRFs/tiRNAs

	Primer
iRNA-1:34-Glu-TTC-2	F: 5ʹAGGCTCCCACATGGTCTAGC3’
	R: 5ʹCAGTGCAGGGTCCGAGGTAT3’
tRF-60:76-Arg-ACG-1-M2	F: 5ʹATTATCCTGGCTGGCTCGCC3’
	R: 5ʹCCGAGGTATTCGCACTGGA3’
tRF-1:23-Ala-AGC-1-M3	F: 5ʹTAAACCGGGGGTATAGCTCAGT3’
	R:5ʹCAGTGCAGGGTCCGAGGTAT3’
tRF-1:16-Lys-TTT-3-M2	F: 5ʹCTTTGCCCGGATAGCTCAGT3’
	R: 5ʹCAGTGCAGGGTCCGAGGTAT3’
tRF-1:22-Leu-AAG-4-M2	F: 5ʹATTGGTAGCGTGGCCGAGT3’
	R: 5ʹCGCAGGGTCCGAGGTATTC3’
tRF-+1:T14-Gly-TCC-1	F: 5ʹACCCTGCGGTACCACTTTGT3’
	R: 5ʹCGCAGGGTCCGAGGTATTC3’
U6	F: 5ʹGCTTCGGCAGCACATATACTAAAAT3’
	R: 5ʹCGCTTCACGAATTTGCGTGTCAT3’

Library Preparation and Sequencing

Firstly, total RNA samples were pretreated to remove RNA modifications that interfere with small RNA-seq library construction. The libraries then were denatured as single-stranded DNA molecules, captured on Illumina flow cells, amplified in situ as sequencing clusters and sequenced for 50 cycles on Illumina NextSeq 500 system following the manufacturer’s instructions. Image analysis and base calling were conducted by Solexa pipeline v1.8 (Off-Line Base Caller software, v1.8).

Sequencing quality was examined by FastQC. The abundance of tRFs/tiRNAs was evaluated by sequencing counts and normalized as counts per million of total aligned reads (CPM). R package edgeR was used to calculate the differentially expressed tRFs/tiRNAs. Hierarchical clustering, Volcano plots, Pie plots and Venn plots were constructed by R or Perl.

Bioinformatic Analysis of tRFs/tiRNAs

We used TargetScan algorithms to predict gene targets of tRFs/tiRNAs.18 Gene Ontology (GO) and Kyoto Encyclopaedia of Genes and Genomes (KEGG) function enrichment analysis were conducted using DAVID Bioinformatics Resources 6.8.19 Graphing was performed with Cytoscape 3.7.2 and GraphPad Prism 8.0.1.

Statistical Analysis

Analysis of data was performed with two-tail unpaired Student t-test using SPSS 23.0 software. qPCR validation was presented as mean ± standard deviation. The statistical significance level was set at P value < 0.05.

Results

General Features of tRFs/tiRNAs Expression in MM

We performed tRFs/tiRNAs sequencing on bone marrow specimens from five MM patients and five normal donors. Clinical characteristics of MM patients and healthy donors are shown in Table 2. After Illumina quality control and 5ʹ, 3ʹ-adaptor trimmed, reads with length 14nt ~ 40nt were retained. tRNA isodecoders share the same anticodon but have differences in their body sequence. All the mature tRNAs and pre‐tRNAs were found to produce tDRs in our profiles. We counted the number of significantly differentially expressed subtype tRFs/tiRNAs against tRNA isodecoders which the CPM of the sample is not less than 20 (Figure 1A and B). tRFs/tiRNAs may be produced from different tRNAs by cleavage into the fragments with identical sequences. The number of tRF-3a and tRF-3b, which mainly occur in the cytoplasm, decreased in myeloma patients. The pie charts depict the distribution of the number for each subtype tRFs/tiRNAs (Figure 1C and D).

Table 2

Baseline Characteristics of MM Patients and Healthy Donors

	MM (n=20)	Healthy Donors (n=18)
Age
Median - yr	59	55
Range - yr	38–75	31–68
Sex
Male	11(55.0)	11(61.1)
Female	9(45.0)	7(38.9)
ECOG performance status
0	7(35.0)	–
1	11(55.0)	–
2	2(10.0)	–
R-ISS stage
I	3(15.0)	–
II	6(30.0)	–
III	11(55.0)	–
Cytogenetics
High risk	6(30.0)	–
Standard risk	12(60.0)	–
Unknown/missing	2(10.0)	–

Figure 1

Categories of tRFs/tiRNAs in MM patients and healthy donors. (A) The number of up-regulated tRFs/tiRNAs against tRNA isodecoders. (B) The number of down-regulated tRFs/tiRNAs against tRNA isodecoders. (C) The distribution of up-regulated tRFs/tiRNAs subtypes. (D) The distribution of down-regulated tRFs/tiRNAs subtypes.

Differentially Expressed tRFs/tiRNAs Between MM and Normal Groups

Venn diagram shows the commonly expressed and specifically expressed tRFs/tiRNAs (Figure 2A). The “commonly expression” means CPM values ≥ 20 in both two groups, and the “specifically expression” means CPM values ≥ 20 in one group and <20 in the other group. We identified 354 upregulated tRFs/tiRNAs and 515 downregulated tRFs/tiRNAs in this project in total. After limiting the condition with fold change ≥1.5 and P < 0.05, there remained 33 up-regulated and 22 down-regulated tRFs/tiRNAs. The heatmap shows hierarchical clustering of significantly differentially expressed tRFs/tiRNAs (Figure 2B). We present the overall tRFs/tiRNAs expression variation between the two groups in the volcano plots (Figure 2C).

Figure 2

Expression profiles of tRFs/tiRNAs in paired healthy donors and MM patients. (A) Venn diagram based on number of commonly expressed and specifically expressed tRFs/tiRNAs between healthy donors and MM patients. (B) The hierarchical clustering heatmap for significantly differentially expressed tRFs/tiRNAs. (C) The volcano plot of tRFs/tiRNAs (MM vs healthy donors).

tiRNA-1:34-Glu-TTC-2 Was Upregulated and tRF-60:76-Arg-ACG-1-M2 Was Downregulated in MM

We performed qPCR to detect the expression of three significantly upregulated or downregulated tRFs/tiRNAs in paired healthy donors and MM samples. The expression of all these six tRFs/tiRNAs was aligned to the sequencing data, while the most significantly upregulated and downregulated tRFs/tiRNAs were tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2, respectively (Figure 3A). Then, we expanded sample size for further verification of tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 expression in 18 healthy donors and 20 MM patients. At last, we confirmed that tiRNA-1:34-Glu-TTC-2 was upregulated and tRF-60:76-Arg-ACG-1-M2 was downregulated in MM (Figure 3B).

Figure 3

qPCR validation of tRFs/tiRNAs expression. (A) Relative expression of three significantly upregulated tRFs/tiRNAs and downregulated tRFs/tiRNAs in paired healthy donors and MM samples. (B) Relative expression of tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 18 healthy donors and 20 MM patients. **P<0.01, ***P<0.001.

Target gene prediction of tRFs/tiRNAs. (A) Target genes of tiRNA-1:34-Glu-TTC-2. (B) Target genes of tRF-60:76-Arg-ACG-1-M2.

Target Gene Prediction with Bioinformatics Tool

We used TargetScan algorithms to explore the putative roles of tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 in MM. Through this strategy, we predicted 193 conserved targets and 310 conserved targets, respectively. The target genes were presented in the network diagrams constructed by Cytoscape (Figure 4A and B).

Figure 4

Continue.

GO/KEGG Enrichment Analysis

We performed GO and KEGG enrichment analysis for target genes of tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 using DAVID database. GO analysis showed that 37 and 58 GO terms were enriched (P < 0.05) for target genes of tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2, respectively. The most enriched terms of target genes of tiRNA-1:34-Glu-TTC-2 were “transcription, DNA-templated” and “regulation of transcription, DNA-templated” in biological process (BP) category, “nucleus” and “intracellular” in cellular component (CC) category, “protein binding” and “metal ion binding” in molecular function (MF) category (Figure 5A). Whereas, the most enriched terms of target genes of tRF-60:76-Arg-ACG-1-M2 were “apoptotic process” and “regulation of transcription from RNA polymerase II promoter” in BP category, “cytoplasm” and “plasma membrane” in CC category, “protein binding” and “zinc ion binding” in MF category (Figure 5B).

Figure 5

Functional analysis of tRFs/tiRNAs. (A) GO enrichment analysis of tiRNA-1:34-Glu-TTC-2. (B) GO enrichment analysis of tRF-60:76-Arg-ACG-1-M2. (C) KEGG pathway analysis of tiRNA-1:34-Glu-TTC-2. (D) KEGG pathway analysis of tRF-60:76-Arg-ACG-1-M2.

KEGG pathway analysis indicated that tiRNA-1:34-Glu-TTC-2 might participate in cancer-related pathways such as Pathways in cancer, Wnt signaling pathway and Hippo signaling pathway (P < 0.05) (Figure 5C). While the most enriched pathways of tRF-60:76-Arg-ACG-1-M2 were cAMP signaling pathway, Endocytosis, Oxytocin signaling pathway and Rap1 signaling pathway (Figure 5D).

Discussion

Although the clinical outcome of MM has been considerably improved over the last several decades, MM is still an incurable type of hematological malignancy. The pathogenesis of MM is complex and multifactorial. Many genes or ncRNAs have been determined to be involved in the pathological process. However, whether tRFs/tiRNAs participate in MM initiation has not been elucidated. Therefore, we performed a preliminary evaluation of the role of tRFs/tiRNAs in MM.

First, we screened the expression profile of tRFs/tiRNAs in MM by RNA-sequencing. All tRFs/tiRNAs were named using the standardized tDR (tRNA-derived RNA) naming system (tDRnamer) from the Lowe Lab at the University of California Santa Cruz. The four parts in the name sequentially represent prefix, position, source tRNA and matching tRNA transcripts.20 tDRnamer provides a consistent, stable name and gives additional annotation for each submitted tDR sequence compared to tRFdb ID and MINTbase ID. The biogenesis of most of tRF-3s is unknown. The decrease in the number of tRF-3 may affect DNA replication, repair and damage response.10 Next, the two most differentially expressed tRFs/tiRNAs were verified by qPCR. To explore the possible biological functions of these RNAs, we predicted their target genes and performed GO/KEGG analysis. For tiRNA-1:34-Glu-TTC-2, its predicted target genes were primarily localized in the nucleus and implicated in transcription regulation and protein binding. However, the target genes of tRF-60:76-Arg-ACG-1-M2 were primarily located in the cytoplasm and participated in apoptotic processes and protein binding. This result is consistent with the common molecular functions of tRFs/tiRNAs in a previous report.12

tRFs/tiRNAs may regulate transcription and apoptotic processes through multiple molecular mechanisms. For example, Ivanov et al found that 5′-tiRNAAla and 5′-tiRNACys could repress global translation by displacing eIF4G and eIF4A (translation eukaryotic initiation factors) from mRNAs and displacing eIF4F from isolated m7G caps.21 Gkatza et al22 showed that the loss of NSUN2, a cytosine-5 RNA methyltransferase, could lead to impaired regulation of protein synthesis by altering the biogenesis of tiRNAs in response to stress.23 However, the mechanisms underlying how tRFs/tiRNAs regulate carcinogenesis have not been elucidated.

To explore whether tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 are involved in the initiation and development of MM, we performed KEGG pathway analysis. Among the top ten enriched pathways, pathways in cancer, the Wnt signaling pathway, the Hippo signaling pathway and the Rap1 signaling pathway are widely known to be related to carcinogenesis in previous reports.23–25 These pathways are also associated with MM initiation and development. Harmen et al reported that aberrant Wnt signaling mediated the proliferation, migration, and drug resistance of MM cells and MM cells secreted Wnt antagonists that promoted the development of osteolytic lesions by impairing osteoblast differentiation.26 The Hippo pathway has also emerged as an important mediator of MM oncogenesis and is implicated in the pathogenesis of MM-induced bone disease.27,28 There are relatively few studies investigating the role played by the Rap1 pathway in MM. Taken together, tiRNA-1:34-Glu-TTC-2 and tRF-60:76-Arg-ACG-1-M2 may be involved in the pathological process of MM but these hypotheses required further validation.

Several study limitations warrant consideration. First, due to the short follow-up time, we have not obtained the efficacy and prognosis data of myeloma patients to date. Second, in this preliminary experiment, we did not collect specimens from cases of different stages of plasma cell disease (monoclonal gammopathy of undetermined significance, smoldering myeloma, multiple myeloma and plasma cell leukemia). The evaluation of the expression of tRFs/tiRNAs in these four stages requires further study. Third, although tRFs/tiRNAs are also involved in the silencing of microRNAs, there are few related data reports. The current bioinformatics databases do not provide information on the network of tRFs/tiRNAs-microRNA-mRNA. We look forward to realize this analysis basing on the in-depth research on tRFs/tiRNAs in the future.

In conclusion, we assessed the potential function of tRFs/tiRNAs in MM and the results establish a foundation for continued research. tRFs/tiRNAs may play important roles in the pathogenesis of myeloma and serve as clinical biomarkers and therapeutic targets in the future. Further experiments are warranted to explore how tRFs/tiRNAs regulate target genes and miRNAs in MM. We look forward to further research on tRFs/tiRNAs in MM.