Spatiotemporal genomic architecture informs precision oncology in glioblastoma

Abstract

Precision medicine in cancer proposes that genomic characterization of tumors can inform personalized targeted therapies. However, this proposition is complicated by spatial and temporal heterogeneity. Here we study genomic and expression profiles across 127 multisector or longitudinal specimens from 52 individuals with glioblastoma (GBM). Using bulk and single-cell data, we find that samples from the same tumor mass share genomic and expression signatures, whereas geographically separated, multifocal tumors and/or long-term recurrent tumors are seeded from different clones. Chemical screening of patient-derived glioma cells (PDCs) shows that therapeutic response is associated with genetic similarity, and multifocal tumors that are enriched with PIK3CA mutations have a heterogeneous drug-response pattern. We show that targeting truncal events is more efficacious than targeting private events in reducing the tumor burden. In summary, this work demonstrates that evolutionary inference from integrated genomic analysis in multisector biopsies can inform targeted therapeutic interventions for patients with GBM.

Conflict of interest statement

Disclosure of Potential Conflicts of Interest

The authors declare no competing financial interests.

Figures

Figure 1. Mutational landscape of multi-region malignant glioma samples

(a) A schematic representation of glioma genomic heterogeneity and differential drug response analysis. Human glioma specimens were acquired based on their spatial order, or longitudinal pairs and subjected for genomic analysis for identification of tumor-initiating (truncal) events. (b) Somatic mutations including Single Nucleotide Variants (SNVs) and small Insertions/Deletions, copy number alterations, and gene fusions of 83 glioma multi-region or multisector-longitudinal specimens from 30 patients are demonstrated. 34 locally adjacent tumor fragments were from 14 patients, 13 multifocal/multicentric (referred as multiple) tissues from 5 patients, and a longitudinal pair GBM14 with leptomeningeal seeding were collected from Samsung Medical Center (SMC). We also curated 34 multisector-longitudinal tumor exomes and/or RNA sequencing from 10 patients in TCGA cohort. All somatic mutations called by SAVI with allele frequency >5% were demonstrated. For each gene we calculated the copy number (CN) based on Excavator. Clonal alterations were determined using ABSOLUTE with cancer cell fraction >80%. (Methods).

Figure 2. Comparison of genetic heterogeneity across glioma multisector / longitudinal samples

Patient samples were classified into five groups: Local, Multiple Lesion, S.T. (Short-Term) Longitudinal Local, L.T. (Long-Term) Longitudinal Local and Longitudinal Distant for comparative analyses. (a) Nei’s genetic distance of the indicated groups are shown. Q-values were calculated by Wilcoxon Rank-Sum Test and corrected for false discovery rate using Benjamini-Hochberg method. S.T. and L.T. Local indicates short-time (<18 months surgical interval) and long-time recurrent tumors ( ≥18 months), respectively. (b) Illustration of leave-one-out results from multinomial logistic regression. Each point indicates one pair of samples, and their coordinates are the probabilities to be local, multiple lesion/longitudinal distant, or longitudinal local. Long-Time recurrent samples were classified together with multiple lesion/longitudinal distant samples, indicating they might follow the same evolutionary model. (c) Tumor evolution behind Big Bang and Multiverse models. The Big Bang model is indicated as a mixture of tumor cells that share many clonal and subclonal alterations. The Multiverse model is indicated with a greater proportion of private events at a clonal level. (d) Pie charts demonstrate the frequencies of PIK3CA mutations in multifocal/multicentric (M-) GBMs (30%, 9/30) and solitary (S-) GBMs (10%, 13/130). The p-value was calculated using Fisher’s exact test.

Fig. 3. Single cell transcriptome from multi-region samples

(a) Expression profile of individual tumor cells from three samples of GBM9 (left initial, right initial, and relapse), according to expression subtypes. For each cell, the subtype with the highest expression is marked with an asterisk. Several EGFR genomic alterations can be identified in the single cell expression data (yellow), despite the abundance of missing data (gray). (b) Topological representation of the expression data of individual tumor cells from GBM9, labeled by sample of origin. Each node represents a set of cells with similar transcriptional profile. A cell can appear in several nodes, and two nodes are connected by an edge if they have at least one cell in common. Topological representation of GBM9 labeled by expression of EGFR(c) and mitotic genes (d). (e) Expression profile of individual tumor cells from GBM10 (two samples: 5-ALA (+) and 5-ALA (-)). The p-value between proneural and 5-ALA was obtained using Fisher’s exact test. The p-value between mesenchymal and LTBP4 expression was calculated based on Spearman’s correlation. ATRX fusion was validated by RT-PCR assays (Supplementary Fig. 7b–c). Topological representation of expression data of individual tumor cells from patients GBM10 (f) and GBM2 (g), labeled by the sample of origin. (h) Expression profile of individual tumor cells from patient GBM2 according to the GBM expression subtypes.

Fig 4. Chemical screening of multi-region patient-derived cells (PDCs)

(a) PDCs were treated with 40 chemical agents, targeting oncogenic signaling pathways in diluted series from 20 μM to 4.88 nM. X-axis indicates Nei’s genetic distances between fragments from the same patient, while y-axis indicates Spearman’s correlation coefficient (SCC) of corresponding fragments based on drug sensitivities measured by Area Under the Curve (AUC). (b) A violin plot for SCC of drug responses of the groups described in Fig. 4a. (c) Mean values of AUCs for six PI3K/AKT/mTOR (PAM) inhibitors (BEZ235, BKM120, BYL719, AZD5363, AZD2014 and Everolimus) of PDCs isolated from M- (n=9) or S-GBMs (n=22) were plotted. (d) The normalized Z-score in each PDC was plotted when the corresponding tissues harbored associated genetic alterations, designated as “shared” or “private”. The private group was determined when the drug response-associated genetic alteration (i.e. EGFR mutations-EGFR inhibitors; PTEN mutations-PI3K/AKT/mTOR pathway inhibitors) was private, and vice versa for the shared group. (e) Preoperative T1-weighted contrast-enhanced MR image and key genomic alterations found in the corresponding tumors and its derivative cells from a multicentric patient (GBM9). Right-side ‘R’ tumors encompassed the right frontal lobes and corpus callosum (CC). ‘L’ indicates the left frontal tumor. Preoperative MRI showed a multifocal infiltrative lesion in both the frontal lobe and CC. (f) Scatterplot of AUC for 40 cancer-targeting compounds on GBM9 PDCs derived from the left and right side tumor. The R was obtained as Pearson’s correlation coefficient. All p-values in this figure were obtained using Wilcoxson Rank-Sum test.