ABSTRACT

INTRODUCTION

Multiple myeloma (MM) is a hematological neoplasm that appears after transformation and uncontrolled proliferation of plasma cells. MM is characterized by a heterogeneous genetic aberrations and very different clinical outcomes (Avet-Loiseau, . Boyd KD). Although treatments for MM have increased by develop new and more sophisticated therapeutics drugs such an immunomodulatory drugs (IMiDs) and proteasome inhibitors (PIs)( Shaji K. Kumar, Raghupathy, Antonio Palumbo) the marked variability of responses indicate that larger studies will be required.

Currently, identification of cytogenetic abnormalities is performed by conventional karyotyping and fluorescence in situ hybridization. However, these disorders are not sufficient to explain the malignant phenotype given that are also observed in premalignant states of MM such a monoclonal gammopathy of unknown significance (MGUS) or smoldering myeloma(SM)(Bergsagel, Hideshima). This finding justifies the need for a comprehensive screening of genetic abnormalities in MM patients, which has not been incorporated in the medical workup yet. Recently, the introduction of massive sequencing of patient genome using next-generation sequencing (NGS) technologies has considerably increased the understanding of the biological features of MM. Many works have described the complex and heterogeneous mutational profile of MM patients(bolli nuevo, walker). Whole Exome Sequency (WES) studies in newly diagnosis MM patients have confirmed that more than 50 genes are mutated in the first manifestation of disease (walker). However, only few genes have been detected recurrently mutated at diagnosis, including KRAS, NRAS, BRAF, DIS3, TP53 and FAM46C, and only TP53 mutations have been recurentelly associated with poor survival. In addition, other studies have assessed the clonal evolution over time, pointing out that systemic treatment with chemotherapy may affect the livelihood of some subclones more than others, and thus may influence the tumor evolution over time(Egan JB, bolli, Keats JJ, ) The introduction of targeted studies allows the detection of mutations even with very low allele frequencies at an affordable price, allowing the incorporation of extensive genetic studies to the clinical workup. In the last years, several groups have applied this approach in order to achieve a better patients stratification and prognosis prediction. Although many studies have highlighted the importance of the subclonal landscape in MM and many efforts have been undertaken to stratify patients and predicts their responses, there is no clear relation between sensitives or refractories clones to treatment, and more information about the prognostic impact of this subclonal profile in series of homogeneously treated MM patients is needed. A large number of clinical trial are being carried out with this aim, unify treatments in order to study more effectively the impact of genetic alterations in prognosis.

In this work, newly diagnosis MM patient homogenously treated have been genetically characterized using a combination of the most recent techniques, including FISH and ultra-deep targeted sequencing. We applied the highest read depth described in the literature with the aim to detect minority subclones ignored to date. We also integrated these data with the clinical features to find out new patterns of behavior, relate them with survival and reveal new insight into the complexity of clonal and subclonal architecture of MM.

Patient samples

Samples were taken from the available 79 newly diagnosed MM patient enrolled in the clinical trial GEM10MAS65 (registered at www.clinicaltrials.gov as #NCT01237249). This is a phase III trial where patients older than 65 years were randomly assigned between two treatments arms: sequential melphalan/prednisone/Velcade (MPV) followed by Revlimid/low dose dexamethasone (Rd) versus alternating melphalan/prednisone/Velcade (MPV) with Revlimid/low dose dexamethasone (Rd). Progression free survival (PFS) and overall survival (OS) were measured from the starting date of the treatment. The median time to progression was 26.4 months with a median follow up of 31.5 months.

Targeted sequencing and mutation calling

Positive plasma cells CD138 were isolated from bone marrow aspirates and DNA was extracted using AllPrep DNA/RNA mini kit (Quiagen). Only 20 ng of DNA were used to prepare libraries using Ampliseq Library Kit 2.0. We also sequenced DNA from the 15 available CD138 negative fractions in order to filter out potential artifacts and corroborate detection sentivity. Samples were sequenced using Ion Torrent platform (IonProton, Thermofisher, Carlsbad, CA, USA) using the M3P gene panel (Mayo Clinic, Arizona). This panel spreads out over 77 genes frequently mutated in MM, which are related to critical pathophysiological pathways, associated to drug resistance or targetable with molecular drugs [paper mayo kortum etal]. Quality filter and alignments was performed using Torrent Suit software (Life Technologies) Single nucleotide variants were calling and annotated using Ion Reporter software applying in-house modifications in call variants process. Variants listed in Single Nucleotide Polimorphism database (dbSNP, http://www.ncbi.nlm.nih.gov/SNP/) were excluded from samples without germline available, as well as variants that were detected in multiple samples. In addition, to test the ability of the workflow previously described, we applied a novel bioinformatics pipeline developed by Spanish National Cancer Research Centre (CNIO). All reported mutations were detected by both bioinformatics approaches. The integrative genomic viewer (IGV) was used to visualize the read alignments, single variants and correct sequencing errors due to homopolymer regions.

Statistical analysis

All statistical analysis was performed using the statistical environment R. Correlation coefficients between mutated genes and cytogenetic aberrations was assessed and plotted using corrplot (https://cran.r-project.org/web/packages/corrplot/). Differences in survival were tested using the log-rank test. Cox proportional hazard regression was employed to obtain hazards ratios (HR) and evaluated at 5% of significance level. A second approach called LASSO (least absolute shrinkage and selection operator) was performed to detect relevant variables among clinical, cytogenetic and mutated genes.  Further details can be found in Data Supplement.

RESULTS

Mutated genes and altered pathways (cambiar título)

We sequenced 79 tumor samples with a mean coverage depth of 1600X. The minimum coverage of the detected variants was 60X and the average coverage 370X. We identified 170 nonsynonymous missense/nonsense/stoploss single variants, 81 of them (48%) were predicted pathogenic by Sift and Polyphen and 61 (36%)  have been described in COSMIC data base.  85% of patients harbored at least 1 mutation with a median of 2.1 mutations per patient. We detected mutations in 53 genes (Figure1), although 6 genes accounted the 39% of the total number of mutations; KRAS= 21.5%, DIS3= 19%, NRAS= 16.5%, BRAF= 10.1%, TP53= 8.8% and ATM= 7.6% of the patients.   48% of patients (38/79) presented at least one mutation in genes envolved in RAS/MAPK pathway, being the most frequently mutated pathway. 72 and 100% of variants in KRAS and NRAS respectively were detected in the hotspot codons 12, 13 and 61, and the targetable V600E BRAF mutation was detected in 1 patient. (Figure 2 supplemental?). NFKB pathway was the second most frequently mutated in our cohort, accounting for the 15% of all mutations distributed in 25 % of the patients (19/79). This pathway included TRAF3 (5 mutations, one nonsense and 4 missense) and TRAF2 (3 missense mutation) all of them predicted pathogenic by Sift and Polyphen. Other pathways importantly altered in the cohort were MYC in 11% of patients (9/79), cereblon and ciclyn both in 9% of patients (7/79).

Multiple mutations within the same gene were observed in 11 patients: 4 of whom harbored 2 and 3 mutations within DIS3 (patient 1-24= Glu501Lys and Phe120Leu at 8 and 53 % of VRF, patient 2-236= Asp487His and Asp479Glu at 4 and 21 of VRF respectively, patient 3 – 321= Tyr753Asn and Glu126Lys at 2 and 58 % of VRF respectively and patient 4-42 = Arg820Trp, Gly249Glu and at 14, 20 and 24% of VRF respectively). Other 2 patients showed 2 and 3 mutation in KRAS (patient 5-168= Gly13Asp and Gln61His at 9 and 13 %of VRF and patient 6-269 = Tyr71Asp, UTR3′ in exon 6 and Gln61Glu at 3, 15 and 37 % of VRF respectively; 2 patients with 2 mutation in NRAS(patient 7-177= .Gln61Lys  and Gly12Ala at 5 and 12% of VRF respectively and patient 8-257=  Gln61Glu  and Gly12Ala at 5 and 6% of VRF respectively), one patient with 3 mutations in MAX (patient 9-190= Arg36Lys , Arg35Leu and Glu32Val at 10, 20 and 26 %of VRF), one patient with 3 mutations in TRAF3 (patient 10-40 = Lys453Asn, His136Tyr  and Phe445Leu 8, 11 and1 3% of VRF) and one patient with 2 mutation in TP53 (patient 11-40 = Asp208Val  and Glu204Ter at 35 and 36% of VRF respectively).

Variant Read Frequency study

The VRF found in our cohort were diverse.  We detected 50% of variants (85/170) below 25% of VRF and 27% (46/170) below 10% (Figure 2). KRAS (n=5) DIS3 (n=5), BRAF (n=4), NRAS (n=4), and TP53 (n=3) were the genes that harbored the largest number of low frequency mutations (Figure 2).

KRAS, NRAS, BRAF and TP53 mutations were, in all cases, lower than 50% of VRF while DIS3 showed mutations in a broad range (from 2 to 85%)(Fig3). ). Most of DIS3 mutations with VRF <10% (4/5) are located in the activity domain RNB, all considered pathogenic by Sift and Polyphen.  TRAF3 (n=5) was the only gene where only mutations lower than 25% of VRF were detected. Furthermore, 4 of the 5 mutations were missense and predicted damaging by Sift and Polyphen (Shit value=0, Polyphen value=1) and the remaining mutation was nonsense. No gene showed low frequency mutations restricted to a particular region, but rather spread along different exons and codons (Figuras genes4).  2 mutations at very low VRF were detected in IRF4 (3 and 4%) and were predicted damaging by Sift and Poliphen (Sift value=0, Polyphen value=1). Surprisingly, these patients underwent an aggressive progression with a premature death.

The large number of mutations detected at very low VRF shows the need to perform depth targeted sequencing, since WES is not powerful enough to detect minority subclones and crucial genetic information could be ignored[JM1].

We should talk about the clinical significance and prognostic implication of low frequency variants.  Survival impact.

Survival analysis and correlation study of clinical variables and the mutational profile

As might be expected, demographics and clinical variables did not show any impact in survival. We also evaluated the impact on survival of each mutated gene individually and grouping according to their pathways. Tree genes, TP53, TRAF3 and CDKN1B, had a significant and negative impact in OS (P-value = 0.00451, 0.00751 and 0.0377 respectively) but not in PFS. All mutation in TP53 are located in the DNA-binding domain except one, which means that  ….  Regarding to altered pathways, neither RAS nor NFKB pathways had any impact in survival, as have been previously reported by other authors. However, mutations in genes of Cereblon pathway (IRF4= 4 and CUL4B = 1) were associated to short OS ( P-value = 0.0309). We also evaluated whether harboring more than one mutation in the same gen effect patient survival and we did not found any significant impact.  Regarding to cytogenetic impact on survival, we found a trend toward shorter SLP (P-value = 0.0922) in patients with at least one high-risk cytogenetic aberration (del17p, t4; 14 and t 14; 16). When mutation in TP53 were additionally considered for this group, the tendency became clearer (P-value = 0.0628). We also observed a negative impact on survival according with the number of mutations; patients with more than 3 mutations died prematurely, with a negative impact in OS (P-value = 0.00148). In these patients with more than 3 mutations, the 30% of mutations are presented in a very low VRF, lower than 10% (15/51). The identifications of these high-risk patients subgroup had not been possible by the conventional targeted sequencing depth. Correlation study between mutated genes and cytogenetic aberration did not show any significant relation (Figure 5).

DISCUSSION

In this study, we report a genetic characterization of 79 homogenously treated patients at diagnosis using a 77 genes customized MM panel (M3P) applying the highest deep read to date. The most frequently mutated genes and pathways found in this work were consistent with previously reported (poner walker y bolli Nuevo) (KRAS 21 vs 21%, NRAS 18 vs 16%, BRAF 7 vs 8% and TP53 8 vs 8%) except for DIS3, where an increase in its mutational incidence was found (9 vs 19%). This fact may be explained by the high number of sublonal mutations with <10% of VRF found for this gene, which have been not detected in other studies. We again confirmed in our cohort the crucial role of alteration of the RAS/MAPK and NF-kB pathways in myeloma cells, which account for the 66% of total mutations. Despite of out homogeneous cohort of patients, and also, he high depth read applied to samples in this study which allowed to detected minority subclones never observed before, mutations in these pathways had no prognostic significance.

Interestingly, patients with non-synonymous mutations in genes related to Cereblon pathways (IRF4 and CUL4B) showed a shorter significant survival in terms of SLP as well as OS. Since recently studies have identified Cereblon pathway as essential for anti-MM activity of immunomodulatory drugs (IMIDs)(11-18 del paper refractory), this finding is especially relevant due to patients enrolled in this clinical trial were treated with Lenalidomide, suggesting a direct and negative influence of cereblon alteration in patient response. There is a low incidence of mutations reported for this pathway in the literature and our knowledge of IMIDs resistance mechanisms is still restricted. Therefore, extensive studies in larger cohorts are required in order to confirm the negative impact in survival of crbn mutations in patients treated with IMIDs.  Previously publications pointed out a clear relationship between high cytogenetics risk (t(4;14), t(14;16), del(17p), gain(1q21)) and impaired clinical outcomes (2 de walker), especially in patients with accumulation of adverse cytogenetic lesions.  In our cohort we found a trend toward shorter PFS in high risk-patients and the tendency became clearer when TP53 mutations were additionally considered (p= 0.09 vs 0.06). Another author reported a correlation between cytogenetics aberration and mutated genes, such as t(11;14) with KRAS and IRF4 mutations. However, due to the lower number of samples and also, the low incidence of mutations in our cohort, we did not find any significant correlation.

We have identified mutation in the 85% of patients, half of them at <25 % of VRF. These values demonstrate the high ability of this methodology in comparison with WES to describe the subclonal landscape of MM. Although low frequency mutations did not impact in survival after a first line of treatment, we unknown their role in subsequent relapses as well and in patients responses to treatment. This approach allows to identify minority subclones at diagnosis and tracks them along patients disease with an affordable price and low sample requirement, essential properties for introducing personalize medicine protocols in clinical workout. The monitoring of subclonal mutations could be crucial for revealing new insight into molecular evolution of myeloma clones with treatment.