Circulating microRNAs Correlate with Multiple Myeloma and Skeletal Osteolytic Lesions

Simple Summary Molecular biomarkers for the diagnosis of multiple myeloma and for the early detection of the associated osteolytic lesions are needed. MicroRNAs are a class of small non-coding RNAs that regulate gene expression post-transcriptionally and have been explored as circulating (extracellular) biomarkers for distinct diseases. Results show that miR-16-5p, miR-20a-5p, and miR-21-5p levels are differently expressed in the plasma of multiple myeloma patients compared with the control group and suggest that their combined expression could be used as a potential circulating biomarker. Furthermore, the expression of plasma microRNAs significantly correlates with myeloma bone disease and with bone lesions in the spine. Abstract Multiple myeloma (MM) is the second most frequent hematological disease and can cause skeletal osteolytic lesions. This study aims to evaluate the expression of circulating microRNAs (miRNAs) in MM patients and to correlate those levels with clinicopathological features, including bone lesions. A panel of miRNAs associated with MM onset and progression, or with bone remodeling, was analyzed in the plasma of 82 subjects (47 MM patients; 35 healthy controls). Results show that miR-16-5p, miR-20a-5p, and miR-21-5p are differently expressed between MM patients and healthy controls. Receiver operating characteristic analyses indicate that their combined expression has potential as a molecular marker (Area Under the Curve, AUC of 0.8249). Furthermore, significant correlations were found between the analyzed miRNAs and disease stage, treatment, β2 microglobulin, serum albumin and creatinine levels, but not with calcium levels or genetic alterations. In this cohort, 65.96% of MM patients had bone lesions, the majority of which were in the vertebrae. Additionally, miR-29c-3p was decreased in patients with osteolytic lesions compared with patients without bone disease. Interestingly, circulating levels of miR-29b-3p correlated with cervical and thoracic vertebral lesions, while miR-195-5p correlated with thoracic lesions. Our findings suggest circulating miRNAs can be promising biomarkers for MM diagnosis and that their levels correlate with myeloma bone disease and osteolytic lesions.


Introduction
Multiple myeloma (MM) is a malignant hematological disease characterized by excessive proliferation of abnormal plasma cells (PCs) in the bone marrow [1]. It accounts for more than 10% of all hematologic cancers and it remains incurable in the majority of cases, with a 45% 5-year survival rate [2,3]. Up to 80% of the newly diagnosed MM patients have osteolytic lesions that develop as a consequence of the exacerbated proliferation of malignant PCs in the bone marrow, causing an unbalance between the activity of boneforming cells (osteoblasts) and bone-resorbing cells (osteoclasts) [1,4]. This leads to the disruption of the normal bone remodeling process and to the development of myeloma bone disease (MBD) [5,6]. Consequently, MM patients are at high risk of skeletal-related incidents such as pathological fractures, spinal cord compression, vertebral collapse, severe bone pain and hypercalcemia [7]. Notably, the spine is the most frequent location of the bone lesions caused by MM and vertebral crush fractures are common [8]. Bisphosphonates are pyrophosphate analogues that cause apoptosis of osteoclasts and are currently the gold standard pharmacological approach to prevent and treat MBD [7,9]. More recently, denosumab, a humanized monoclonal antibody against RANKL, has also been used to decrease osteoclastogenesis and resorption by osteoclasts [7]. Imaging tools, such as whole-body low-dose computed tomography (CT), positron emission tomography-CT, and whole-body magnetic resonance imaging, are used to detect osteolytic lesions [10], and there are no currently available molecular biomarkers to detect MBD. Therefore, new tools for the early diagnosis of MM, as well as for the rapid identification of MBD, are needed. This will contribute to initiate the treatment at an early stage, avoiding further complications.
MicroRNAs (miRNAs) are a class of small RNA molecules that is transcribed from the genome but do not encode proteins [11]. In their mature form, miRNAs are about 19-22 nucleotides in length and are able to regulate the expression of messenger RNAs [11] and long non-coding RNAs [12]. Importantly, miRNA expression is disrupted in a vast number of diseases and hematological malignancies, such as MM [13,14]. Specifically, miRNAs are involved in MM onset and progression, which triggered their study as therapeutic targets [14]. Since the last two decades, circulating (extracellular) miRNAs have emerged as potential biomarkers for a wide range of diseases [11,15], including for MM [13]. Considering miRNAs can be detected in human blood serum/plasma, urine, saliva and other body fluids [11], these molecules are promising diagnosis tools through non-invasive or minimally invasive methods [15], leading to the development of several miRNA-based products for disease diagnosis. The stability of endogenous miRNAs in plasma samples is an important prerequisite for biomarkers [16]. It has been shown that miRNAs can be packed into extracellular vesicles [17] or can be associated with protein complexes [18], namely AGO proteins, protecting them from degradation by nucleases. Still, one of the main current challenges is to increase the sensibility and specificity of these circulating biomarkers [15].
The aim of this study is to identify circulating (extracellular) miRNAs to be used as biomarkers for the diagnosis of MM or MBD. For this purpose, a panel of miRNAs associated with MM onset and progression or with MBD pathophysiology, including regulation of osteoblasts and osteoclasts-related mechanisms, was evaluated. Levels of circulating miRNAs were analyzed in plasma samples and correlated with MM clinicopathological characteristics and with bone lesions, including in the spine. Results show that circulating miRNAs are promising candidates for the diagnosis of MM and for the MM-associated skeletal osteolytic lesions.

Patients and Plasma Samples
Human peripheral blood from patients diagnosed with MM (n = 47) was collected at Oncology Department at Hospital de Braga (Braga, Portugal) between July 2018 and January 2021. Patients with other tumor types or inflammatory diseases were excluded from the study. Peripheral blood from healthy donors (n = 35) was collected at Immunohemotherapy Department of Centro Hospitalar de São João (CHSJ, Porto, Portugal). The blood donors used in this study do not have history of cancer or inflammatory diseases, had not been diagnosed with infectious diseases and were not under the effect of anti-inflammatory or anti-depressive drugs. A written informed consent was obtained from all subjects before blood collection. The patients' clinical data was pseudonymized. The study was approved by Hospital de Braga and CHSJ Ethics Committee for Health, and the protocol conforms to the declaration of Helsinki. Clinical data was collected from medical records. Stage of the disease was determined following the International Staging System (ISS) for multiple myeloma, genetic alterations were assessed by cytogenetics/fluorescence in situ hybridization (FISH), and bone osteolytic lesions were assessed by magnetic resonance imaging, computerized tomography or conventional X-ray.
The blood samples were transferred from the EDTA collection tubes to 15 mL tubes, centrifuged at 1200× g, for 20 min, at room temperature (RT), and the plasma was collected. All samples were processed within three hours following blood collection to prevent hemolysis. No signs of colors from faint pink to bright red were detected in plasma samples following blood centrifugation. Plasma was further centrifuged at 2500× g, for 10 min, at 4 • C, to remove cell debris and stored at −80 • C in microtubes until further use.
RNA was converted to cDNA through TaqMan ® Advanced miRNA cDNA Synthesis Kit (Applied Biosystems), according to the manufacturer s protocol. Firstly, for the poly(A) tailing reaction, RNA samples were mixed with Poly(A) Buffer, ATP, Poly(A) Enzyme and RNase-free water, briefly vortexed and incubated in a thermocycler (MyCycler Thermo Cycler, Bio-Rad, Hercules, CA, USA) at 37 • C for 45 min followed by 10 min at 65 • C. For the adaptor ligation reaction, DNA Ligase Buffer, 50% PEG 8000, Ligation Adaptor, RNA Ligase and RNase-free water were mixed with the previous reaction and incubated at 16 • C for 1 h. Then, the reverse transcription was performed by adding RT Buffer, dNTP, Universal RT Primers, RT Enzyme Mix and RNase-free water and incubated at 42 • C for 15 min and 5 min at 85 • C. Finally, for cDNA amplification, miR-Amp Master Mix, miR-Amp Primer Mix and RNase-free water was added to 5 µL of the RT reaction product and incubated with the following settings: enzyme activation for 5 min, at 95 • C; denaturation at 95 • C for 3 s and annealing/extension at 60 • C for 30 s, for 14 cycles; stop reaction at 99 • C for 10 min. Then, cDNA was diluted in 50 µL of RNase-free water and stored at −20 • C until further use. For the RT-qPCR reaction, cDNA was mixed with TaqMan Fast Advanced Master Mix (Applied Biosystems), TaqMan Advanced miRNA Assay for miR-16-5p, miR-20a-5p, miR-21-5p, miR-29a-3p, miR-29b-3p, miR-29c-3p, miR-93-5p, miR-99a-5p, miR-146a-5p or miR-195-5p, and nuclease-free water in a 7500 Fast Real-Time PCR System (Applied Biosystems) and incubated as follows: 95 • C for 20 s; 40 cycles of 95 • C for 3 s and 60 • C for 30 s. RT-qPCR reactions were performed in duplicates. The results were analyzed in the 7500 Software v2.0.6 (Applied Biosystems). Exogenous cel-miR-39-3p and cel-miR-54-3p were used as normalizers. The relative gene expression was calculated using the ∆Cq method, in accordance with MIQE guidelines [26].

Statistical Analysis
Statistical data analysis was performed using the GraphPad Prism v8.1.2 software (GraphPad Software Inc, San Diego, CA, USA) and SPSS Statistics v26 (IBM Corp, Armonk, NY, USA). D'Agostino and Pearson omnibus normality tests were performed to evaluate if data followed a Gaussian distribution. Considering that data did not follow a normal distribution, the non-parametric Mann-Whitney U test was used to evaluate significant differences between groups. Categorical data was tested using the Fisher's exact test.
Relative expression of miRNAs was compared across groups using a linear regression, adjusted for the baseline characteristic identified as a possible confounder (age). The diagnostic value of the circulating miRNAs in plasma samples was calculated using Receiver Operating Characteristics (ROC) curve. To obtain the value of combined miRNAs expression a binary logistic regression was performed, considering the disease as the dependent variable and the miRNAs levels as covariates. For the ROC curve, the p-value tests the null hypothesis that the area under the curve equals 0.5. The cut-off points with the highest sensitivity and specificity were determined. Considering the data follows a non-parametric distribution, Spearman's rank correlation was used to calculate correlation coefficient and significance between miRNA expression and clinicopathological parameters. Recursive partitioning (RP) was used to split the dataset iteratively based on certain statistical criteria. The RP classification trees were performed using rpart package from R software and Gini index criteria. Statistical significance was considered for p < 0.05.

Patients Characteristics
Demographic, clinical and biochemical characteristics, treatment regimen, and genetic alterations in the patients' cohort are shown in Table 1. Patients with MM were significantly older than healthy donors (p < 0.001), but no differences were found between groups regarding gender. The majority of the MM patients (53%) are at stage II, according to the International Staging System (ISS). In addition, IgG was the immunoglobulin most commonly produced by these patients (70%). Genetic alterations are present in MM patients, mainly P53 and RB1 deletion in 23% and 35% of the patients, respectively (Table 1).

Differences in the Expression of Circulating miRNAs between Multiple Myeloma Patients and Control Group
The expression of the selected 10 miRNAs was evaluated in the plasma of 47 MM patients and 35 healthy donors. Results show that miR-16-5p and miR-20a-5p were significantly downregulated, while levels of miR-21-5p and miR-29c-3p were significantly upregulated, in MM patients compared with the control group ( Figure 1). Other members of the miR-29 family, namely miR-29a-3p and miR-29b-3p, were not statistically different between groups (Figure 1). Levels of four additional miRNAs, namely miR-93-3p, miR-99a-5p, miR-146a-5p and miR-195-5p, were also tested but no statistically significant differences were found between groups ( Figure 1). Notably, miR-16-5p and miR-93-3p were among the most highly expressed miRNAs in plasma samples with median CT values of 13.96 and 19.84 cycles, respectively. However, miR-93-3p shows high variability, both in the MM and in the control group. On the other hand, miR-146a-5p was the miRNA with the lowest expression in the plasma, with a median Ct value of 28.77 cycles. Importantly, all tested miRNAs show Ct values below the threshold of 35 cycles. Considering age as a variable that is significantly different between MM samples and the control group (Table 1), miRNAs miR-16-5p, miR-20a-5p, miR-21-5p and miR-29c-3p levels were tested again, adjusted for the age variable. The results from the multivariate linear regression analysis show that miR-16-5p, miR-20a-5p, and miR-21-5p maintain statistical significance, while miR-29c-3p differences could not be proven as age-independent (Table 2). Moreover, results from RP analysis show that the worst prognostic group refers to patients who have low expression of miR-20a-5p (miR-20a-5p < 0.011). The RP was subsequently used to maximize the variance across different nodes and minimize the variance within each node, allowing the identification of potential interactions between miRNAs biomarkers and MM disease. Using combined data, RP analysis identified different risk groups for MM disease (Figure 2). Only age-independent miRNAs differentially expressed between healthy individuals and MM patients were used for this analysis. The worst prognostic group referred to individuals who had low values of miR-20a-5p (94% out of miR-20a-5p < 0.011 individuals have MM disease). On the other hand, individuals with high levels of miR-20a-5p (miR-20a-5p ≥ 0.011), low levels of miR-21-5p (miR-21-5p < 0.0023) and high levels of miR-16-5p (miR-16-5p ≥ 4) or with high levels of miR-20a-5p (miR-20a-5p ≥ 0.011), low levels of miR-21-5p (miR-21-5p < 0.0023), low levels of miR-16-5p (miR-16-5p < 4) and low levels of miR-21-5p (miR-21-5p < 970 × 10 −6 ) had the best prognosis (90% of these individuals are healthy) (Figure 2). Herein, among the tested miRNAs, miR-16-5p, miR-20a-5p, and miR-21-5p were identified as differently expressed in the plasma of MM patients compared with healthy donors. < 0.0023) and high levels of miR-16-5p (miR-16-5p ≥ 4) or with high levels of miR-20a-5p (miR-20a-5p ≥ 0.011), low levels of miR-21-5p (miR-21-5p < 0.0023), low levels of miR-16-5p (miR-16-5p < 4) and low levels of miR-21-5p (miR-21-5p < 970 × 10 −6 ) had the best prognosis (90% of these individuals are healthy) (Figure 2). Herein, among the tested miRNAs, miR-16-5p, miR-20a-5p, and miR-21-5p were identified as differently expressed in the plasma of MM patients compared with healthy donors.    2. Decision tree created with 4 splits of the data. The root node (Node 1) contains all the data in the dataset. The algorithm splits the data based on a defined statistic, creating two new nodes (Node 2 and Node 3). Using the same statistic, it splits the data again at Node 2, creating two more leaf nodes (Nodes 4 and 5). For instance, Node 4 splits into Node 8 and 9; and Node 5 splits in two leaf nodes (Node 10 and 11). Lastly, it splits the data at Node 9, creating two more nodes (Nodes 18 and 19). The decision tree makes its prediction for each data row by traversing to the leaf nodes (one of the terminal nodes: Node 8, 18, 19, 10, 11 or 3).

Specificity and Sensibility of Circulting microRNAs as Multiple Myeloma Biomarkers
Next, attempting to find a miRNA or a combination of miRNAs that could successfully be used as a biomarker for MM, ROC analysis was performed for the miRNAs that were previously found differentially expressed between the MM group and the control group.
ROC curve analysis of miR-16-5p expression showed 78.72% sensitivity and 51.43% specificity in discriminating between MM and healthy donors, corresponding to an Area Under the Curve (AUC) of 0.665 (p < 0.05, 95% CI: 0.5458-0.7843). The miR-21-5p has an AUC of 0.654 (p < 0.05; sensitivity = 57.45%; specificity = 71.43%, 95% CI: 0.5365-0.7717). Considering the analysis of individual miRNAs, the highest AUC value was found for miR-20a-5p expression (0.707, p < 0.01; sensitivity = 36.17%; specificity = 97.14; 95% CI: 0.5961-0.8179). Nevertheless, combination of these three miRNAs improved the efficacy Figure 2. Decision tree created with 4 splits of the data. The root node (Node 1) contains all the data in the dataset. The algorithm splits the data based on a defined statistic, creating two new nodes (Node 2 and Node 3). Using the same statistic, it splits the data again at Node 2, creating two more leaf nodes (Nodes 4 and 5). For instance, Node 4 splits into Node 8 and 9; and Node 5 splits in two leaf nodes (Node 10 and 11). Lastly, it splits the data at Node 9, creating two more nodes (Nodes 18 and 19). The decision tree makes its prediction for each data row by traversing to the leaf nodes (one of the terminal nodes: Node 8, 18, 19, 10, 11 or 3).

Specificity and Sensibility of Circulting microRNAs as Multiple Myeloma Biomarkers
Next, attempting to find a miRNA or a combination of miRNAs that could successfully be used as a biomarker for MM, ROC analysis was performed for the miRNAs that were previously found differentially expressed between the MM group and the control group.

Correlation between Circulating microRNAs and Clinicopathological Variables in Multiple Myeloma Patients
Correlations of the miRNA levels in MM plasma samples with clinical and biochemical properties, treatment or chromosomal abnormalities was also tested. The miRNAs miR-20a-5p, miR-99a-5p and miR-195-5p negatively correlate with stage of the disease, following ISS classification (Table 3). Moreover, negative correlations were found between circulating miRNA levels and serum proteins. Among those, miR-16-5p, miR-20a-5p, miR-29b-3p, miR-29c-3p, miR-99a-5p and miR-195-5p were negatively correlated with β2M levels, a serum marker of tumor burden in hematologic malignancies [27,28], while miR-29a-3p was the only extracellular miRNA that negatively correlated with serum albumin (Table 3), a protein with prognostic value in MM patients [29,30]. Interestingly, besides miR-20a-5p, the three members of the miR-29 family showed a negative correlation with creatinine serum levels (Table 3), commonly used to measure the degree of insufficiency or renal failure [31]. No associations were found between the expression of circulating miRNAs and calcium levels ( Table 3).

Correlation between Circulating microRNAs and Clinicopathological Variables in Multiple Myeloma Patients
Correlations of the miRNA levels in MM plasma samples with clinical and biochemical properties, treatment or chromosomal abnormalities was also tested. The miRNAs miR-20a-5p, miR-99a-5p and miR-195-5p negatively correlate with stage of the disease, following ISS classification (Table 3). Moreover, negative correlations were found between circulating miRNA levels and serum proteins. Among those, miR-16-5p, miR-20a-5p, miR-29b-3p, miR-29c-3p, miR-99a-5p and miR-195-5p were negatively correlated with β2M levels, a serum marker of tumor burden in hematologic malignancies [27,28], while miR-29a-3p was the only extracellular miRNA that negatively correlated with serum albumin (Table 3), a protein with prognostic value in MM patients [29,30]. Interestingly, besides miR-20a-5p, the three members of the miR-29 family showed a negative correlation with creatinine serum levels (Table 3), commonly used to measure the degree of insufficiency or renal failure [31]. No associations were found between the expression of circulating miRNAs and calcium levels ( Table 3).
Regarding treatment of MM patients, a negative correlation was found between two members of the miR-29 family (namely miR-29b-3p and miR-29c-3p) and the use of chemotherapy or VTD regimen, while miR-93-5p and miR-146a-5p were negatively correlated with bone therapy (Table 4). Regarding genetic alterations, specifically del17p, del13q, translocation t (11;14), and translocation t (4;14), no significant correlations were found for any of the miRNAs (Table 4).

Correlation between Circulating microRNAs and Skeletal Osteolytic Lesions
Thirty one, out of the 47 MM patients included in the cohort, have bone lesions (65.96%). Of those, 55% have osteolytic lesions in the vertebrae, with or without lesions in other bone regions of the skeleton ( Figure 4A). The most affected regions are the thoracic (80.77%; T1-T12) and the lumbar vertebrae (69.23%; L1-L5), while sacrum and cervical lesions account for 30.77% and 19.23%, respectively ( Figure 4B). Considering lesions in spinal regions, 61.52% of the patients had lesions in at least two spinal regions, while only 7.69% had injuries in all the regions ( Figure 4C).
Our results show that circulating miR-29c-3p was the only miRNA significantly different between MM patients with and without bone lesions. Patients with bone lesions have a downregulation of about 36% in the expression of miR-29c-3p compared with patients without bone lesions, independent of the location of the lesions ( Figure 4D, Supplemental Table S2). In agreement, when analyzing all the miRNAs of the panel by RP, only miRNA-29c-3p was identified and the results show that the worst prognostic group refers to patients who have low values of miR-29c-3p (miR-29c-3p < 0.0074) ( Figure 4E). Importantly, age and gender are not significantly different between these two groups, miR-29c-3p being an independent marker. Furthermore, ROC curve analysis revealed that miR-29c-3p has discriminatory ability as diagnostic biomarker (AUC = 0.695; p < 0.05) ( Figure 4F). In agreement, miR-29c-3p is statistically correlated with bone lesions (r = −0.308; p = 0.042) ( Table 5). When considering lesions in specific spinal regions, miR-29b-3p correlates with cervical and thoracic vertebral lesions, while miR-195-5p correlates with thoracic lesions only (Table 5). No correlations were found between miRNAs expression and lesions in the lumbar or sacrum regions.

Discussion
Early detection of MM and of MM-related bone disease may lead to faster interventions and more effective treatments. Circulating miRNAs are excellent candidates to be used as biomarkers for the early diagnosis of diseases through non-invasive or minimally invasive methods [11]. In this study, plasma samples were screened for the expression of 10 miRNAs that had been selected due to their role in MM onset/progression or in the bone remodeling process. Herein, miR-16-5p, miR-20a-5p and miR-21-5p were identified as potential biomarkers for MM, particularly when combined together. Interestingly, Rocci et al. identified miR-16 and miR-21 as consistently expressed in at least 20% of MM samples, and further associated low miR-16 expression with poor overall survival [32], while Qu et al. showed miR-20a as downregulated in 33 MM patients versus 20 control samples in a Chinese cohort [33]. An advantage of using a panel of miRNAs compared with individual miRNA levels is that it increases accuracy and diagnosis value, as shown for other tumor types, including breast [34,35], lung [36], and bladder [37], among others. Importantly, differences in circulating miRNAs between MM patients and healthy controls do not necessarily translate the differential expression found in the malignant PCs versus normal PCs. Interestingly, the secretion of miRNAs to the plasma may be caused by the tumor cells themselves, by the cells in the tumor microenvironment, including immune cells, or by cells that have their function indirectly dysregulated by the tumor. The miR-16-5p is a tumor suppressor downregulated in MM-PCs, while miR-21-5p is an oncogene overexpressed in malignant PCs [25]. The same expression profile was identified in our cohort of plasma/blood samples from MM patients, compared with healthy controls, for these two miRNAs. However, this is not the case for miR-20a-5p. In vitro and in vivo studies showed that miR-20a-5p has an oncogenic role, being upregulated in MM-PCs [38,39]. Nevertheless, levels of this transcript were found to be decreased in the plasma samples from MM patients, suggesting alterations in miR-20a-5p secretion to the plasma may not be directly related with the changes in the malignant cells themselves, but rather related with changes in cells from the tumor microenvironment. This is in line with previous findings that miR-20a-5p levels are decreased in both the bone marrow microenvironment and in the peripheral blood of MM patients [40].
In this study, the association between circulating miRNAs and osteolytic lesions was also investigated, with a particular focus in the spine. The extracellular levels of miR-29c-3p were decreased in MM patients with bone disease compared with those without bone lesions. Recently, Papanota et al. identified let-7b-5p, miR-143-3p, miR-17-5p, miR-214-3p and miR-335-5p as associated with osteolytic bone disease in MM. Both studies are complementary as they use distinct miRNA panels, and contribute to reach the goal of finding novel non-invasive markers for osteolytic lesions [41]. Interestingly, Kapinas et al. showed that miR-29c-3p is highly expressed in fully mature and mineralized mousederived osteoblasts [42]. This miRNA is a member of the miR-29 family, known to directly target 16 extracellular matrix genes [43], including genes expressed in the bone. Although miR-29c-3p correlates with bone disease, no significant correlations were found when considering the distinct spinal regions. On the other hand, miR-29b-3p negatively correlated with vertebral lesions in the cervical and thoracic regions. This miRNA was shown to promote osteoblasts differentiation at early stages by directly targeting inhibitors of osteogenesis, while it reduces collagen levels at late stages to facilitate maturation of the collagen fibrillar matrix and, consequently, mineral deposition [44]. Circulating levels of miR-29b-3p were previously found to correlate with histomorphometric parameters in the bone [45]. However, to the best of our knowledge, this is the first study correlating extracellular miRNA levels and osteolytic lesions in the spine caused by MM.

Conclusions
Globally, we identified a promising combination of miRNAs in the plasma as biomarkers for MM. In addition, we specifically correlated circulating miRNA levels with the bone