Comprehensive genomic analysis of refractory multiple myeloma reveals a complex mutational landscape associated with drug resistance and novel therapeutic vulnerabilities

Abstract

The outcomes of patients with multiple myeloma (MM) refractory to immunomodulatory agents (IMiDs) and proteasome inhibitors (PIs) remain poor. In this study, we performed whole genome and transcriptome sequencing of 39 heavily pretreated relapsed/refractory MM (RRMM) patients to identify mechanisms of resistance and potential therapeutic targets. We observed a high mutational load and indications of increased genomic instability. Recurrently mutated genes in RRMM, which had not been previously reported or only observed at a lower frequency in newly diagnosed MM, included NRAS, BRAF, TP53, SLC4A7, MLLT4, EWSR1, HCFC2, and COPS3. We found multiple genomic regions with bi-allelic events affecting tumor suppressor genes and demonstrated a significant adverse impact of bi-allelic TP53 alterations on survival. With regard to potentially resistance conferring mutations, recurrently mutated gene networks included genes with relevance for PI and IMiD activity; the latter particularly affecting members of the Cereblon and the COP9 signalosome complex. We observed a major impact of signatures associated with exposure to melphalan or impaired DNA double-strand break homologous recombination repair in RRMM. The latter coincided with mutations in genes associated with PARP inhibitor sensitivity in 49% of RRMM patients; a finding with potential therapeutic implications. In conclusion, this comprehensive genomic characterization revealed a complex mutational and structural landscape in RRMM and highlights potential implications for therapeutic strategies.

Introduction

The introduction of novel therapies such as immunomodulatory agents (IMiDs) and proteasome inhibitors (PIs) has improved the outcomes of patients with multiple myeloma (MM), including those with relapsed MM (RMM) following second or third lines of therapy.¹ However, survival remains short if the disease becomes resistant to the major drug classes despite the advent of monoclonal antibodies.^2,3 Relapsed/refractory MM (RRMM) therefore represents a patient population of particularly high unmet medical need.⁴ Thus, a better understanding of the pathophysiology of RRMM is key to improving outcome of these patients.

In recent years, advances have been made in elucidating the genomic landscape of newly diagnosed MM (NDMM). These studies have revealed marked clonal heterogeneity with recurrently mutated genes each only affecting a minority of patients.^5-7 Clonal evolution over the course of the disease caused both by therapeutic interventions and ongoing genetic instability leads to more resistant clones in RMM, and ultimately to RRMM.^8-11 Accordingly, RRMM represents, in many respects, a fundamentally different biological disease entity. However, in contrast to NDMM, genomic data on RRMM is still limited.

Data from targeted sequencing in RRMM has identified an evolved set of mutated genes with enrichment for certain oncogenic drivers such as KRAS, NRAS, and TP53 mutations and the development of what are presumed to be resistance-associated mutations.¹¹ Mutations conferring resistance to PIs and IMiDs have been described in proteasomal subunits and in the IMiD target gene Cereblon (CRBN), respectively.^11,12 It is noteworthy that most of these mutations only occur in a minority of patients and at low allele frequency.¹³ Additional mechanisms of resistance may, therefore, be important in RRMM. An unbiased and comprehensive molecular study is, consequently, required to fully dissect the biology underlying RRMM.

Here, we used whole genome sequencing (WGS) and RNA sequencing (RNA-Seq) to comprehensively analyze a carefully selected cohort of 39 heavily pretreated RRMM patients with at least double refractory disease, revealing a complex mutational and structural landscape and highlighting potential implications for personalized therapeutic strategies.

Methods

Patient characteristics

WGS and RNA-Seq was performed on samples from 39 highly refractory MM patients. Nine of these samples have previously been analyzed and reported by our group using targeted sequencing.¹¹ The median number of prior therapy lines was five (2-13), all had relapsed following IMiDs and PIs, and 34 (87%) had had an autologous transplant. All patients were at least double-refractory, 62% were at least triple-refractory, and 38% were quadruple-refractory to IMiDs and PIs (Online Supplementary Table S1). In addition, 8% were also refractory to anti-CD38 monoclonal antibodies. Median progression-free survival (PFS) calculated from the time of sampling was 3.5 months and median overall survival (OS) was 7.4 months. Using fluorescence in situ hybridization (FISH)-based methods, 21 (54%) RRMM patients harbored high-risk cytogenetic aberrations; specifically, 16 (41%) had deletion 17p. For a comparison with NDMM, we further analyzed a WGS data set of 21 NDMM patients. For comparative analyses of mutational signatures, we included WGS data of 15 RMM patients with a median of two (1-4) prior therapy lines, thus less heavily pretreated than our RRMM cohort. Both additional WGS data sets are publicly available and had a median coverage of 40x.^14,15

Sample acquisition and preparation

Between March 2014 and October 2017, tumor samples from 39 RRMM patients were collected at Heidelberg University Hospital. Written informed consent was obtained prior to sampling in accordance with the Declaration of Helsinki. This study was approved by the Institutional Ethics Committee. CD138+ plasma cells were isolated, as described previously.¹¹ DNA and RNA were extracted using the AllPrep kit (QIAGEN, Hilden, Germany). Saliva, buccal swabs or bone marrow stroma cells obtained from cultured CD138- cells were used as germline controls. Saliva was collected in Oragene-Dx tubes and DNA was extracted using prepIT-L2P (DNA Genotek, Ottawa, Canada). DNA from buccal swabs was extracted using the blackPREP Swab DNA kit (Analytik Jena, Jena, Germany). DNA from stroma cells was extracted using the QIAamp DNA Mini kit (QIAGEN, Hilden, Germany).

Whole genome sequencing

As described previously,¹⁶ DNA libraries were prepared following the Illumina TruSeq Nano DNA Library protocol using the TruSeq DNA Nano kit (Illumina, Hayward, CA) and then sequenced on two lanes on the HiSeq X (2×151 bp) using the HiSeq X Ten Reagent Kit v2.5 (Illumina, Hayward, CA) to a median coverage of 77x.

Alignment and small variant calling

The raw reads were mapped to the human reference genome (build 37, version hs37d5), using BWA mem¹⁷ (version 0.7.8). To assess the effect of differing sequencing depths between the NDMM (median coverage 40x), RMM (40x), and RRMM (77x) samples on variant calling, subsampling of RRMM samples was performed for comparative analyses using Sambamba¹⁸ (version 0.6.6) to achieve a 50% lower coverage.

Small nucleotide variants (SNVs) were called using SAM-tools mpileup (version 0.1.19) and bcftools view. Indels were called using Platypus¹⁹ (version 0.8.1). Variants were annotated with Gencode²⁰ (version 19) and ANNOVAR.²¹ Splicing SNVs or SNVs resulting in nonsynonymous coding were called ‘functional SNVs’. For the prediction of the functional relevance of SNVs, we calculated Combined Annotation Dependent Depletion (CADD) scores (version 1.3) and used a cut-off of 20. Driver genes were identified using IntOGen (version 3.0.5).²²

Analysis of mutational signatures and Kataegis clusters

A supervised analysis of mutational signatures was performed with the R package YAPSA,²³ based on the mutational catalogue of the 30 known signatures from COSMIC v2 (https://cancer.sanger.ac.uk/signatures/signatures_v2), supplemented with the MM1 signature which was recently linked to melphalan exposure.^15,24 We defined Kataegis-like clusters to be regions of increased SNV density with at least five SNVs with, at most, a 1000 bp inter-mutational distance in one sample.²⁵

Calculation of HRDetect scores

HRDetect is a weighted model used to predict BRCA1/2 deficient tumors.²⁶ We used an implementation of HRDetect available at https://github.com/eyzhao/hrdetect-pipeline.²⁷

Detection of copy number variation and structural variants

Copy number states were called, as described previously,^16,25 and estimation of tumor purity and ploidy was performed using ACEseq (allele-specific copy number estimation from sequencing; https://www.biorxiv.org/content/early/2017/10/29/210807).

Structural variants (SVs) were detected using the DKFZ SOPHIA workflow (version 2.0.2, https://github.com/DKFZODCF/SophiaWorkflow).^16,25 SV candidate detection is a process of split-read and discordant mate evidence collection across each breakpoint as precursors for an SV. SV candidates (pairs of breakpoints) are filtered by a complex decision tree trained by expert assessment of orthogonal FISH data.

RNA sequencing

RNA-Seq libraries were prepared using the Illumina TruSeq stranded mRNA kit and were sequenced on the Illumina HiSeq 2000 V4 platform. The paired-end reads were mapped to the STAR index-generated reference genome (build 37, version hs37d5) with gencode (version 19) using STAR²⁸ (version 2.5.2b). The gene expressions were quantified using featureCounts (Subread version 1.5.1). Gene fusions were detected using Arriba version 1.0.0 (https://github.com/suhrig/arriba), as described previously.²⁹

A detailed description of the bioinformatics workflow and subsequent analyses is provided in the Online Supplementary Methods.

Results

High mutational load and genomic instability in RRMM

RRMM displayed a complex mutational landscape and an increase of both chromosomal and nucleotide aberrations compared to NDMM (Figure 1A, Online Supplementary Figure S1). With a median of 67 (range 7-496), the overall load of SVs was significantly higher in RRMM (P=0.002) (Figure 1B, Online Supplementary Figure S1). Complex structural rearrangements and catastrophic events were a frequent finding in RRMM patients, notably chromoplexy (n=10) and chromothripsis (n=5) (Figure 1C, D). Numerical chromosomal aberrations, which occurred significantly more frequently in RRMM compared to NDMM, included gain(1q) and deletions of 1p, 13q, and 17p (all P<0.05, Online Supplementary Figure S1).

The overall mutational load in RRMM was also significantly higher than in NDMM (P<10^-5) (Figure 1E, Online Supplementary Figure S1). We observed a median number of 116 (range 42-237) functional SNVs and five (range 1-15) functional indels per patient and, overall, a median of 3.94 somatic small variants per megabase in RRMM. RRMM patients showed a much higher prevalence of SNVs outside of Kataegis-like clusters than NDMM patients (P<10^-5), indicating activity of mutational mechanisms on a broader scale in RRMM (Online Supplementary Figure S1).

Genomic instability, assessed by calculation of the unbiased sum of homologous recombination deficiency (HRD), large-scale transition (LST) and telomeric allelic imbalance (TAI) scores, was significantly increased in RRMM compared to NDMM (P=0.004, Figure 1F, Online Supplementary Figure S1).

Driver gene aberrations in RRMM

Significantly mutated driver genes in RRMM featured prominently members of the mitogen-activated protein kinase (MAPK) pathway and TP53 (Figure 2, Online Supplementary Table S2) and showed a large overlap with known drivers in NDMM. However, the prevalence of BRAF (P<0.001) and TP53 (P<0.0001) mutations was significantly higher in RRMM compared to data from NDMM patients published by Walker et al.;⁷ a trend was also seen for NRAS mutations (P=0.05). In addition, we identified the following genes as significantly mutated drivers in RRMM: the sodium bicarbonate co-transporter SLC4A7 (13%); the Ras target MLLT4 (5%); the RNA binding protein EWSR1 (5%); the MLL complex member HCFC2 (5%); the COP9 signalosome subunit COPS3 (5%) (Figure 2). RNA-Seq data confirmed expression of the vast majority of variants and virtually all variants were predicted to be functionally relevant by CADD score (Online Supplementary Figure S2).

We found multiple bi-allelic events or ‘double hits’ in RRMM affecting known tumor suppressor genes (TSGs) (Online Supplementary Figure S3). In total, 25/39 RRMM patients presented with at least one double hit in TSGs, significantly increased compared to NDMM (3/21, P<0.001). The most frequently affected TSGs in RRMM were TP53 (n=8), RB1 (n=7), and TRAF3 (n=6). RRMM patients with a double hit of TP53 had significantly inferior PFS (P=0.004) and a trend for inferior OS (P=0.07); both calculated from time of sampling in univariate log-rank tests compared to patients with no or a single hit of TP53 (Figure 3). While treatment of patients before and after sampling was heterogeneous, the adverse impact of bi-allelic TP53 alterations remained significant in multivariate analysis including age, number of prior therapies and ISS as possible confounders with regard to both PFS (hazard ratio (HR) 4.02; P=0.01) and OS (HR 4.77; P=0.02). ‘Double hits’ in other single TSGs, or a combination thereof, did not show a significant impact on survival which would have been independent of TP53 events.

Figure 1.High genomic complexity and mutational load in RRMM vs. NDMM. (A) SV and SNV load per patient in RRMM vs. NDMM. SNV and SV counts are plotted for each patient in both cohorts showing a higher overall mutational load in RRMM (blue) vs. NDMM (red). Each dot represents an individual patient. The example cases shown in panels C and D are annotated as RRMM_16 and RRMM_15, respectively. (B) Differences in SV types in RRMM vs. NDMM. Median and range of number of overall SVs per patient in RRMM (blue) vs. NDMM (red) are shown as well as deletions (DEL), duplications (DUP), inversions (INV), and translocations (TRA). (C) Example case of RRMM displaying chromoplexy. Green lines represent translocations, blue lines deletions, red lines duplications, and black lines inversions. Transparency of lines is based on estimated SV clonality. Variant existence is represented by bar plots. The outer layer represents copy number variations and displays the copy-neutral nature of the chromoplexy event. (D) Example case of RRMM displaying chromothripsis. Green lines represent translocations, blue lines deletions, red lines duplications, and black lines inversions. Transparency of lines is based on estimated SV clonality. Variant existence is represented by bar plots. The outer layer represents copy number variations. (E) Genome-wide small variant mutational load in RRMM vs. NDMM. The number of mutations per patient and length of genome in megabases (MB) is shown in RRMM (blue) vs. NDMM (red). (F) Genomic instability scores in RRMM vs. NDMM. The unbiased sum of HRD, LST, and TAI scores is shown for RRMM (blue) vs. NDMM (red), illustrating a higher genomic instability in RRMM. To compensate for differing sequencing depths in both cohorts, the RRMM dataset was subsampled for these analyses.

Oncogenic networks in RRMM

Next, we analyzed whether gene groups, resistance mechanisms or signaling networks were recurrently affected by small variants in our RRMM cohort (Online Supplementary Table S4). We first addressed genes linked to the mechanisms of action of PIs and IMiDs (Online Supplementary Figure S3-5, Online Supplementary Table S5). We found that, overall, 21% of RRMM patients as opposed to only 5% of NDMM patients harbored mutations in genes considered relevant to PI activity. These were mainly protea-some subunits (PSMB5, PSMC2, PSMC6, PSMD2, PSMD11, and PSME3). Recurrent mutations were also detected in TJP1, previously found to modulate PI sensitivity in MM.³⁰ However, statistical significance was not reached (P=0.14), likely due to the limited sample size. Regarding IMiD resistance, candidate mutations were found to be significantly more frequent in RRMM patients compared to NDMM (31% vs. 5%, P=0.02). We observed mutations in all four members of the CRBN complex (CRBN, CUL4B, RBX1, DDB1), IZKF1 as downstream target, as well as in three members of the COP9 signalosome complex (COPS3, COPS4, COPS8) and in CAND1, which are regulators of cul-lin-RING ligase neddylation and which were recently identified in a CRISPR-Cas9 screen as being essential for the IMiD mechanism of action.³¹

Figure 2.Significantly mutated driver genes in RRMM. Significantly mutated drivers and their prevalence in the RRMM cohort are shown as well as copy number aberrations (CNAs) of chromosome arms 13q, 1q, 17p, 1p, presence or absence of hyperdiploid karyotype, and the genomic instability score.

Figure 3.Impact of TP53 alteration on PFS (A) and OS (B) in RRMM. PFS and OS are shown for RRMM patients with a bi-allelic event involving the TP53 locus (red), a mono-allelic event (blue) or no event (black), illustrating the inferior outcome of patients with bi-allelic TP53 aberrations. PFS and OS were calculated from time of sampling.

The MAPK pathway harbored mutations in 77% of RRMM patients (52% in NDMM), mostly due to mutations in NRAS, KRAS and BRAF. NFκB signaling was affected in 23% of such patients (5% in NDMM) with recurrent mutations in NFKB1, NFKB2, TRAF3 and CHUK. Both pathways showed a trend to being more frequently affected in RRMM, as did the large functional group of epigenetic modifiers which was affected in 77% of RRMM patients (52% in NDMM, all P=0.08). Of potential therapeutic interest, mutations in genes associated with sensitivity to PARP inhibitors were found in 49% of RRMM compared to 29% of NDMM patients, however, not reaching statistical significance (P=0.17), with recurrent mutations in AT M, NBN, and TOP3A (Online Supplementary Figure S3, S6, Online Supplementary Table S5).

Enrichment of mutational signatures of impaired DNA damage repair in RRMM

To identify mechanisms contributing to high mutational load and high genomic instability observed in RRMM as compared to NDMM, we analyzed mutational signatures (Figure 4A, Online Supplementary Figure S7). For further evaluation of different disease stages, we also included a cohort of relapsed but less heavily pretreated RMM patients as an intermediate stage. This cohort had previously been analyzed by Maura et al.¹⁵ In supervised fitting with signature-specific cutoffs, we found a significantly higher contribution from COSMIC signatures AC3 (associated with deficiency in homologous recombination repair) and the melphalan signature MM1 (P=0.006 and P<0.001, respectively) in RRMM compared to NDMM (Figure 4B) at the cost of the clock-like signatures AC1 (spontaneous deamination) and AC5 (clock-like but unknown) with significantly lower contributions (P<0.001 and P=0.005, respectively). In RMM, we detected signature MM1 but to a lesser extent when compared to our RRMM cohort (P=0.004), and signature AC3 was found with exposures in-between the NDMM and RRMM cohorts (Figure 4A, Online Supplementary Figure S7). A comparison of the clinical information on the published RMM cohort with our RRMM patients confirmed a significantly higher overall number of prior therapies as well as more extensive exposure to both novel agents and to high-dose melphalan therapy in our RRMM patients. Information on all mutational signatures detected in this manuscript, including asserted mutational mechanisms, is summarized in Online Supplementary Table S6.

To further assess the finding of enriched deficiency in homologous recombination repair as indicated by signature AC3, we applied HRDetect, which is a weighted model used to predict BRCA1/2 deficient tumors.²⁶ HRDetect scores were significantly higher in RRMM than in both NDMM and RMM (P<0.0001 and P=0.01, respectively, Figure 4C, Online Supplementary Figure S7). This was in line with the observation of higher exposure to mutational signature AC3 in RRMM than in NDMM and RMM. In 6/39 RRMM samples, the HRDetect scores exceeded the value of 0.7, a cutoff which was recently established to identify tumors with a high level of BRCA1/BRCA2 deficiency.²⁶ Taken together, these observations indicate a shift in the activities of different mutational mechanisms during the course of the disease. This might be a consequence of mechanisms intrinsic to the tumor cells (e.g., acquisition of DNA repair deficiencies) or may result from exposure to therapeutic agents.

Figure 4.Exposure to mutational signatures in RRMM vs. RMM vs. NDMM patients. (A) Absolute exposure to mutational signatures in RRMM vs. RMM vs. NDMM patients. Exposure to mutational signatures based on the Alexandrov COSMIC (AC) catalogue with the addition of the MM1 signature recently linked to melphalan exposure is shown for RRMM (blue) vs. RMM (green) vs. NDMM (red) patients. Most notable is an increased impact of signatures AC3 (light brown) and MM1 (black) in RRMM. (B) Relative exposure to mutational signatures in RRMM vs. RMM vs. NDMM patients. Exposure to mutational signatures is shown for RRMM (blue) vs. RMM (green) vs. NDMM (red) patient cohorts. Significant differences in exposure are indicated with the following P-values: * <0.05, ** <0.01, *** < 0.001. Most notable is a significant increase of signature AC3 in RRMM compared to NDMM and RMM as well as a significant step-wise increase of signature MM1 in RRMM vs. RMM vs. NDMM. (C) HRDetect scores in RRMM vs. RMM vs. NDMM patients. HRDetect scores in RRMM (blue) vs. RMM (green) vs. NDMM (red) indicate a significant increase in impaired homologous recombination repair features in RRMM. To compensate for differing sequencing depths in both cohorts, the RRMM dataset was subsampled for this analysis.

Enriched chromosomal translocations in RRMM

Next, we sought to identify translocations, which were enriched or unique to RRMM (Figure 5A, Online Supplementary Figure S8). While MYC rearrangements per se were not enriched in RRMM, FAM46C-MYC translocations (n=4) and local rearrangements of chr8q24.21 (n=4) were both found exclusively in RRMM patients (Figure 5B). The super-enhancer of FAM46C showed rearrangements involving a number of target genes such as the transcription factor and presumed oncogene LMO4. An analysis of expression data confirmed enhancer hijacking (Online Supplementary Figure S9). IgH translocations involving MYCN were observed in two cases as part of composite t(4;14)-t(2;4) translocations also involving the MMSET locus, co-activating MMSET and MYCN. In these cases, MYCN was highly expressed, while MYC expression was completely suppressed (Online Supplementary Figure S10). Both cases showed extramedullary disease, secondary plasma cell leukemia and a distinct plasmablastic cytological appearance. In one case, retrospective analysis of earlier samples revealed the absence of the secondary t(2;4) translocation and an absence of concomitant MYCN activation and MYC suppression, suggesting the emergence of MYCN overexpression with disease progression.

Discussion

This comprehensive study of extensively pretreated, highly refractory MM by WGS and RNA-Seq suggests that the pathophysiology of NDMM and the biology of refractory disease are strikingly different. In RRMM, we observed marked genomic instability with impaired DNA repair mechanisms, in particular homologous recombination repair (HRR), which we confirmed using the well-established HRDetect algorithm. Mutational signatures indicative of impaired DNA repair, such as AC3, have been reported in NDMM although not to the extent observed here in this highly refractory setting.^32,33 A recent study has reported the absence of mutational signature AC3 and HRR deficiency in NDMM and has classified the appearance of signature AC3 in NDMM as a false positive effect of the applied deconvolution method.¹⁵ A study on the genomic make-up of relapsed but less heavily pretreated RMM using whole-exome sequencing did not report mutational signature AC3, but rather a novel signature associated with alkylator therapies.¹³ In our work, we find both mutational signatures and thereby provide evidence that they are not reflective of the same underlying mutational process and can be deconvoluted with sufficient statistical power in the underlying data.

Figure 5.Immunoglobulin translocations and MYC rearrangements in RRMM. (A) Immunoglobulin translocations in RRMM. Translocations involving the immunoglobulin heavy chain (IGH) locus, the lambda light chain (IGL), and the kappa light chain (IGK) locus are shown. The number of patients with involvement of the respective partner genes are given in brackets. Patient RRMM_34, harboring a very complex IGL translocation, was excluded from this graph for reasons of readability and is shown separately in Online Supplementary Figure S8. Cytobands, chromosome arms and chromosomes were also stretched and compressed to emphasize targets of immunoglobulin translocations and to improve readability. (B) MYC rearrangements in RRMM. Most notable are local rearrangements as well as those involving the IGL locus or FAM46C. Orange lines represent immunoglobulin locus translocations; purple lines secondary immunoglobulin locus related translocations. Secondary translocations of the immunoglobulin loci were defined as secondary events with one of the breakpoints of the SV not further away than 2MB from the target breakpoint of a given primary immunoglobulin translocation (i.e., IG —PrimarySV—> PrimaryTarget —SecondarySV—>SecondaryTarget). Green lines indicate non-IG locus translocations, blue lines deletions, red lines duplications, and black lines inversions.

Walker et al. recently used the term double hit MM to indicate a subset of NDMM with a very poor prognosis, including patients with bi-allelic TP53 inactivation.³⁴ In our study, patients with a TP53 double hit experienced the worst outcomes in RRMM, confirming recent observations in the relapsed setting.^10,35 Further genes with recurrent biallelic aberrations in RRMM include RB1 and TRAF3, though they did not independently affect the prognosis in our cohort. This was likely attributable to the limited number of patients who had ‘double hits’ in these TSGs but not concomitant bi-allelic TP53 events.

BRAF mutations occur at a higher prevalence in RRMM. The frequency seen in our cohort seems to exceed even those previously reported in RRMM,³⁶ a finding which may reflect the heavily pretreated nature of the patient cohort presented here. However, a possible selection bias has to be borne in mind as a strong driver mutation such as BRAF V600E might be associated with higher tumor load, thus potentially resulting in a higher success rate of plasma cell purification for WGS. Further, the limited sample size might have contributed to this finding. Nonetheless, the eight cases of the druggable mutation BRAF V600E are of particular therapeutic interest.^37,38

Among potential resistance-conferring SNVs, individual genes were only affected at low frequencies. However, at the level of functional networks, recurrence was actually seen: we found mutations in several proteasomal subunits as well as in TJP1 which modulates PI sensitivity in MM.^12,30,39 While the functional impact of most of these mutations remains to be proven, their enrichment in RRMM supports an association with PI resistance. The same holds true for mutations in genes presumably associated with IMiD resistance, such as CRBN, CUL4B, and IZKF1. Furthermore, we detected mutations in three members of the COP9 signalosome complex (COPS3, COPS4, and COPS8) and CAND1, further supporting their functional impact on IMiD activity.³¹

The increased mutational load might explain both the higher capacity of MM cells to adapt to treatment and facilitate the emergence of resistance. One of the mechanisms contributing to the higher mutational load in RRMM appears to be impaired DNA double-strand break HRR with the potential therapeutic implication of synthetic lethality to pharmacological inhibitors of DNA damage response, such as ATR inhibitors. At the same time, we found mutations in genes associated with sensitivity to PARP inhibitors in 49% of RRMM patients. This further strengthens the rationale for assessing the therapeutic efficacy of ATR and PARP inhibition in RRMM as has been shown, for example, in solid cancers with BRCAness characteristics.⁴⁰ In fact, there is plenty of pre-clinical evidence for synthetic lethality conferred by such inhibitors in MM cell lines with high rates of ongoing DNA damage.^41-43 As PI treatment has been suggested to induce a BRCAness-like state in MM cells via impairment of DNA repair pathways,⁴⁴ there may also be a rationale for combining PI and PARP inhibitors in MM. However, these concepts need to be confirmed in clinical trials. One major limitation of our study is that comparative analyses were performed between independent NDMM, RMM, and RRMM patient cohorts, limiting the ability to draw conclusions regarding the tumor evolution under treatment. Thus, one can only speculate as to whether genomic instability was pre-existing in these cases, e.g., restricted to certain focal lesions,⁴⁵ and was then selected for or whether it was newly acquired following treatment with e.g., DNA damaging drugs. Longitudinal analyses that also address spatial heterogeneity will, therefore, be of particular interest. While this real-world RRMM cohort represents similarly ultra-refractory MM patients, the route taken to this end-stage disease differed greatly between individual patients. It is well conceivable that important biological differences exist between RRMM patients reaching refractoriness after multiple lines of treatment and those progressing quickly through a limited number of therapies. Our RRMM cohort size might, therefore, be too small to overcome the heterogeneity of this patient population with sufficient statistical power. Further analyses on larger or more homogeneous RRMM patient cohorts will help elucidate these issues.

In conclusion, based on our observation that RRMM is characterized by marked genomic instability, which enables MM cells to rapidly adapt to selective therapeutic pressure, treatment strategies focused on exploiting impaired HRR should be evaluated within prospective clinical trials. Such strategies might be particularly useful in the current era of novel immunotherapies in MM as recent reports suggest genomic instability as a mechanism of resistance to CAR T cell treatment.^46,47 Targeting impaired DNA repair mechanisms may, therefore, help to improve the outcomes of patients with RRMM.

Footnotes

• Received June 7, 2021
• Accepted January 7, 2022

Correspondence

*NG, NP, UHT, and DH contributed equally as co-first authors.

^#NW, MS, and MSR contributed equally as co-senior authors

Marc S. Raab m.raab@dkfz-heidelberg.de

Nicola Giesen nicola.giesen@med.uni-heidelberg.de

Disclosures

NG received honoraria from MSD and advisory board from Pfizer; EKM received honoraria from Janssen, Celgene, Takeda, has a consulting or advisory role with Janssen, Celgene, Takeda, received research funding from Takeda and travel support from Janssen, Takeda, Celgene and Mundipharma; HG received grants and/or a provision of Investigational Medicinal Products from Amgen, BMS, Celgene, Chugai, Dietmar-Hopp-Stiftung, Janssen, John Hopkins University, Sanofi and research support from Amgen, BMS, Celgene, Chugai, Janssen, Molecular Partners, MSD, Sanofi, Mundipharma, Takeda, Novartis, advisory boards from Adaptive Biotechnology, Amgen, BMS, Celgene, Janssen, Sanofi, Takeda and honoraria from ArtTempi, BMS, Celgene, Chugai, Janssen, Novartis and Sanofi; MSR received honoraria, consultancy, research grants and travel support from Novartis and Amgen. The other authors declare no competing interests.

Contributions

NG, NW, MS, and MSR designed the study, collected and analyzed data, and wrote the paper; NP, UHT, and DH analyzed the data and wrote the paper; JX, SU, SB, MF, SSM, AJ, and BB analyzed data; EKM, CMT, NM, and HG collected data. All authors critically reviewed and approved the final version of the manuscript.

Data-Sharing statement

Sequence data has been made available at the European Genome-phenome Archive (EGAS00001004363).