ABSTRACT
Selinexor is the first selective inhibitor of nuclear export (SINE) to be approved for treatment of relapsed or refractory multiple myeloma (MM). There are currently no known genomic biomarkers or assays to help select MM patients at higher likelihood of response to selinexor. Here, we aim to characterize transcriptomic correlates of response to selinexor-based therapy, and present a novel, three-gene expression signature that predicts selinexor response in MM. We analyzed RNA sequencing of CD138+ tumor cells from bone marrow of 100 MM patients who participated in the BOSTON study and correlation with clinical outcomes. We validated this gene signature in 64 patients from the STORM cohort of triple-class refractory MM, and additionally in an external cohort of 35 patients treated in a real world setting outside of clinical trials. We also found that the signature tracked with response in a cohort of 57 patients with recurrent glioblastoma treated with selinexor. Furthermore, the genes involved in the signature, WNT10A, DUSP1, and ETV7, reveal a potential mechanism through upregulated interferon-mediated apoptotic signaling that may prime tumors to respond to selinexor-based therapy. This signature has important clinical relevance as it could identify cancer patients that are most likely to benefit from treatment with selinexor-based therapy.
INTRODUCTION
Selinexor is the first selective inhibitor of nuclear export (SINE) approved for treatment of relapsed or refractory multiple myeloma (MM) and diffuse large b-cell lymphoma (DLBCL).1, 2 This approval is well supported by recent clinical trial data, most notably the STORM and BOSTON studies. In the STORM phase II clinical trial, oral selinexor and dexamethasone were administered twice weekly in patients with triple-class refractory MM, with 26% of patients exhibiting a partial response (PR) or better, and a median progression-free survival (PFS) at 3.7 months.1, 3 The phase III BOSTON trial compared the efficacy of once weekly selinexor in combination with once weekly bortezomib and dexamethasone with standard twice weekly bortezomib and dexamethasone in patients with previously treated MM and found an overall response rate (ORR) of 76.4%, and a median PFS of 13.9 months for patients receiving the novel treatment regimen, in comparison to 9.5 months for those receiving the standard treatment.1
Mechanistically, selinexor binds to the karyopherin exportin 1 (XPO1), which is responsible for shuttling more than 200 oncogenic and tumor suppressor proteins and mRNA transcripts to the cytoplasm.4 This inhibition of nuclear export, the sequestration of tumor suppressor proteins in the nucleus, and the prevention of select oncogene mRNA translation into oncoproteins ultimately induces cancer cell death while permitting the survival of non-malignant cells.5–7 Selinexor has a unique adverse event profile due to its novel mechanism of action, and adverse events (AE) are common, with severe AE reported in 52% of patients in the BOSTON trial, and 80% of patients in the STORM trial.1, 3 Identifying biomarkers that help predict treatment responses and toxicity is essential to targeted selinexor-based therapeutic intervention.
Systematic approaches to biomarker discovery for selinexor response that leverage next-generation sequencing are generally lacking in the literature. While MM cells tend to over-express XPO1 compared to normal plasma cells, XPO1 alterations have not correlated significantly with response to treatment with SINE compounds with the exception of a few in-vitro studies.4–6, 8 The STORM study presented a brief transcriptomic analysis that identified a potential four-gene signature based on imputed protein activity.3 A few candidate biomarkers identified through bioinformatics analyses have also been reported in conference proceedings.9, 10 However, limitations of these studies include small sample sizes and lack rigorous validation with external cohorts of selinexor-treated patients.
Here, we rigorously analyzed RNA sequencing data from 213 patients from multiple studies of selinexor-based therapy to characterize transcriptomic correlates of response to selinexor. We discovered a novel three-gene signature predicting response to selinexor using RNA-seq data from the BOSTON study and validated it on data from the STORM clinical trial and on an additional independent cohort of patients treated with selinexor at the Mount Sinai Hospital in “real-world” conditions outside of a clinical trial. The three-gene signature is biologically interpretable and opens a path for evaluating a mechanism of response in future studies. More importantly, this signature has the potential to identify patients most likely to benefit from treatment with selinexor-based therapy, ultimately reducing toxicities and improving outcomes.
RESULTS
Patient characteristics and transcriptomic profiling of selinexor-treated MM tumors
BOSTON
We performed RNA-seq on CD138+ bone marrow plasma cells from 100 patients who participated in the BOSTON study (Table 1). We included 53 samples from the selinexor, bortezomib, and dexamethasone (XVd) arm of the study, and 47 samples from the bortezomib and dexamethasone (Vd) arm. In the XVd arm, most patients were male (60.3%), compared with 46.8% in the Vd arm. In both study arms, and across the entire cohort, the median age was 67 years. The ORR in the XVd arm, defined as partial response (PR) or better, was 75.47%, compared to an ORR of 65.95% in the Vd Arm (Table 1). These observations suggest that the addition of selinexor to a regimen of bortezomib and dexamethasone in MM offers a clinical benefit over bortezomib and dexamethasone alone, and are consistent with the findings of previous studies, including the main BOSTON trial.1
Cohort Clinical Characteristics.
For validation purposes, RNA-sequencing on CD138+ selected cells from two external cohorts were leveraged: MM patients treated with selinexor-dexamethasone who participated in the STORM trial (N = 64), and MM patients treated post-FDA approval at Mount SInai Hospital outside of a clinical trial setting (MSSM cohort, N = 35).3 Patients in the STORM cohort had failed at least five prior lines of therapy.3 The median PFS was 62 days, and 18 patients (28.1%) achieved a response of PR or better (Table 1). Patients in the MSSM cohort received a median of 7 lines of therapy (range = 1-18) prior to starting selinexor based therapy, which was often administered in combination with a variety of other agents. The median PFS was 66 days, and 10 patients (28.5%) achieved a response of PR or better (Table 1). The most common combination drug strategy in addition to the selinexor backbone was carfilzomib and dexamethasone, used in 7 of 35 (20%) patients.
We performed principal component analysis (PCA) on the normalized, batch-corrected expression matrix to understand sources of variation in the BOSTON, STORM, and MSSM datasets, which did not identify any specific bias (Fig S1). Overall, we did not find any significant or precise predictors of response to selinexor from clinical, demographic, or cytogenetic markers.
Differential expression analysis identifies genes associated with selinexor response
The strategy used for differential expression (DE) analysis is summarized in Fig 1A. To identify genes whose expression is associated with longer PFS or better depth of response to selinexor, we performed a series of DE comparisons across 9 unique different PFS or depth of response (defined according to International Myeloma Working Group [IMWG] criteria) thresholds in both the XVd and Vd arms (Fig 1A-B).11 Definitions of each BOSTON comparison are detailed in Table S1. Across all response thresholds, there were a total of 107 unique significant DE genes between better and worse responders in the XVd Arm (27 upregulated, 70 downregulated, FDR < 0.05) and 560 unique DE genes in the Vd Arm (398 upregulated, 162 downregulated, FDR < 0.05). To identify genes whose association with PFS- or IMWG response category in the XVd arm was not also associated with bortezomib/dexamethasone therapy, we retained for further analysis genes that were significant in the XVd arm but not in the Vd arm for each corresponding cutoff as shown in Fig 2B. In total, there was a moderate overlap (up to 31%) across genes identified through the different comparisons, with 6 out of 24 (25.0%) uniquely downregulated genes and 12 of 33 (36.26%) uniquely upregulated genes overlapping across at least two different cutoffs.
Differentially expressed genes in selinexor responders. A) Experimental design for identifying DE genes specific to selinexor. B) Bar plot indicating number of selinexor specific DE genes across various cutoffs of PFS, depth of response, or both. C) Volcano plot of the BOSTON XVd DE comparison with a PFS ≥ 120 days response cutoff showing the three significantly upregulated genes used to generate the three-gene response signature.
Upregulation of the signature is associated with longer PFS in pre-selinexor-treated MM. A) Kaplan-Meier curve of low versus high expression of the three gene GSVA signature in the XVd arm of the BOSTON cohort shows a significant association with longer PFS in patients with upregulation of the signature. B) Using the Vd arm as a negative control, the three gene signature is not associated with PFS. C) Three gene signature validates in the STORM cohort of triple-class refractory MM patients who received selinexor as a single-agent treatment. D) The signature also validates in 35 patients that received selinexor in a retrospective cohort at MSSM that is not part of a clinical trial and is a more heterogeneous patient population reflective of “real world” treatment settings. E) Signature does not validate in non-selinexor, standard of care treated MMRF-COMMPASS samples.
A three-gene signature predicts selinexor response
Next, we identified a gene expression signature for selinexor response. Each DE gene set identified through the 9 unique comparisons was split into upregulated and downregulated subsets based on logFC in the XVd Arm, resulting in 17 candidate gene signatures for further analysis (9 upregulated, 8 downregulated). We performed Gene Set Variation Analysis (GSVA) on the normalized and batch-corrected expression matrix to calculate a unique score for each of the 17 candidate gene signatures. For each signature, a univariate proportional hazard model was generated in the XVd arm. The gene signature with the best performance in model training was identified based on a ranking procedure using repeated four-fold cross validation and the spearman correlation of the signature with PFS (Fig 2A; see Methods).
The best performing signature was composed of three genes, WNT10A, DUSP1, and ETV7. It was correlated with PFS (Spearman rho = 0.46, P = 0.0007) and was upregulated in XVd patients with PFS ≥ 120 days. The signature was predictive when using a proportional hazard model (FDR=0.047, HR=0.36 [95% CI = 0.14-0.84]; log-rank P = 0.017, Fig 2A). To further ensure that the predictive effect of the signature was not due to random variations, we also performed a permutation test and found that the signature was more predictive than a GSVA score composed of three randomly selected genes (10,000 permutations, P = 0.03; see methods). In patients who achieved a partial response (PR) or better, higher signature expression was significantly associated with longer duration of response (DOR), defined as the number of days from IMWG partial response to progression (Fig S2A). Furthermore, we also found that the signature was significantly higher in patients who achieved IMWG response category of very good partial response (VGPR) or better, suggesting that the signature is associated with both duration and depth of response (Wilcoxon P = 0.014, Fig 3A-C).
Higher signature expression is associated with better depth of response. A) IMWG response category of the BOSTON XVd arm versus three-gene signature score on the y-axis, B) Wilcoxon test of signature expression in BOSTON XVd responders defined as a PR or better, and C) as VGPR or better. D) STORM Xd signature expression across IMWG best overall response categories, as well as comparison with E) expression in responders versus non-responders defined as PR or better, and F) as VGPR or better. G) MSSM non-trial signature expression across IMWG best overall response categories, as well as comparison with H) expression in responders versus non-responders defined as PR or better, and I) as VGPR or better.
The same analysis was applied to the Vd arm as a negative control under the rationale that a signature specific to selinexor response would not accurately distinguish long or short PFS in cohorts treated with non-selinexor based therapy. The three-gene signature did not track with PFS or DOR in the Vd arm, as shown in Fig 2B and Fig S2B.
The three-gene signature is validated in external cohorts
STORM
A log-rank test comparing signature expression higher or lower than the cutoff of zero in the STORM dataset also validated the finding from BOSTON that higher signature expression is associated with PFS (log rank P = 0.039; N=64). We found that the linear association with PFS and the signature performed nominally well, despite not reaching statistical significance (N=64; P = 0.056, HR=0.621 [95% CI = 0.483-0.389], Fig 2C). The poorer signature performance in a cox regression may be explained by lower PFS, as patients in the STORM cohort were triple-class refractory with more advanced disease to begin with. We also tested whether the association between signature expression and depth of response also validated in the STORM cohort and found that patients with a VGPR or better had significantly higher expression of the signature than those with a worse depth of response (Wilcoxon P = 0.035, Fig 3D-F).
MSSM
Using cox regression, we found that higher expression of the three-gene signature was significantly associated with survival (P = 0.0145, HR = 0.41 [95% CI = -0.467-0.551]). This result was also replicated via log rank testing (P = 0.0023, Fig 2D). Furthermore, we found a statistically significant correlation of the gene expression signature with PFS via Spearman correlation (rho = 0.4, P = 0.01). When compared with depth of response, we found that higher expression of the signature trends towards a better IMWG response category (Fig 3G-I).
MMRF-CoMMpass
We additionally used the MMRF-CoMMpass dataset (N=767) as a negative control and found that the signature was not predictive of PFS in patients who were treated with non-selinexor based, standard of care therapies. The signature did not show a significant association with progression free survival via cox regression (P = 0.14, HR = 0.89 [95% CI = 0.73-1.04]) or via log rank test (Wilcoxon P = 0.32, Fig 2E). Taken together, these results suggest that the signature is specific to selinexor treatment response and is not reflective of overall prognosis.
KING Study in Recurrent Glioblastoma
Lastly, we sought to test whether the signature would be predictive in cohorts of patients with other cancer types treated with selinexor from other cancer types. To test this hypothesis, we obtained RNA-seq of tumor samples from 57 patients with recurrent glioblastoma who were treated with selinexor monotherapy as part of the phase II KING trial (NCT01986348).12 Overall, we found that patients with higher expression of the signature experienced improved PFS, although statistical significance was not achieved (log rank P = 0.078, Cox PH P = 0.0734, HR = 0.549, Fig S3A). However, patients who achieved a clinical benefit with a partial response (PR) or better had significantly higher signature expression (Wilcoxon P = 0.0034, Fig S3B).
Gene-set enrichment analysis reveals response-associated activation in wnt, apoptosis, and mapk signaling pathways
We next identified pathways associated with response to selinexor in the BOSTON XVd cohort using gene set enrichment analysis (GSEA) on the MSigDB Hallmark gene sets.13 Notably, we found downregulation of MYC targets, as well as upregulation in KRAS, apoptotic, and interferon signaling among significantly enriched pathways in XVd responders (Fig S4).
Knockdown of signature genes in-vitro induces resistance in previously sensitive cells
To further validate the signature genes, we generated siRNA knockdowns of ETV7, WNT10A, and DUSP1 in vitro in JJN3 cell lines with two replicates (Fig 4). In JJN3 cells, the ETV7 knockdown showed a mean 22.7% (N=2) increase in percent cell viability compared to non-target control after 72h of treatment with 100nM of selinexor (Fig 4B). These results demonstrate that knockdown of ETV7 can cause previously sensitive MM cells to achieve resistance against selinexor.
Knockdown of signature genes in sensitive MM cell lines. A) qPCR efficacy of siRNA knockdown of ETV7, WNT10A, and DUSP1 shown as the mean of two replicates in JJN3 cells. B) Percent cell viability of JJN3 cells for each knockout condition after 72H of treatment with 100nM of selinexor.
METHODS
Patient Selection & RNA Sequencing
BOSTON and STORM
CD138+ cells were purified from bone marrow aspirates obtained from 100 patients who participated in the BOSTON study and 64 patients who participated in the STORM study. RNA extraction was performed using Qiagen Allprep RNA mini kit and library preparation was performed with the TruSeq Stranded mRNA (non-FFPE compatible) kit. For samples where the RNA quality was low, the Smart-Seq V4 Ultra Low Input Nextera XT kit was used. Total RNA sequencing was performed with 100 bp reads using an Illumina HiSeq 2500 instrument.
MSSM
The 35 patients were physician-referred as part of the multiple myeloma banking protocol approved by the Mount Sinai institutional review board (IRB). Informed written consent was obtained from each patient. CD138+ cells were isolated from bone marrow aspirates obtained prior to the start of selinexor-based therapy, and RNA isolation and sequencing was performed as previously described.14
MMRF-CoMMpass
Gene-level counts for 767 RNA-seq samples from the MMRF-CoMMpass dataset were downloaded from dbGaP (accession #phs00748).
KING
Gene-level read counts were obtained from pretreatment tumors from 57 patients with recurrent glioblastoma that were enrolled in the phase 2 KING study.12
Bioinformatics Processing
Raw reads were aligned to the GRCh38 human reference genome using STAR.15 Gene-level counts, obtained through featureCounts, were filtered to remove immunoglobulin and ribosomal transcripts, and to remove genes whose counts across all samples had zero variance.16, 17 Counts were converted to Log2CPM, normalized with voom, and corrected for batch effects or covariates identified through variancePartition analysis using the sva ComBat package.17–19
Differential Expression
DE genes were identified in the BOSTON XVd and the Vd arm separately using DESeq2.20 Within each arm, DE genes relative to responders versus non-responders were generated across a total of 9 response cutoffs based on PFS, ORR, or a combination of both. For each response comparison, genes were designated as DE in selinexor responders if they were significantly DE in the XVd Arm (PAdj < 0.05) and were not significantly DE in the corresponding comparison within the Vd arm (P > 0.05). DESeq2 analysis was performed on unfiltered raw counts per the tool’s requirements.20 Gene set enrichment analysis was applied to the selinexor-specific DESeq2 results in the XVd arm using the curated Hallmark and Reactome signatures publicly available from MSigDB using the fGSEA package.13, 21
Survival Analysis
Each DE gene set identified through the 9 comparisons was split into upregulated and downregulated subsets based on logFC in the XVd Arm, resulting in 17 candidate gene signatures for further analysis (9 upregulated, 8 downregulated). Gene set variation analysis (GSVA) scores were calculated for all samples for the 17 candidate gene signatures from the covariate-normalized expression matrix.22 For each candidate gene signature, a univariate Cox proportional hazard model was generated in the BOSTON XVd arm. The gene signature with the best model performance in the discovery cohort was selected by ranking the highest Somer’s Dxy after repeated four-fold cross validation (nrepeats = 1000). The gene signature was also evaluated with a Cox proportional hazard model in the Vd arm as a negative control to ensure that performance was specific to selinexor. A permutation test was performed to evaluate whether the selected gene signature is more significant via log-rank Kaplan-Meier testing than a GSVA score composed of three randomly selected genes.
All validation tests were executed by first performing quantile-normalization of the expression matrix with the distribution of expression in the BOSTON dataset using feature-specific quantile normalization, followed by calculation of a GSVA score for the gene signature.22, 23 Survival was tested between groups that had low expression versus high expression, using a cut-off of zero, with a log-rank test and univariate cox regression carried out with the same procedures used for the discovery cohort.
siRNA Knockdown Validation
Cell lines, drug, and culture conditions JJN3 and MM1.S cell lines were cultured in RPMI1640 medium (Mod.) 1X with l-Glutamine (Corning), supplemented with 10% FBS (Gemini Bio-Products) and 1% penicillin-streptomycin 100× solution (Corning). Selinexor (KPT-330, cat #S7252) was purchased from Selleck Chemicals (storage at −20°C). All cells were propagated in standard cell culture conditions (5% CO2, 37°C) in cell culture–treated T75/T150 flasks (Falcon). All cell lines were authenticated by short tandem repeat DNA profiling and routinely tested for Mycoplasma infection (MycoAlert, Lonza Bioscience).
siRNA knockdown and Selinexor Treatment
For transfections, ON-TARGETplus SMARTpool siRNA construct targeting ETV7 (#L-017938-00-0005), WNT10A (#L-008232-00-0005) and DUSP1 (#L-003484-02-0005) and a control non-targeting (NT) siRNA no. 1 (#D-001810-01-20; all siRNA from Dharmacon, Horizon discovery) were used. JJN3 cells were transiently transfected with siRNA for ETV7, WNT10A, DUSP1, and NT control siRNA using the cell line 4D Nucleofector Solution SG program EO100 (Lonza, Basel, Switzerland). Cells were then treated with 5nM, 10nM, 50nM and 100nM Selinexor 24 hours after transfection, followed by harvesting and analysis using quantitative RT-PCR and cell viability assays.
Cell Viability
Cell viability was assessed using fluorometric resazurin reduction method (CellTiter-Blue; Promega) according to the manufacturer’s protocol. After 24h of siRNA knockdown, JJN3 cells (100ul; 1 × 10 5 cells/well) were cultured (quadruplet per condition) in a round-bottomed 96-well microtiter plates without and with increasing concentration of Selinexor for 72 hours in standard cell culture conditions (5% CO2, 37°C). Then, the number of viable cells in each treated well was calculated by adding 20 uL of CellTiter-Blue Reagent to each well and further incubation with the dye for 1 hour. The fluorescence (560Ex/590Em) was measured with the FLUOstar Microplate Reader (BMG LABTECH). The linear least-squares regression of the standard curve was used to calculate the number of viable cells in each treated well. DMSO treated cells were used to normalize the cell viability in drug treated cells. Cell counts were confirmed on the Countess Automated Cell Counter (Life Technologies). Viability was calculated as % out of DMSO control (Number of live cells in Selinexor treatment/number of live cells in DMSO control x 100 %) in each transfection group.
Quantitative RT-PCR
To assess the knockdown efficacy after 24 and 48-hour siRNA transfection, total RNA was purified from 1 million cells using the RNeasy Plus Mini Kit (Qiagen) per manufacturer protocol. cDNA was prepared from mRNA using SuperScript VILO Master Mix (Thermo Fisher Scientific) and quantified using SsoFast EvaGreen Supermix (Bio-Rad Laboratories) on a CFX96 Touch Real-Time PCR Detection System (Bio-Rad Laboratories).
Samples and controls were run in triplicate. Gene expression was normalized to hypoxanthine phosphoribosyltransferase (HPRT) and expressed relative to cells transfected with control non-target siRNA using the 2−ΔΔCT formula. Thermal cycler conditions were: initial step of 30 seconds at 95°C followed by 40 cycles of 5 seconds at 95°C (denature) and 5 seconds at 60°C (anneal/extend) followed by 5 minutes at 55°C to 95°C in increments of 0.5°C. Primers used for quantitative RT-PCR are listed below:
HPRT Forward Primer: 5′-AAAGGACCCCACGAAGTGTT-3′
HPRT Reverse Primer: 5′-TCAAGGGCATATCCTACAACAA-3′
ETV7 Forward Primer: 5′- CTGCTGTGGGATTACGTGTATC -3′
ETV7 Reverse Primer: 5′- GTTCTTGTGATTTCCCCAGAGTC-3′
DUSP1 Forward Primer: 5′- AGTACCCCACTCTACGATCAGG -3′
DUSP1 Reverse Primer: 5′- GAAGCGTGATACGCACTGC -3′
WNT10A Forward Primer: 5′- GGAGACTCGCAACAAGATCCC -3′
WNT10A Reverse Primer: 5′- CGATGGCGTAGGCAAAAGC -3′
DISCUSSION
Selinexor is approved as a second-line therapy for MM and its efficacy is well supported by clinical trials. However, there are currently no known biomarkers to better guide selection of patients whose tumors are more sensitive to selinexor-based therapy. Furthermore, while the mechanism of action is well characterized for selinexor, little is known about the correlates of response or resistance to selinexor-based therapy. Here, we describe the transcriptomic characteristics of patients who respond to selinexor therapy. Further, we report the discovery of a robustly validated three-gene signature that is predictive of response to selinexor-based therapy in MM.
There is little literature to date on patient populations describing correlates to response in the context of selinexor therapy. Some studies have shown candidate biomarker activity in miRNAs as regulators of XPO1 and its targets, and certain mutations in XPO1 and XPO5 as prognostic markers for survival, but they have not been correlated to SINE drug sensitivity.8, 24–27 There are few studies that have explored biomarkers in selinexor therapy, and fewer with validation in external cohorts. One notable study found a signature based on imputed activity comprising four master regulator proteins, IRF3, ARL2BP, ZBTB17, and ATRX, but did not provide any external cohort validation data.3
We report here a gene expression signature composed of the combined upregulation of three genes, ETV7, WNT10A, and DUSP1, that precisely and accurately predicts both depth and duration of response to selinexor-based therapy in patients with MM. Furthermore, we present robust external validation and negative controls. Finally, the signature genes are validated in-vitro. This is the first robustly validated signature for selinexor response to date. Given that more than 50% of patients treated with selinexor in the BOSTON and STORM studies had adverse events, the discovery of a predictive gene signature holds tremendous potential for biomarker-guided selection of candidates for selinexor-based therapy. Since the signature was validated in multiple heterogeneous patient cohorts, including selinexor-dexamethasone monotherapy (STORM), and in various combination drug scenarios outside of clinical trials (MSSM), it is both flexible and applicable to a wide variety of real-world scenarios where selinexor may be given in combination with other drugs. One limitation of the signature is that, since it relies on a GSVA score, its accuracy is dependent on a cohort with multiple samples and prediction improves with greater sample sizes. However, since it is composed of just three genes, it also holds potential for fast and simple implementation in clinical settings, potentially via qPCR-based assay development strategies or in combination with ex-vivo drug sensitivity assays. Furthermore, since the signature uses gene expression directly, it is more interpretable and easier to implement than previously cited signatures, which used imputed protein activity and are not as thoroughly validated.3
Interferon signaling is responsible for the anti-viral immune response and has anti-proliferative properties. It has been shown to play an important role in apoptotic signaling in MM, mostly through IRF4, MYC, and BCL6.28–30 Here, we found a strong enrichment for upregulated interferon signaling and genes involved in interferon signaling in patients who responded to XVd therapy. We also found a corresponding enrichment of dysregulated MYC and apoptotic signaling. Interferon signaling has been found to modulate response to XPO1 inhibition by eltanexor to treat viral infection.31 Additionally, ETV7 whose depletion with RNA interference induced resistance in in vitro experiments, has been identified as an interferon-stimulated gene.32 Based on these results and existing literature, we hypothesize that there may be a potential biological rationale for further evaluating the role of interferon signaling in MM selinexor response.
In conclusion, we report a novel gene expression signature for response to selinexor-based therapy in patients with MM. We have validated our findings in several external transcriptomic datasets of MM patients treated with selinexor-based regimens. Ongoing in vitro and mechanistic studies will help determine whether they are causative of response or simply a correlative readout of some other selinexor response mechanism. This signature has important clinical significance as it could identify patients most likely to benefit from treatment with selinexor-based therapy, especially in earlier lines of therapy.
Data Availability
Raw sequencing data will be deposited in SRA upon publication. The code necessary to calculate selinexor signatures is provided in a GitHub repository at https://github.com/Lagana-Lab/selinescore.
AUTHOR CONTRIBUTIONS
PR performed most of the downstream bioinformatics analysis, with assistance from SB and AL. PR, DTM, and DB processed the raw sequencing data. SA, AA, and VL, compiled clinical data for analysis. YGP performed the cell line experiments with assistance from AA and VL. CW and YL provided data from the KING study. AC, SJ, SR, SP, DM, JR, HJC, YL provided patient samples. VL performed library preparation for sequencing. PR wrote the initial manuscript with assistance from JJ and YGP. AL, SP, CW and YL conceptualized and supervised the study.
CONFLICT OF INTEREST
This study was supported by Karyopharm Therapeutics. YL reports being employed by and owning stock in Karyopharm Therapeutics. CW reports being employed by Karyopharm Therapeutics. AC reports receiving grant support and consulting fees from Millennium/Takeda, grant support, advisory board fees, and consulting fees from Celgene, Novartis Pharmaceuticals, Amgen, and Janssen, consulting fees from Bristol-Myers Squibb, advisory board fees from Sanofi and Oncopeptides, grant support from Pharmacyclics, and grant support and advisory board fees from Seattle Genetics. DM reports current employment by Janssen. HJC reports research support from Bristol-Myers Squibb and Takeda. JR reports receiving consulting fees and fees for serving on a speakers bureau from Amgen, advisory board fees and fees for serving on a speakers bureau from Celgene, Takeda, and Janssen, and advisory board fees from Sanofi, Karyopharm Therapeutics, Oncopeptides, Adaptive Biotechnologies, and Bristol-Myers Squibb. SR reports receiving consulting fees from Karyopharm Therapeutics, Janssen Pharmaceuticals, and Bristol-Myers Squibb. SJ reports receiving advisory board fees and consulting fees from Celgene, Bristol-Myers Squibb, Janssen Pharmaceuticals, and Merck. SP reports receiving grant support from Karyopharm Therapeutics.
DATA AND CODE AVAILABILITY
Raw sequencing data will be deposited in SRA upon publication. The code necessary to calculate selinexor signatures is provided in a GitHub repository at https://github.com/Lagana-Lab/selinescore.
SUPPLEMENTARY FIGURES
A) PCA of BOSTON cohort normalized expression matrix colored by patient PFS in days. B) Spearman correlation of principal components with clinical variables, annotated with P values for correlations with P < 0.05, in BOSTON XVd and C) BOSTON Vd. D) PCA of STORM cohort normalized expression matrix, and E) MSSM cohort samples colored by patient PFS in days.
A) Kaplan Meier curve of higher versus lower signature expression as a function of duration of response (DOR) in patients who achieved a partial response (PR) or better in A) the BOSTON XVd cohort shows a significant association, whereas there is no significant difference in B) patients in the BOSTON Vd cohort of non-selinexor based therapy.
A) Survival analysis of the KING recurrent glioblastoma cohort stratified by three-gene signature GSVA scores. B) Three gene signature GSVA scores in the KING cohort in responders (purple) versus non-responders (gray).
GSEA enrichments in XVd responders across the various response cutoffs.