Identification of aberrantly expressed lncRNAs and ceRNA networks in multiple myeloma: a combined high-throughput sequencing and microarray analysis

Background: This study aimed to explore the potential effects of long non-coding RNAs (lncRNAs) in multiple myeloma (MM) patients using two detection methods: high-throughput sequencing and microarray.

Methods: In this study, lncRNAs were detected in 20 newly diagnosed MM patients, with 10 patients analyzed by whole transcriptome-specific RNA sequencing and 10 patients analyzed by microarray (Affymetrix Human Clariom D). The expression levels of lncRNAs, microRNAs, and messenger RNAs (mRNAs) were analyzed, and the differentially expressed lncRNAs identified by both methods were selected. The significant differentially expressed lncRNAs were further validated using PCR.

Results: This study established the aberrant expression of certain lncRNAs involved in the occurrence of MM, with AC007278.2 and FAM157C showing the most significant differences. The top 5 common pathways identified by the Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis were the chemokine signaling pathway, inflammatory mediator regulation, Th17 cell differentiation, apoptosis, and NF-kappa B signaling pathway. Furthermore, three microRNAs (miRNAs) (miR-4772-3p, miR-617, and miR-618) were found to constitute competing endogenous RNA (ceRNA) networks in both sequencing and microarray analyses.

Conclusions: By the combination analysis, our understanding of lncRNAs in MM will be increased significantly. More overlapping differentially expressed lncRNAs were found to predict therapeutic targets precisely.

1 Introduction

Multiple myeloma (MM) is a plasma cell neoplasm characterized by the clonal proliferation of malignant plasma cells in the bone marrow microenvironment (1). Proteasome inhibitors, immunomodulatory drugs, and CD38-targeting antibodies have been extending survival (2, 3). Despite the recent advances in clinical and experimental oncology, many patients will relapse and die from MM. Disease progression and subsequent relapses are characterized by increasingly resistant disease and subclonal evolution. Risk stratification based on disease-specific is important for prognosis and the treatment strategy (4). Many current research strategies were used to guide the therapy and refine the immunotherapeutic approaches (5, 6).

Accurate prognostic evaluation and therapeutic targets of multiple myeloma depend on gene detection. Thus, a detailed understanding of the mechanisms is very important. Most long non-coding RNAs (lncRNAs) regulate tumor proto and suppressor oncogene pathways, as transcribers record regulatory factors (7, 8). In addition, epigenetic modification and interaction with RNA-binding proteins had also been confirmed to be the main mechanisms and functions of lncRNAs (9). Through a competing endogenous RNA (ceRNA) network, lncRNAs regulate the function of microRNA (miRNA) as a miRNA sponge (10, 11). LncRNAs are involved in the occurrence and development of MM. The action modes of different lncRNAs are very different, and upregulated or downregulated of lncRNAs can lead to various carcinogenic tendencies in plasma cells (12, 13). The aberrant expression of lncRNAs in more aggressive stages of MM has been found (14). Ronchetti et al. (14) provided a special view of the expression of lncRNA in MM. In previous studies, through high-throughput sequencing, some lncRNAs were found abnormally expressed in newly diagnosed MM, including maternally expressed gene 3 (MEG3), colon cancer-associated transcript 1 (CCAT1), and coiled-coil domain-containing 26 (CCDC26) (15–17). In addition to known lncRNAs, some unpublished lncRNAs were found, and new RNAs related to MM pathogenic genes were identified, providing new targeted therapies for MM.

RNA-Seq can be applied to investigate different populations of RNA, including messenger RNA (mRNA), circular RNA (circRNA), and lncRNA. High-throughput sequencing technology can be called a milestone in the development of sequencing technology. In addition to high-throughput sequencing, gene microarray has also become an important method of gene detection, and gene chip has the analytical ability to fully display the genetic spectrum of target samples through the mechanism of hybridization (18). Compared with sequencing, gene microarray is not possible to discover novel genes, but it will focus on lncRNAs that have been found to play a role in disease. Affymetrix Human Clariom D Array is a new generation of gene microarray chips, which detects mRNA, lncRNA, miRNA, and circRNA information of known sequences. Whether sequencing or gene microarray, software tools can perform hierarchical clustering, Gene Ontology (GO) enrichment analysis, pathway enrichment analysis, co-expression network (lncRNA–circRNA–mRNA), lncRNA target pathway network, and ceRNA network (lncRNA–miRNA–mRNA or circRNA–miRNA–mRNA). Through the analysis of this information, many unknown functions of lncRNA in MM can be inferred.

In our study, we conducted a comparative analysis of high-throughput sequencing data and gene chip data for lncRNAs. The aim was to identify overlapping differentially expressed genes, thereby allowing for the identification of more precise target genes that are significantly associated with the occurrence and progression of MM.

2 Materials and methods

2.1 Patients and cell samples of bone marrow

The plasma cells of bone marrow were collected from 20 MM patients, and normal individuals were taken as the control. The patients met the updated diagnostic criteria of the International Myeloma Working Group (IMWG) for MM (19). Inclusion criteria: 1) p reviously untreated patients diagnosed with multiple myeloma according to the IMWG criteria. 2) According to IMWG standard disease, there are evaluable blood and urine M protein indexes (serum m-protein ≥1 g/dl (≥10 g/L) and serum m-protein ≥0.5 g/dl (≥5 g/L) in IgA, IgD, IgE, or IgM multiple myeloma subjects; or urine m-protein levels ≥200 mg/24 h; or light chain, no measurable lesions in serum or urine, serum immunoglobulin free light chain (FLC) ≥10 mg/dl and serum immunoglobulin Kappa Lambda FLC ratio abnormalities). There may be a measurable extramedullary mass. Exclusion criteria: 1) p atients did not meet the MM diagnostic criteria updated by the International Myeloma Working Group, 2) p atients with relapsed and refractory multiple myeloma, 3) MM combined with other tumors, 4) MM prepared for autologous stem cell transplantation, and 5) s olitary plasmacytoma. This study met the standards of medical ethics and was approved by the hospital ethics committee. All participants provided signed informed consent, and this study was approved by the institutional ethical review board of the Chao-Yang Hospital, Capital Medical University (ethical approval code 2016-science-84, 2019-science-118). The median age of the patients was 63 (53–74) in the sequencing group and 61 (47-82) in the chip group.

Mononuclear cells were extracted from 5 ml of bone marrow. Plasma cells were sorted by CD138 magnetic beads. Total RNA was extracted using a TRIzol reagent. Non-coding RNAs of 10 patients were detected by sequencing and those of 10 patients by gene chip. Total RNA was purified using miRNeasy Serum/Plasma kit (Qiagen, Hilden, Germany).

2.2 RNA sequencing and identification of lncRNAs by high-throughput sequencing

Raw sequencing data were filtered out. The clean data were obtained after filtering. To ensure the quality of the data used for information analysis, the raw data (raw reads) obtained from sequencing were initially filtered using in-house programs. The filtering conditions included the following: 1) removal of reads contaminated with adapter sequences (reads with adapter contamination > 5 bp were removed; for paired-end sequencing, if one end was contaminated, both ends of the read were removed); 2) removal of low-quality reads (reads with base quality values of Q ≤ 19 accounting for > 15% of the total bases; for paired-end sequencing, if one end was of low quality, both ends of the read were removed); 3) removal of reads with an N ratio > 5% (for paired-end sequencing, if one end had an N ratio > 5%, both ends of the read were removed). The filtered data, referred to as clean data, were subjected to quality and quantity analyses, including Q30 calculation, data volume calculation, and base composition analysis. Reference genomes were downloaded from Ensembl (http://www.ensembl.org/index.html), and the clean data were mapped to the reference genome through HISAT2 (http://ccb.jhu.edu/software/hisat2/index.shtml) (20–23). FPKM (Fragments Per Kilobase Million Mapped Reads) were calculated for the gene expression of all samples. DESeq (http://www.bioconductor.org/packages/release/bioc/html/DESeq.html) was used for differential expression analysis (21, 22). The mRNAs with Spearman’s correlation coefficient (p ≥ 0.9) were selected as the trans-targets, and mRNAs with a distance of less than 50 kb were selected as the cis-targets. To identify potential lncRNAs from the sequencing data, a series of steps were employed. First, transcripts with coding potential were filtered out based on transcript length, exon number, and other relevant criteria. Transcripts with low read coverage were also eliminated across all samples, as well as known coding RNAs and non-coding RNAs such as rRNAs, tRNAs, snoRNAs, and snRNAs. Next, the class_code information from gffcompare was used to categorize the remaining transcripts into lincRNA, intronic lncRNA, or anti-sense lncRNA. Finally, a comprehensive analysis was performed using various coding potential analysis tools, including CNCI, CPC, PFAM protein domain analysis, and CPAT, to confirm the non-coding nature of the remaining transcripts. These transcripts were then identified as novel lncRNAs.

2.3 Identification of DE lncRNAs for the microarray

The human Clariom™ D array (Affymetrix, Santa Clara, CA, USA) was used to analyze the expression of lncRNAs, miRNAs, and mRNAs. The ceRNA network was generated using GeneChip Command Console Software (AGCC version 4.0; Affymetrix) (24). One-Cycle Target Labeling and Control Reagents (Affymetrix, Santa Clara, CA, USA) were used for cDNA synthesis, and cRNA was generated using GeneChip^® IVT Labeling Kit (Affymetrix). cRNA was fragmented and then hybridized to Affymetrix Human Clariom D array. After being washed and stained in the Affymetrix Fluidics Station 450, GeneChips were scanned by using Affymetrix GeneChip Command Console (AGCC version 4.0) installed in GeneChip^® Scanner 3000 7G. The original data (CEL.files) were analyzed by Robust Multichip Analysis (RMA) algorithm. Summarizing Affymetrix default analysis settings and global scaling as a quantile normalization method, Log2 RMA signal intensity was obtained.

Hierarchical clustering is a powerful and visualized tool for analyzing the status of all genomic datasets. The matrix distance between the gene expression data was calculated by Cluster_Treeview software (Palo Alto, CA, USA). The statistical method was a t-test for the mean of two samples, the differentially expressed genes were selected when the p- value < 0.05 and the threshold of fold change (FC) value ≥2.0. The heatmap of all differential genes and all samples was illustrated in visualization by the Tree view (25, 26).

2.4 Enrichment analysis

GO (http://geneontology.org/) provides a standard vocabulary of gene function, including molecular function, cell component, and biological process. Kyoto Encyclopedia of Genes and Genomes (KEGG; http://www.kegg.jp/) contains the molecular interaction and reaction networks of the genes. Functional enrichment analysis based on hypergeometric distribution, the significantly enriched GO terms, and KEGG pathways are identified.

2.5 CeRNA analysis

CeRNA is a deep analysis of mRNA, lncRNA, circRNA, and microRNA (27, 28). Target genes of the miRNAs were predicted by MiRanda, PITA, and TargetScan, and the target genes were the cumulative outcome of the prediction using at least two programs.

Pearson’s correlation coefficient (PCC) was calculated according to the expression of mRNA and lncRNA. The number and binding score of miRNA response elements (MREs) can evaluate the ability of ceRNA molecules (27, 28).

2.6 Quantitative real-time reverse transcriptase–PCR

The expression of two downregulated overlapping differentially expressed (DE) lncRNAs in MM patients was examined using high-throughput sequencing and microarray using the quantitative reverse transcriptase– polymerase chain reaction (RT-PCR). Total RNA was isolated by RNA Extraction Kit DNase I (GenePool, Beijing, China; Cat# GPQ1801). The cDNA was synthesized using lncRNA cDNA Synthesis Kit with GenePool (Cat# GPQ1806), and quantitative PCR was performed on the Line Gene 9600 Plus (Bioer Technology, Zhejiang, China). Briefly, PCRs were performed in a triplicate mixture (20 µl) containing dNTP Mix 2.5 mM each 4 µl, Primer Mix 2 µl, RNA Template 2 µl, 5× RT Buffer 4 µl, DTT 0.1M 2 μl, HiFiScript 200 U/μl, RNase-Free Water up to 20 μl. GADPH was used as an internal control. After amplification, a melting curve analysis was performed to confirm reaction specificity. The amplification procedure was as follows: 95°C for 35 s, (95°C for 5 s, 60°C for 30 s) × 45 cycles, and conduct melting curve analysis simultaneously at 60-95°C. The upstream primer of AC007278.2 is CCATTCATCCACAACGCAAGG, and the downstream primer is TGGCAGACACAGAGTATCTTCAC. The upstream primer of FAM157C is TCAGGACCTACTCGGGAGAC, and the downstream primer is CTGAGGAGCTGAGGAGAAGGA.

3 Results

3.1 Clinical characteristics of patients with MM in the two groups

The clinical and genetic characteristics of the 20 patients with MM in the two groups were recorded and shown in Table 1. The high-risk MM patients in the sequencing group and gene chip group were 70% and 80%, respectively.

TABLE 1

Table 1 Clinical characteristics of patients with MM.

3.2 DE lncRNAs in sequencing group

Compared to the control group, MM cells exhibited a significant difference in the expression of lncRNAs, with a total of 4,419 DE lncRNAs identified. Among these DE lncRNAs, 1,559 were known lncRNAs that showed differential expression, with 910 upregulated and 649 downregulated. The volcanic map of differentially expressed known lncRNAs in the sequencing group is shown in Figure 1. The clustering results of the known lncRNAs are shown in Figure 2. The top 20 DE known lncRNAs are shown in Table 2.

FIGURE 1

Figure 1 Volcanic map of differentially expressed known lncRNAs in the sequencing group: the abscissa is the expression of multiple changes in different experimental groups, the ordinate is the statistical significance of the expression change, red and green respectively represent differentially expressed upregulated and downregulated genes, and gray represents non- differentially expressed genes (910 upregulated and 649 downregulated). Different colors represent different classifications. lncRNAs, long non-coding RNAs.

FIGURE 2

Figure 2 Heatmaps of DE lncRNAs in the sequencing group. The various colors represent different expression levels. Blue color indicates the lower expression level, and yellow color indicates higher expression levels. DE lncRNAs, differentially expressed long non-coding RNAs.

TABLE 2

Table 2 Analysis of top 20 known DE lncRNAs in sequencing group.

3.3 D ifferentially expressed lncRNAs in gene microarray group

A total of 1,263 DE lncRNAs were identified in microarray (748 upregulated and 515 downregulated). The clustering results of the lncRNAs of the samples in the gene microarray group according to the expression of differentially expressed genes are shown in Figure 3.

FIGURE 3

Figure 3 Heatmaps of DE lncRNAs in gene microarray group. Each column represents specimens with different serial numbers. C represents the control group, and MM represents the experimental group. The red box indicates the upregulated lncRNA molecule. The blue box indicates the down regulated lncRNA molecule. DE lncRNAs, differentially expressed long non-coding RNAs.

3.4 Overlapping DE lncRNA in sequencing group and gene microarray group

A total of 64 overlapping DE lncRNAs were detected by the two methods (26 genes were upregulated, as shown in Table 3, and 38 were downregulated, as shown in Table 4). The correlation between the two methods is shown in Figure 4 (R = 0.698, p < 0.001). Among the 26 upregulated genes, IFNG-AS1, AC093818.1, AF127936.5, and OLMALINC were expressed in different diseases, especially IFNG-AS1, which has been reported in many malignant tumors (29–31). Among the down regulated genes, DLEU2, LINC01270, LINC01270, GK-IT1, LINC01001, DLEU2, CHRM3-AS2, AC007278.2, PLBD1-AS1, FAM157C, KCNJ2-AS1, and LINC01127 were abnormally expressed in other diseases in previous studies (32, 33). AC007278.2 and FAM157C were identified as the most significantly differentially expressed lncRNAs between the MM cells and the control group.

TABLE 3

Table 3 Overlapping upregulated DE lncRNAs in sequencing group and gene microarray group.

TABLE 4

Table 4 Overlapping downregulated DE lncRNAs in sequencing group and gene microarray group.

FIGURE 4

Figure 4 A scatter plot of the expression of the microarray group versus the sequencing group of the common molecules. The horizontal coordinates represent the fold-change expression of sequencing, and the vertical coordinates represent the fold-change expression of microarray. The dots indicate the common molecules between differentially expressed (p < 0.05) molecules based on sequencing and microarray. The line indicates regression line. p-Value and R-value were calculated by Spearman’s correlation.

3.5 Functional prediction of DE lncRNAs in sequencing

According to co-location and co-expression of DE lncRNAs, the most significant enrichment in biological process was due to biological regulation, regulation of cellular process, and regulation of biological process; the most significant enrichment in cellular component was due to membrane, membrane part, and vesicle. Protein binding and activity were most significantly enriched in molecular function. The results are shown in Figures 5A–C.

FIGURE 5

Figure 5 GO analysis of lncRNAs in the sequencing group. (A) Enrichment of biological process. (B) Enrichment of molecular function. (C) Cellular component changed mRNA-targeted ncRNAs in MM. GO, Gene Ontology; lncRNAs, long non-coding RNAs; mRNA, messenger RNA; ncRNAs, non-coding RNAs; MM, multiple myeloma.

The KEGG pathway is shown in Figure 6; the 15 most significant enrichment pathways included osteoclast differentiation, chemokine signaling pathway, alcoholism, malaria, systemic lupus erythematosus, leishmaniasis, leukocyte transendothelial migration, inflammatory mediator regulation of TRP channels, amoebiasis, relaxin signaling pathway, viral carcinogenesis, phagosome, pertussis, Th17 cell differentiation, and apoptosis (Figure 6).

FIGURE 6

Figure 6 The KEGG pathway of lncRNAs in sequencing group: KEGG pathway scatterplots in the plasma cells of MM. The KO enrichment of all samples was combined, and the distribution map was drawn according to the significant Q value of the enrichment of the samples. Each point indicates the degree of enrichment of the KO entry, and the closer the color approaches red, the higher the degree of enrichment. The size of each point indicates the number of genes enriched in the KO entry. The larger the point, the more the genes are enriched in the GO entry, and vice versa. KEGG, Kyoto Encyclopedia of Genes and Genomes; lncRNAs, long non-coding RNAs; MM, multiple myeloma; KO, KEGG Orthology.

3.6 Functional prediction of DE lncRNAs in microarray

The upregulated gene functions mainly focus on protein synthesis, most of which have been confirmed to be related to cancer replication, such as rRNA processing, protein translation, and SRP- dependent co- translation proteins targeting membranes. Inhibiting neutrophil degranulation, inflammatory response, immune response, and other immune functions of killing myeloma cells while reducing hematopoietic capacity and reducing cell apoptosis is consistent with the pathogenesis of MM. KEGG analysis had similar results in gene pathway enrichment, as shown in Figures 7A, B.

FIGURE 7

Figure 7 KEGG enrichment analysis of differentially expressed in microarray group. (A) Upregulated regulatory pathways. (B) Downregulated regulatory pathways. KEGG, Kyoto Encyclopedia of Genes and Genomes.

3.7 Common functional prediction of DE lncRNAs in sequencing and microarray

Comparing the two lncRNA detection methods, it was found that the common pathway in the top 5 KEGG was the chemokine signaling pathway, inflammatory mediator regulation, Th17 cell differentiation, apoptosis, and NF-kappa B signaling pathway. After further analysis, the differential lncRNAs of the corresponding pathway were found. Among them, there are 20 differential genes in the inflammatory mediator regulation pathway, 18 differential genes in the Th17 cell differentiation pathway, 29 differential genes in the apoptosis pathway, and 16 differential genes in the NF- kappa B signaling pathway. The corresponding path of AC007278.2 is cytokine receptor interaction in both the sequencing group and gene microarray group.

The targeted pathway of AC007278.2 in both groups was Cytokine –Cytokine receptor interaction. AC007278.2 (sense_intronic Position:chr2:102433957-102435340:+) is Toll/interleukin-1 receptor (TIR) domain (34). The TIR homology domain is an intracellular signaling domain found in MyD88, interleukin- 1 receptor, and the Toll receptor. LncRNA AC007278.2 promotes T helper cell 1 differentiation by regulating NF-κB import into the nucleus and activating interleukin-1 receptor and interleukin-18 receptor. The KEGG targeted pathway of AC007278.2 is interleukin- 18 receptor 1 and IL1RRP. The KEGG targeted pathway of FAM157C (cRNA Position:chr16:90102271-90186204:+) has potential kinase structure characteristics; the KEGG targeted pathway of FAM157C is adenylate cyclase 4.

3.8 qPCR verification of lncRNA

The top 20 DE lncRNAs were verified by overlapping DE lncRNAs, including AC007278. 2 and FAM157C, especially AC007278.2, which was the most differentially expressed lncRNAs in the sequencing group. The qPCR results were consistent with the sequencing results (Figures 8, 9). AC007278. 2 and FAM157C were downregulated in both the sequencing group and the chip group.

FIGURE 8

Figure 8 The qPCR results of AC007278.2: AC007278.2 was downregulated in sequencing group and microarray group compared with the control group. The ordinate is the relative expression of AC007278. 2.

FIGURE 9

Figure 9 The qPCR results of FAM157C: FAM157C was downregulated in sequencing group and microarray group compared with the control group. The ordinate is the relative expression of FAM157C.

3.9 Analysis of the regulatory network of lncRNAs (ceRNAs)

A total of 19 known lncRNAs and 116 novel lncRNAs were speculated to constitute ceRNA with miRNAs and mRNA in the sequencing group; the proportion of novel lncRNAs was significantly higher than that of known lncRNAs. Only a small number of RNAs were shown to be involved in the ceRNA network in the microarray group, including 9 lncRNAs, 8 miRNAs, and 51 mRNAs (Figure 10).

FIGURE 10

Figure 10 CeRNA network in microarray group: gene molecules with ceRNA relationship are shown in this network. Square, triangle, and diamond represent lncRNA, miRNA, and mRNA, respectively. ceRNA, competing endogenous RNA; lncRNA, long non-coding RNA; miRNA, microRNA; mRNA, messenger RNA.

Among these genes, three miRNAs constituted ceRNA of lncRNA–miRNA–mRNA in both the sequencing group and microarray group, including miR-4772-3p, miR-617, and miR-618. MiR-4772-3p was upregulated in malignant plasma cells, which constituted ceRNA with ZNF850 and novel lncRNA-MSTRG.60361 in the sequencing group. MiR-4772-3p was speculated to constitute ceRNA with BMS1P5, RPL10AP6, RPL4P4, RPSAP19, RPL12P4, and RPS13P2 in microarray. MiR-617 constituted ceRNA with novel lncRNA-MSTRG.255209, MSTRG.42167, and MSTRG.42200 in the sequencing group, while miR-617 constituted ceRNA with BMS1P5 in the microarray. MiR-618 constituted ceRNA with HIC1, MSTRG.42140, and MSTRG.42235 in the sequencing group, the pathway dominated by BTB/POZ. MiR-618 was speculated to constitute ceRNA with BLOC1S5-TXNDC5 and BMS1P5 in the microarray group.

4 Discussion

MM is the second most common hematological malignancy. Recently, microarrays and sequencing have become important tools for identifying treatment targets and developing personalized treatments. In our study, we detected overlapping differential expression of lncRNAs and identified their common functions through high-throughput sequencing and genetic microarray. We also investigated the ceRNA of multiple myeloma and selected two significant overlapping differentially expressed lncRNAs that could be used as potential therapeutic targets for MM.

The study found that lncRNAs were significantly dysregulated in MM cells. A total of 1,263 known lncRNAs were detected to be significantly regulated in malignant plasma cells of newly diagnosed MM (NDMM) by microarray assay, and more DE lncRNAs were found by high- throughput sequencing, in addition to 1,559 known lncRNAs. After comparison, the two-study mode has a common finding. Among them, 64 overlapping lncRNAs were found to be differentially expressed in both the sequencing group and the microarray group. Some lncRNAs up regulated in different diseases, such as IFNG-AS1 (35), AC093818.1, AF127936.5, and OLMALINC.IFNG-AS1, have been reported in many malignant tumors. Among the down regulated lncRNAs, DLEU2, LINC01270, LINC01270, GK-IT1, LINC01001, DLEU2, CHRM3-AS2, AC007278.2, PLBD1-AS1, FAM157C, KCNJ2-AS1, and LINC01127 were abnormally expressed in other diseases in previous studies. Among them, AC007278.2 and FAM157C showed more significant differences in the two groups, which may be important regulatory factors for the development of multiple myeloma. Especially, AC007278.2 is the top 1 DE known lncRNA in the sequencing group, which may become a potential therapeutic target of multiple myeloma.

The study found that lncRNAs were significantly dysregulated in MM cells. Microarray assay detected 1,263 known lncRNAs that were significantly regulated in malignant plasma cells of NDMM, and more DE lncRNAs were found by high- throughput sequencing, including 1,559 known lncRNAs. Comparison between the two methods revealed common findings, with 64 overlapping lncRNAs identified as differentially expressed in both sequencing and microarray groups. Some lncRNAs were found to be upregulated in different diseases, such as IFNG-AS1, AC093818.1, AF127936.5, and OLMALINC; others were downregulated in various diseases, such as DLEU2, LINC01270, GK-IT1, LINC01001, CHRM3-AS2, AC007278.2, PLBD1-AS1, FAM157C, and KCNJ2-AS1. Among them, AC007278.2 and FAM157C showed more significant differences in the two groups, indicating their potential role as important regulatory factors in the development of MM. Specifically, AC007278.2 was identified as the top DE known lncRNA in the sequencing group and may represent a potential therapeutic target for MM.

The comparison of the function of lncRNAs using two detection methods revealed that the top 5 KEGG pathways with commonly dysregulated lncRNAs were the chemokine signaling pathway, inflammatory mediator regulation, Th17 cell differentiation, apoptosis, and NF-kappa B signaling pathway. Moreover, lncRNAs were found to be involved in the regulation of apoptosis pathways and immune function- related pathways. The identification of these pathways and related lncRNAs provides a clearer understanding of the role of lncRNAs in MM. Further in vivo and in vitro validations will be conducted in future studies to confirm these findings.

Comparing ceRNA networks, lncRNA–miRNA–mRNA networks were speculated from the two methods. The analysis revealed that miR-4772-3p, miR-617, and miR-618 were identified as ceRNAs with lncRNAs in malignant plasma cells in both groups. Although miR-4772-3p has been less investigated in previous research on malignant tumors, it was speculated as an important miRNA of ceRNA in both detection methods, requiring further research for confirmation. Our findings also suggest that miR-618 may be involved in MM mechanisms. Previously, miR-618 was shown to play a role in tumor inhibition by targeting XIAP expression (36). In our study, miR-618 constituted ceRNA with HIC1, MSTRG.42140, and MSTRG.42235 in the sequencing group, with the pathway being dominated by BTB/POZ.

The study revealed common differentially expressed lncRNAs, functions, and ceRNA in newly diagnosed MM provided by whole-genome sequencing and whole-genome microarray. Gene chip is a closed detection system that can only detect known gene information, while sequencing is an open detection system that can find new gene data. After nearly 20 years of development, the chip platform has greatly improved quality control standards and analysis mode, enabling faster results. However, sequencing is the mainstream technology, and it is powerful in completing massive data analysis and finding more new genes. The results obtained by the two detection methods can mutually verify the expression and function of lncRNAs in MM. Co-expression of DE lncRNAs was identified, which is associated with the pathogenesis of MM, contributing to the further exploration of different biological characteristics in MM.

In conclusion, both sequencing and microarray are powerful methods for generating genome-wide information, and our comparative study suggests that the two methods produce similar association profiles, enabling the integration of both data sets. However, a limitation of our study is the limited number of relevant studies currently available. A more comprehensive study with larger sample size is needed to clarify the expression of non-coding RNA in MM and validate potential target genes, as well as to explore the underlying mechanisms.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

Ethics statement

All participants signed the informed consent, and this study was approved by institutional ethical review board of the Chao-Yang Hospital, Capital Medical University (Ethical approval code 2016-science-84, 2019-science-118).

Author contributions

Conceptualization: M-QL and W-MC. Data curation: YJ. Formal analysis: WG. Funding acquisition: W-MC. Project administration: LB. Resources: H-XZ. Software: YW. Supervision: LB. Validation: H-XZ. Writing—original draft: Y-QH. All authors contributed to the article and approved the submitted version.

Funding

This study was supported by the Yangfan Project Special Foundation of Beijing Hospital Authority (ZXLY201606), Beijing JST Research Funding (IR-202105), and Special Medicine Innovation Scientific Special Project (2018ZX09733003).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2023.1160342/full#supplementary-material

References

1. Sanoja-Flores L, Flores-Montero J, Garces JJ, Paiva B, Puig N, Garcia-Mateo A, et al. Next generation flow for minimally-invasive blood characterization of mgus and multiple myeloma at diagnosis based on circulating tumor plasma cells (Ctpc). Blood Cancer J (2018) 8(12):117. doi: 10.1038/s41408-018-0153-9

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Gandhi UH, Cornell RF, Lakshman A, Gahvari ZJ, McGehee E, Jagosky MH, et al. Outcomes of patients with multiple myeloma refractory to Cd38-targeted monoclonal antibody therapy. Leukemia (2019) 33(9):2266–75. doi: 10.1038/s41375-019-0435-7

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Franssen LE, Stege CAM, Zweegman S, van de Donk N, Nijhof IS. Resistance mechanisms towards Cd38-directed antibody therapy in multiple myeloma. J Clin Med (2020) 9(4):1195–1214. doi: 10.3390/jcm9041195

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Fotiou D, Gavriatopoulou M, Terpos E. Multiple myeloma and thrombosis: prophylaxis and risk prediction tools. Cancers (Basel) (2020) 12(1):191–207. doi: 10.3390/cancers12010191

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Ferguson ID, Patiño-Escobar B, Tuomivaara ST, Lin YT, Nix MA, Leung KK, et al. The surfaceome of multiple myeloma cells suggests potential immunotherapeutic strategies and protein markers of drug resistance. Nat Commun (2022) 13(1):4121. doi: 10.1038/s41467-022-31810-6

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Bal S, Schmidt TM, Costa LJ, Callander NS. Clinical implications of measurable residual disease assessment in multiple myeloma in the era of quadruplet therapy. Leuk Lymphoma (2022) 14:141–11. doi: 10.1080/10428194.2022.2123231

CrossRef Full Text | Google Scholar

7. Carrasco-Leon A, Ezponda T, Meydan C, Valcárcel LV, Ordoñez R, Kulis M, et al. Characterization of complete lncrnas transcriptome reveals the functional and clinical impact of lncrnas in multiple myeloma. Leukemia (2021) 35(5):1438–50. doi: 10.1038/s41375-021-01147-y

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Drak Alsibai K, Meseure D. Tumor microenvironment and noncoding rnas as Co-drivers of epithelial-mesenchymal transition and cancer metastasis. Dev Dyn (2018) 247(3):405–31. doi: 10.1002/dvdy.24548

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Guttman M, Donaghey J, Carey BW, Garber M, Grenier JK, Munson G, et al. Lincrnas act in the circuitry controlling pluripotency and differentiation. Nature (2011) 477(7364):295–300. doi: 10.1038/nature10398

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Shen X, Kong S, Yang Q, Yin Q, Cong H, Wang X, et al. Pcat-1 promotes cell growth by sponging mir-129 Via Map3k7/Nf-Kb pathway in multiple myeloma. J Cell Mol Med (2020) 24(6):3492–503. doi: 10.1111/jcmm.15035

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Karreth FA, Reschke M, Ruocco A, Ng C, Chapuy B, Léopold V, et al. The braf pseudogene functions as a competitive endogenous rna and induces lymphoma in vivo. Cell (2015) 161(2):319–32. doi: 10.1016/j.cell.2015.02.043

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Butova R, Vychytilova-Faltejskova P, Souckova A, Sevcikova S, Hajek R. Long non-coding rnas in multiple myeloma. Noncoding RNA (2019) 5(1):13–27. doi: 10.3390/ncrna5010013

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Pan Y, Zhang Y, Liu W, Huang Y, Shen X, Jing R, et al. Lncrna H19 overexpression induces bortezomib resistance in multiple myeloma by targeting mcl-1 Via mir-29b-3p. Cell Death Dis (2019) 10(2):106. doi: 10.1038/s41419-018-1219-0

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Ronchetti D, Agnelli L, Taiana E, Galletti S, Manzoni M, Todoerti K, et al. Distinct lncrna transcriptional fingerprints characterize progressive stages of multiple myeloma. Oncotarget (2016) 7(12):14814–30. doi: 10.18632/oncotarget.7442

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Liu Z, Chen Q, Hann SS. The functions and oncogenic roles of Ccat1 in human cancer. BioMed Pharmacother (2019) 115:108943. doi: 10.1016/j.biopha.2019.108943

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Lu M, Wu Y, Gao W, Tian Y, Wang G, Liu A, et al. Novel non-coding rna analysis in multiple myeloma identified through high-throughput sequencing. Front Genet (2021) 12:625019. doi: 10.3389/fgene.2021.625019

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Ghafouri-Fard S, Esmaeili M, Taheri M. Expression of non-coding rnas in hematological malignancies. Eur J Pharmacol (2020) 875:172976. doi: 10.1016/j.ejphar.2020.172976

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Zhang H, Zhang Q, Liao Z. Microarray data analysis of molecular mechanism associated with stroke progression. J Mol Neurosci (2019) 67(3):424–33. doi: 10.1007/s12031-018-1247-3

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Rajkumar SV, Dimopoulos MA, Palumbo A, Blade J, Merlini G, Mateos MV, et al. International myeloma working group updated criteria for the diagnosis of multiple myeloma. Lancet Oncol (2014) 15(12):e538–48. doi: 10.1016/s1470-2045(14)70442-5

PubMed Abstract | CrossRef Full Text | Google Scholar

20. IEEE international symposium on biomedical imaging. IEEE Pulse (2017) 8(5):65. doi: 10.1109/mpul.2017.2746243

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Kim D, Langmead B, Salzberg SL. Hisat: a fast spliced aligner with low memory requirements. Nat Methods (2015) 12(4):357–60. doi: 10.1038/nmeth.3317

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Langmead B, Salzberg SL. Fast gapped-read alignment with bowtie 2. Nat Methods (2012) 9(4):357–9. doi: 10.1038/nmeth.1923

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by rna-seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol (2010) 28(5):511–5. doi: 10.1038/nbt.1621

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Garzon R, Marcucci G, Croce CM. Targeting micrornas in cancer: rationale, strategies and challenges. Nat Rev Drug Discovery (2010) 9(10):775–89. doi: 10.1038/nrd3179

CrossRef Full Text | Google Scholar

25. Sioson AA, Mane SP, Li P, Sha W, Heath LS, Bohnert HJ, et al. The statistics of identifying differentially expressed genes in expresso and Tm4: a comparison. BMC Bioinf (2006) 7:215. doi: 10.1186/1471-2105-7-215

CrossRef Full Text | Google Scholar

26. Schlicker A, Domingues FS, Rahnenführer J, Lengauer T. A new measure for functional similarity of gene products based on gene ontology. BMC Bioinf (2006) 7:302. doi: 10.1186/1471-2105-7-302

CrossRef Full Text | Google Scholar

27. Hornakova A, List M, Vreeken J, Schulz MH. Jami: fast computation of conditional mutual information for cerna network analysis. Bioinformatics (2018) 34(17):3050–1. doi: 10.1093/bioinformatics/bty221

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Zhu Y, Bian Y, Zhang Q, Hu J, Li L, Yang M, et al. Construction and analysis of dysregulated lncrna-associated cerna network in colorectal cancer. J Cell Biochem (2019) 120(6):9250–63. doi: 10.1002/jcb.28201

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Padua D, Mahurkar-Joshi S, Law IK, Polytarchou C, Vu JP, Pisegna JR, et al. A long noncoding rna signature for ulcerative colitis identifies ifng-As1 as an enhancer of inflammation. Am J Physiol Gastrointest Liver Physiol (2016) 311(3):G446–57. doi: 10.1152/ajpgi.00212.2016

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Wang Z, Cao Z, Wang Z. Significance of long non-coding rna ifng-As1 in the progression and clinical prognosis in colon adenocarcinoma. Bioengineered (2021) 12(2):11342–50. doi: 10.1080/21655979.2021.2003944

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Jafari L, Izadirad M, Vatanmakanian M, Ghaedi H, Farsiani MA, Mohammadi MH, et al. Ifng-As1 and Maf4 long non-coding rnas are upregulated in acute leukemia patients who underwent bone marrow transplantation. Curr Res Transl Med (2021) 69(4):103307. doi: 10.1016/j.retram.2021.103307

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Wang S, Wang E, Chen Q, Yang Y, Xu L, Zhang X, et al. Uncovering potential lncrnas and mrnas in the progression from acute myocardial infarction to myocardial fibrosis to heart failure. Front Cardiovasc Med (2021) 8:664044. doi: 10.3389/fcvm.2021.664044

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Zamaraev AV, Volik PI, Sukhikh GT, Kopeina GS, Zhivotovsky B. Long non-coding rnas: a view to kill ovarian cancer. Biochim Biophys Acta Rev Cancer (2021) 1876(1):188584. doi: 10.1016/j.bbcan.2021.188584

PubMed Abstract | CrossRef Full Text | Google Scholar

34. You Y, Zhao X, Wu Y, Mao J, Ge L, Guo J, et al. Integrated transcriptome profiling revealed that elevated long non-coding rna-Ac007278.2 expression repressed Ccr7 transcription in systemic lupus erythematosus. Front Immunol (2021) 12:615859. doi: 10.3389/fimmu.2021.615859

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Stein N, Berhani O, Schmiedel D, Duev-Cohen A, Seidel E, Kol I, et al. Ifng-As1 enhances interferon gamma production in human natural killer cells. iScience (2019) 11:466–73. doi: 10.1016/j.isci.2018.12.034

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Das PK, Asha SY, Abe I, Islam F, Lam AK. Roles of non-coding rnas on anaplastic thyroid carcinomas. Cancers (Basel) (2020) 12(11):3159–77. doi: 10.3390/cancers12113159

PubMed Abstract | CrossRef Full Text | Google Scholar