Machine Learning Sequence Modeling Identifies Gene Regulatory Responses to Bone Marrow Stromal Interactions in Multiple Myeloma

Introduction Multiple myeloma (MM) is a hematological malignancy of plasma cells which affects different organs including the bone marrow. The interaction between MM and the bone marrow stromal cells (BMSC) has been shown to affect the disease progression and response to treatment (Kumar et al. Multiple myeloma. Nat. Rev. Dis. Primers, 2017). An in vitro system to isolate the signaling mechanisms between BMSC and MM cells has been developed by (Dziadowicz et al., Cancers, 2022), who co-cultured MM cell lines in transwell media (TSW) exposed to BMSC secreted factors and compared their ATAC-seq profiles to MM cell lines (MONO) in isolation without the exposure to BMSC. While some of the gene expression changes induced by the BMSC-MM crosstalk were reported in (Dziadowicz et al., Cancers, 2022), here we model the entire chromatin accessibility profiles with a DNA-sequence based machine learning model to detect the altered Transcription Factor (TF) and enhancer activity responsible for driving the gene expression changes. We used sequence-based machine learning models (Beer et al., 2020, Annu. Rev. Genom. Hum. Genet.) and applied a systematic approach to model the chromatin accessibility changes induced by the interaction of BMSC and MM cells (MM-TSW) vs. MM cell lines cultured in isolation (MM-MONO). Methods We trained gapped-kmer SVM (gkm-SVM) machine learning models on MM1S and RPMI8226 MM cell line ATAC-seq data generated by Dziadowicz et al. as well as 1270 DNase-seq chromatin accessibility profiles of various human primary tissues and cell lines (including MM1S and RPMI8226) from the ENCODE consortium. We trained gkm-SVM on the top 2000 differentially-accessible distal ATAC peaks (>2000bp from transcription start site, e.g. enhancers) to find the differential TF activity between MM-TSW and MM-MONO (Figure 1) and found the TF binding site motifs explaining the gkmSVM output kmer weight distribution (Figure 2). We also identified the set of differentially-expressed TFs and performed gene set enrichment using MSigDB. Results We found that RUNX and Ebox-binding (5‘-CACCTG-3‘) TF family members are the most active TFs in MM1S-MONO and RPMI8226-MONO enhancers (gkmSVM AUROCs~0.92). The contribution of the Ebox TFs to distal chromatin accessibility in the MM cell lines is greater than in any other ENCODE sample (n=1270). Principal Component Analysis (PCA) on distal ATAC peaks compared to promoters showed that the MM transcriptional response to BMSC is more clear