NGS-Based Molecular Karyotyping of Multiple Myeloma: Results from the GEM12 Clinical Trial

Simple Summary Multiple Myeloma (MM) is considered an incurable chronic disease, which prognosis depends on the presence of different genomic alterations. To accomplish a complete molecular diagnosis in a single essay, we have designed and validated a capture-based NGS approach to reliably identify pathogenic mutations (SNVs and indels), genomic alterations (CNVs and chromosomic translocations), and IGH rearrangements. We have observed a good correlation of the results obtained using our capture panel with data obtained by both FISH and WES techniques. In this study, the molecular classification performed using our approach was significantly associated with the stratification and outcome of MM patients. Additionally, this panel has been proven to detect specific IGH rearrangements that could be used as biomarkers in patient follow-ups through minimal residual disease (MRD) assays. In conclusion, we think that MM patients could benefit from the use of this capture-based NGS approach with a more accurate, single-essay molecular diagnosis. Abstract Next-generation sequencing (NGS) has greatly improved our ability to detect the genomic aberrations occurring in multiple myeloma (MM); however, its transfer to routine clinical labs and its validation in clinical trials remains to be established. We designed a capture-based NGS targeted panel to identify, in a single assay, known genetic alterations for the prognostic stratification of MM. The NGS panel was designed for the simultaneous study of single nucleotide and copy number variations, insertions and deletions, chromosomal translocations and V(D)J rearrangements. The panel was validated using a cohort of 149 MM patients enrolled in the GEM2012MENOS65 clinical trial. The results showed great global accuracy, with positive and negative predictive values close to 90% when compared with available data from fluorescence in situ hybridization and whole-exome sequencing. While the treatments used in the clinical trial showed high efficacy, patients defined as high-risk by the panel had shorter progression-free survival (p = 0.0015). As expected, the mutational status of TP53 was significant in predicting patient outcomes (p = 0.021). The NGS panel also efficiently detected clonal IGH rearrangements in 81% of patients. In conclusion, molecular karyotyping using a targeted NGS panel can identify relevant prognostic chromosomal abnormalities and translocations for the clinical management of MM patients.


Introduction
Multiple myeloma (MM) is a heterogeneous disease with a complex clonal and subclonal architecture and with few recurrent mutations [1,2]. The expansion of genetic events, including translocations, copy number variations (CNVs) and point mutations, facilitates tumor progression and is also responsible for the short remissions and numerous relapses along the disease timeline [3]. Some of these genetic abnormalities have a major impact on prognosis [4].
However, only a few translocations and CNVs define outcome risk in patients with MM. Cytogenetic aberrations are routinely evaluated using fluorescence in situ hybridization (FISH) panels and single-nucleotide polymorphism (SNP) microarrays in bone marrow (BM)-isolated plasma cells of MM patients [4,5]. Deletion of 17p (del(17p)) and translocations t(4;14), t(4;16) and t (14;20) were identified as poor prognostic biomarkers; likewise, gain(1q) and del(1p) are considered high-risk features for some groups [3]. The coexistence of different chromosomal aberrations can also modify the prognostic risk--for example, trisomy of chromosome 3 or 5 is associated with a better prognosis in patients with tumors displaying other high-risk genetic features [6,7], whereas del(1p32) in patients with t(4; 14), as well as del(6q) in patients with del(17p) and t (11;14) [8,9], increases the risk. In this context, the development of new molecular technologies has facilitated the classification of very high-risk patients who are characterized by either biallelic inactivation of TP53 or amplification (≥4 copies) of 1q21 in patients presenting ISS-III (International Stage System) [10]. The presence of this "double-hit" high-risk signature is currently considered a key prognostic factor [10][11][12].
The advent of next-generation sequencing (NGS) has ushered in a new era of diagnostic understanding and classification of the genomic landscape of MM [3] and has revealed the presence of unique and shared mutations in coexisting subclones of malignant/pathogenic plasma cells (PPCs) and the absence of a universal driver mutation [2,3,13,14]. While several mutated genes have been described to date in MM, only loss-of-function mutations of TP53 have had a clear prognostic impact. Accordingly, obtaining a complete genomic profile of MM, including single nucleotide variations (SNVs), insertion/deletion mutations (Indels), CNVs and translocations, would be crucial for the accurate diagnostic and prognostic stratification of MM patients. Along this line, there is great interest in the development of an NGS panel to identify--in a single assay--all possible chromosomal alterations (socalled molecular karyotype), which would considerably simplify the genomic diagnosis of patients with MM in a cost-effective manner [15][16][17] and would likely impact future targeted therapy approaches [18]. Although some efforts were made in this direction [15][16][17], the results have not been validated in cohorts of homogeneously treated patients with MM from clinical trials.
In the present study, we designed a capture-based targeted NGS panel to identify relevant genetic aberrations for the prognostic stratification of MM in a single assay. The panel aims to simultaneously detect SNVs, indels, CNVs, chromosomal translocations and V(D)J clonotypes [19] with a higher resolution and a smaller amount of required starting material than conventional techniques. NGS molecular karyotyping was validated and compared with data from FISH and whole-exome sequencing (WES). We also clinically validated the MM-specific NGS panel in samples obtained at diagnosis from patients with MM enrolled in the GEM12 clinical trial [20]. Our findings suggest that this new approach could eventually replace the conventional and standard procedures used in laboratories for the molecular diagnosis of MM.

Patient Cohort and Samples
We studied genomic (g)DNA from CD138 + plasma cells of BM aspirates from 149 patients with newly diagnosed MM enrolled in the GEM2012MENOS65 clinical trial [21]. Patients with active newly diagnosed MM were treated according to the PETHEMA GEM2012MENOS65 and the companion PETHEMA/GEM2014MAIN clinical trials [21]. Each study site's independent ethics committee reviewed and approved the protocol, amendments, and informed consent forms. The study was designed and conducted per the ethical principles of the Declaration of Helsinki and the International Council for Harmonization Guidelines. The study was registered at www.clinicaltrials.gov (accessed on 27 September 2017) as #NCT01916252 and EudraCT as #2012-005683-10. To summarize, patients received VRd plus autologous stem cell transplantation, VRd consolidation and maintenance with Lenalidomide or Lenalidomide plus Ixazomib.
The clinicopathological features of the patients are summarized in Table 1. The average age at diagnosis was 60 years (range: 42-65), with a similar frequency for males (53%) and females (47%). The clinical characteristics of the patient cohort were representative of the patients included in the clinical trial [21]. Plasma cells from BM aspirates were enriched using anti-CD138 + immunomagnetic beads, obtaining a purity >85% in all cases. DNA was extracted using the AllPrep DNA/RNA Mini Kit (Qiagen, Hilden, Germany) and capture libraries were generated with the SureSelectQ XT Reagent Kit (Agilent Technologies, Santa Clara, CA, USA) using 50 ng of genomic DNA. Libraries were multiplexed and sequenced on an Illumina NextSeq 500 platform (High Output Kit v2, Illumina, San Diego, CA, USA) in 150 bp paired-end mode.

Capture Next-Generation Sequencing Panel Design
Two panel versions were used in this study; however, differences between panel versions are mainly focused on the improvement of coverage in the 14q32 region for translocation identification. Design was performed to detect the most common and relevant aberrations known to occur in MM [3], covering the coding regions of 26 genes involved in the development and progression of tumors (ATM, ATR, ATRIP, BCL7A, BRAF, CCND1, CYLD, DIS3, EGR1, FAM46C, FGFR3, HIST1H1E, IRF4, KRAS, LTB, MAX, NRAS, NRM, PRDM1, PRKD2, RB1, TRAF3, ZFHX4, CRBN and NFKB2). Additionally, the complete exonic region was covered for TP53. Genes considered to be treatment targets or candidates for drug resistance in MM (PSMD1, XBP1, PSMB5, PSMC2, PSMC6, DDB1, IKZF1 and IKZF3) [16][17][18] and genes related to new immunotherapy treatments (CD38, CD19 and SLAMF7) were also included [22]. In addition, specific regions within canonical IGH locus (14q32) were included to capture IGH rearrangements, as well as germline polymorphisms distributed throughout the genome, to have a good representation of each chromosome to reliably detect CNVs. The total target size of the final design was 400 kb (SureDesign, https://earray.chem.agilent.com/suredesign/index.htm (accessed on 12 February 2019). Design and order information are available upon request to corresponding authors.

Bioinformatic Analyses
Raw FASTQ files were evaluated using quality control checks from FastQC [23] and Trimmomatic [24] was used to remove low-quality bases, adapters and other technical sequencing artifacts. Each FASTQ file was aligned to the human reference genome (GRCh38/hg38) using Novoalign (http://www.novocraft.com/products/novoalign/ (accessed on 1 March 2020). Optical and PCR duplication elimination was carried out with Sambamba [25]. SNVs and indels were identified through the combination of VarScan2 [26] and bcftools [27]. Variants were annotated using several functional (RefSeq, Pfam), population-related (dbSNP, 1000 Genomes, ESP6500, ExAC, genomeAD), in silico functional impact prediction (dbNSFP, dbscSNV), and cancer-related (COSMIC and ICGC) databases. Only variants described in the Catalog Of Somatic Mutations In Cancer database with prediction of the effect on protein function (codon stop, splicing, frameshift, nonsynonymous and Indels) were selected and filtered considering the following parameters: GMAF population frequency ≤0.01 (obtained from any database), total coverage at position ≥50, frequency of the variant in the sample ≥0.01, absence of annotations indicative of an artifactual-type distribution/quality of reads, and distance to the design bed regions ≤200. Potential CNVs were analyzed through an adaptation of the CONTRA v2.0.8 algorithm [28], using all samples as baseline. Breakpoint analysis was carried out with LUMPY [29] and Delly [30] to detect other potential rearrangements, selecting those events supported by more than 10 reads, followed by inspection of other rearrangement values (homology of consensus sequence, balance in the number of reads supporting each side of the breakpoint, frequency of the breakpoint in the cohort, etc.). IGH rearrangements were evaluated using MiXCR [31] and Vidjil [31], ruling out those supported by <5 reads or presenting with a frequency <0.01.

Whole-Exome Sequencing
WES data were obtained from peripheral blood cells (control) and PPC (tumor) samples from 51 patients to validate the mutations, CNVs and translocations. Exome sequencing was performed using the SureSelect Human All Exon Kit (v5, 51 Mb, Agilent Technologies, Santa Clara, CA, USA) by Macrogen Inc. (Geumcheon-gu, Seoul, Korea), with an output of 0.6 Gb per sample. Bioinformatics analyses were performed as described above, except for total coverage at position that we used ≥10 instead of ≥50. Those variants observed in control exome data were ruled out for consideration as germline variants. Loss of heterozygosity analysis was performed based on a selection of germline SNPs across the genome (minor allele frequency >0.3 and <0.7). Only CNVs observed on chromosomes 1 and 17 were reported.

Statistical Analyses
Statistical analyses were performed using RStudio and specific libraries from Bioconductor (https://bioconductor.org/ (accessed on 25 May 2021). Progression-free survival (PFS) was evaluated using Kaplan-Meier curves and two-sided log-rank tests. Multivariate analysis to determine hazard ratios was performed using the coxph and ggforest functions in R. The results were considered statistically significant at p < 0.05. Cytogenetic high-risk profile was considered if the selected CD138 + cells presented a translocation involving either chr4 or chr16 with chr14 (t(4;14)(p16;q32) and t(14;16)(q32;q23) or the presence of del(17p).

Targeted Capture Next-Generation Sequencing Panel and Exome Statistics
Sequencing statistics were calculated for each panel version (see Supplementary Tables S1A and S2A). Panel 1 had an average of 7.3 million reads (range 3-14 million), which translated into a mean base depth of 1746 (range 567-3272) with an average fraction of bases of 88.5% (range 84-94%) at >200×. The average reads in Panel 2 were 33 million for PPCs, resulting in a mean base depth of 2427× (range 66-5967×) and 71% (range 8.58-88.85%) of bases at >200×. Although the sequencing results were notably impacted by the quality of the sample DNA, the mean base depth was >500× for 91% (136 out of 149) of PPCs and the fraction of bases at >200× was >75% in 88% (131 out of 149) of cases from both panel versions. In addition, we estimated the mean base depth per million reads to be 239 (range 98-300) and 75 (range 2-165) for panel versions 1 and 2, respectively. Regarding the IGH genomic region, we observed an average of 1723 (range 709-3451) and 8691 (range 142-20,634) mapped reads in panels 1 and 2, respectively (see Supplementary Table S1B). Exome statistics are shown in Supplementary Tables S1C and S2. No significant differences were found between tumor and control sequencing data regarding the number of reads, coverage or depth distribution. The mean base depth was 66× (range 63-116×) and 80× (range 60-115×) for tumors and controls, respectively, indicating that no bias would be expected in SNV or CNV analyses.

Molecular Karyotyping Characterization of Patients with MM Using NGS Panels
The global molecular characterization of the 130 patients with at least one molecular alteration is shown in Figure 1. Below, we describe the different types of alterations, including SNVs, indels, CNVs, translocations and IGH rearrangements. SNVs and Indels: In total, 272 somatic mutations were found across 124 patients (see Supplementary Table S3), with an average of 1.7 (0-9) mutations per patient and at least one somatic mutation in approximately 82% of patients. Most of the mutations (253, 93%) were exonic and nonsynonymous SNVs, and 19 (7%) were short indels. As expected, the genes more frequently mutated in terms of percentage of patients carrying at least one mutation were NRAS (26%), KRAS (23%), TP53 (12%), BRAF (11%), DIS3 (8%), and ATM (8%).
IGH VDJ rearrangements: Overall, 185 functional IGH rearrangements were identified in 121 patients (81%) (Supplementary Tables S5 and S6). The number of rearrangements identified with MiXCR was larger than with Vidjil (165 vs. 117) and likewise was the number of patients (119 vs. 107) in which at least one rearrangement was identified. Rearrangements identified by MiXCR showed a higher average frequency (56% vs. 20%), but this could have been due to the final reads used to assemble the rearrangements (188 vs. 533 on average). A literature-based comparative study of major IGH clusters is presented in Supplementary Table S7 and Figure S1 shows a similar clonal detection rate between both techniques. In addition, IGH rearrangements observed in some of the patients were validated through PCR-based sequencing.

Comparative Results for NGS Panel, FISH and WES
Panel versus WES: Paired tumor-control WES data were available for 51 patients. Within the targeted regions, 47 somatic mutations were identified in the 51 patients using the NGS panel, and 37 of them (79%) were identified by WES (Supplementary Table S8). No somatic mutation was found in the WES data that was not present in the NGS panel data. The mean variant allele frequency (VAF) and depth for the NGS panel data were 0.23 (0.01-0.84) and 567× (134-920×), respectively, whereas the VAF and depth for the WES data were 0.26 (0.01-0.88) and 52× (9-182×), respectively. In the case of discordant positions (VAF = 0 in exome data), the mean VAF and depth were 0.04 (0.01-0.14) and 622 (166-891), respectively, for the NGS panel data, whereas in exome data, the average depth was 31 (10-50). The observed frequencies for del(1p), gain(1q) and del(17p) in the exome data were 10%, 29% and 4%, respectively, which were very similar to the frequencies obtained by the targeted approach (8%, 26% and 5%, respectively).  Table 2 (A)). When compared to the CNV results from the WES data (Table 2 (B) and Table  S9), PPV, NPV, sensitivity, specificity and global accuracy were 0.97, 0.96, 0.85, 0.99 and 0.96, respectively. Alteration-specific detailed PPV, NPV, sensitivity, specificity and global accuracy are shown in Table 2 and Table S10.

Conclusions
We show here that our targeted NGS panel based on molecular karyotyping can detect, in a single assay and with high precision, the key genomic alterations for MM stratification, including CNVs, SNVs, IGH clonal rearrangements and translocations. Moreover, the identification of an MM high-risk genetic profile with the NGS panel showed a

Conclusions
We show here that our targeted NGS panel based on molecular karyotyping can detect, in a single assay and with high precision, the key genomic alterations for MM stratification, including CNVs, SNVs, IGH clonal rearrangements and translocations. Moreover, the identification of an MM high-risk genetic profile with the NGS panel showed a suitable cor-relation with current gold standard techniques: FISH for CNVs and translocations and WES for CNVs and SNVs. Other targeted sequencing panels were previously reported [15,16,34] but with limitations in all cases, such as the inability to detect either translocations, CNVs or IGH rearrangements.
The genes often affected in MM are deeply covered by our panel, allowing the detection of both SNVs and indels at a similar frequency to that described in the literature. As expected, the comparative analysis with WES data revealed a greater sensitivity in coding regions targeted by our panel. WES-limited depth across targeted genes could have hampered the detection of mutations, with a VAF of 1-5% [2,35]. According to our results, current standard WES parameters and mean base depth would make it an unreliable tool for the identification of mutations within relevant genes in MM.
The sensitivity shown for translocations reached an average of 60% when compared to FISH data. In this regard, the small number of positive cases (23%) might have negatively affected estimates of the sensitivity to detect translocations. Published results are similar [15,34], pointing out the difficulty in detecting translocations by NGS. Since the nature of the genomic regions prone to produce translocation breakpoint is complex (masked, repetitive homologous regions), the efficiency of the whole process depends on random facts such as where DNA breakpoints occur during the fragmentation step or the efficiency of the probes specific of these regions.
Moreover, we observed no bias in the identification of different events, although del(1p) and t(14;16) presented a lower sensitivity. Taking into account the results obtained from CNVs using WES data, we concluded that this low sensitivity is probably due to FISH limitations in the case of del(1p), derived from the location of the probe; and panel limitations in the case of t (14;16). We are currently working on different improvements in the design of the NGS panel that could increase the sensitivity to detect some of the translocations and CNVs.
On the other hand, we found excellent global accuracy, PPV and NPV (approximately 90%) for the NGS panel when the results were compared with available FISH data. Given that FISH analysis has a limit of detection of approximately 10%, and that FISH is only able to show the presence of a maximum of two concomitant alterations in the same cell clone, we consider that our NGS panel can be a valid diagnostic tool.
The NGS panel was superior to FISH for CNV detection, which was confirmed by WES data. This is likely because FISH probes are designed in very specific locations within chromosomes and partial deletions/gains can occur in more telomeric/centromeric regions that do not affect the FISH signal. An example of this is shown in Supplementary Figure S2 (left), in which partial deletions/gains were identified by WES and the NGS panel but not by FISH. Additionally, as shown in Supplementary Figure S2 (right), alterations spreading a shorter length could be identified more precisely with either FISH or NGS panel data that go beyond the resolution of WES. The importance of these specific cases is shown by the difference observed in the impact on PFS. In our cohort of patients, the impact of gain (1q) was significant (p = 0.027) using panel data (see Figure 2A) and not when using FISH.
Disease progression is a key event in MM and the stratification of high-risk profiles remains an important goal to guide clinical decision-making [5]. Our results show that the accuracy of the NGS panel for the identification of high-risk patients with short PFS is superior to that of the other genetic techniques. Notably, our approach was able to identify patients with poor outcomes even though this series of patients received highly effective treatments [21].
The main limitation was the small number of true positive cases for some specific alterations, leading to a smaller sensitivity than recent studies [36]. Increasing the density of probes in certain regions for CNV and breakpoint detection is the next step. In return, a strength of our study is the fact that all the patients were homogeneously treated in a clinical trial.
In conclusion, molecular karyotyping by our targeted captured-based NGS panel specific for MM identifies SNVs, CNVs, IGH clonal rearrangements and translocations relevant to MM management. Although this approach has the potential to replace conventional cytogenetic approaches for the definition of genomic high-risk genomic profiles, simplifying the genetic study of MM, improvements in panel design will be addressed to improve the sensitivity to detect translocation events.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/cancers14205169/s1, Figure S1: Distribution of major IGH gene clusters involved in the rearrangements identified in 5 different studies [37][38][39][40]; Figure S2: Visual representation of CNV data in chromosomes 1 and 17 from two different patients; Table S1: Sequencing statistics; Table S2: Raw sequencing statistics for each sample. (A) Panel data. (B) WES data; Table  S3: SNVs and indels identified in MM patients; Table S4: Comparative table between FISH and NGS results, related to FISH probe locations; Table S5: IGH arrangement summary of PPCs; Table S6: IG rearrangements identified by using both MiXCR and Vidjil; Table S7: IG rearrangement identification comparative study; Table S8: Comparative mutational results from exome and panel data; Table  S9: CNV events identified by panel and exome analysis across chromosomes 1 and 17; Table S10: Estimation of PPV, NPV, sensitivity, specificity and global accuracy.

Funding:
The study was supported by the Accelerator Grant: next steps: Early detection and intervention: Understanding the mechanisms of transformation and hidden resistance of incurable hematological malignancies, grants from Instituto de Salud Carlos III PI15/01484 and CRIS Foundation 2020/0063 and 2021/0088.

Institutional Review Board Statement:
The study was conducted in accordance with the Declaration of Helsinki, and approved by the Ethics Committee of Instituto de Investigación Hospital 12 de Octubre (i+12)-FIBH12O, project associated code 15/01484, 14 of July 2015.
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data not presented in this article can be shared up on request.