Inter- and intra-tumor heterogeneity of metastatic prostate cancer determined by digital spatial gene expression profiling

Abstract

Metastatic prostate cancer (mPC) comprises a spectrum of diverse phenotypes. However, the extent of inter- and intra-tumor heterogeneity is not established. Here we use digital spatial profiling (DSP) technology to quantitate transcript and protein abundance in spatially-distinct regions of mPCs. By assessing multiple discrete areas across multiple metastases, we find a high level of intra-patient homogeneity with respect to tumor phenotype. However, there are notable exceptions including tumors comprised of regions with high and low androgen receptor (AR) and neuroendocrine activity. While the vast majority of metastases examined are devoid of significant inflammatory infiltrates and lack PD1, PD-L1 and CTLA4, the B7-H3/CD276 immune checkpoint protein is highly expressed, particularly in mPCs with high AR activity. Our results demonstrate the utility of DSP for accurately classifying tumor phenotype, assessing tumor heterogeneity, and identifying aspects of tumor biology involving the immunological composition of metastases.

Fig. 1. Digital spatial profiling of archived formalin-fixed paraffin-embedded prostate cancer metastases.

a Schematic of study participants. N = 27 patients from the UW Rapid Autopsy Program were selected with two sites of metastasis per patient (N = 56). One 500 µM region of interest (ROI) per core was selected for DSP. Of the total 168 ROIs, 141 were utilized for analysis. Reasons for exclusion include missing/100% fat (N = 14), 100% stroma (N = 1), quality control failure (N = 4), or ≤100 genes detected (N = 8). b Schematic of multianalyte DSP workflow. Serial sections of the TMA were run through GeoMx RNA (top) or GeoMx protein (bottom) assays. Both assays were read out by next-generation sequencing. c Split violin plot overlaid on scatterplot of RNA assay counts by patient (N = 3-12 ROIs for most patients). Gene counts are displayed on the left side (green) and the geometric mean of negative probe counts is displayed on the right side (gray). d Histogram of the number of genes above the limit of quantitation (LOQ) per ROI for the RNA assay. LOQ was defined as the geometric mean of the negative probes × geometric standard deviation of negative probes squared. e Split violin plot overlaid on scatterplot of protein assay counts by patient (N = 3-12 ROIs for most patients). Antibody counts are displayed on the left side (green) and the geometric mean of negative control antibodies (mouse and rabbits IgGs) displayed on the right side (gray). f Histogram of the number of proteins above the limit of quantitation (LOQ) per ROI. LOQ was defined as three times the geometric mean of the negative control antibodies.

Fig. 2. DSP classifies mPC subtypes and quantitates the expression of therapeutic targets.

a Heatmap of DSP gene expression correlated with bulk RNAseq across androgen receptor (AR), neuroendocrine (NE), cell cycle progression (CCP), and FGFR/MEK gene signatures (N = 141 ROIs averaged to 53 tumors from 26 patients). Results are expressed as mean gene signature Z-scores and mean log₂ negative-normalized (NN) gene expression and presented according to color scales. RNAseq class and DSP class are the phenotypes assigned to the samples using each dataset. b–e Scatterplots comparing gene expression of RNAseq GSVA scores to mean DSP Z-scores across AR, NE, CCP, and FGRF-MEPK gene signatures (N = 141 ROIs averaged to 53 tumors from 26 patients.) Two-sided test for association using Pearson’s correlation coefficient, r; p value shown on plots. f Multidimensional scaling (MDS) plot of mCRPC phenotypes as defined by DSP using the mean DSP log₂ negative-normalized expression of 23 AR and NE genes (N = 141 ROIs averaged to 53 tumors from 26 patients.). g–j Scatterplot comparison of single genes EZH2, FOLH1, WEE1, and BCL2 mean DSP log₂ negative-normalized expression (N = 141 ROIs averaged to 53 tumors from 26 patients) vs RNAseq log₂ FPKM (N = 53 tumors from 26 patients.) Two-sided test for association using Pearson’s correlation coefficient, r; p value shown on plots.

Fig. 3. DSP identifies intertumoral heterogeneity in prostate cancer phenotypes.

a Heatmap of 141 ROIs averaged to 53 tumor cores grouped by 26 patients highlighting androgen receptor (AR), neuroendocrine (NE), cell cycle progression (CCP), FGFR/MEK, and RB1 gene signatures. Mean gene signature Z-scores are shown according to color scale. Data graphed as boxplots indicating differences in b AR gene signature and c NE gene signature mean Z-score across 138 ROIs averaged from 52 sites from 25 patients with at least two tumors included on the TMAs. Blue circles highlight similar AR signature expression and differential NE signature expression in two different sites of metastasis in patient 15-096. L3 periaortic, M2 diaphragm. Boxes represent the median and interquartile range (IQR) and the upper and lower whiskers extending to the values that are within 1.5 × IQR; data beyond the end of the whiskers are outliers and plotted as points. d Heatmap of 35 ROIs averaged to 14 tumor cores from seven discordant patient samples adapted from a. Results are expressed as mean gene signature Z-scores and presented according to color scale in a. e Volcano plot demonstrating intrapatient heterogeneity in individual 12-011. Two sample sites, bladder and lymph node (LN) with different mCRPC phenotypes were compared (N = 3 regions of interest (ROIs) per site) and genes associated with a NE phenotype are enriched in the sample I8_LN AR+/NE+ when compared to J1_bladder AR+/NE−. f Volcano plot demonstrating intrapatient heterogeneity in individual 15-010. Two sample sites, lymph node and liver with different mCRPC phenotypes were compared (N = 3 ROIs per site). Genes associated with AR signature were enriched in sample K2_LN AR^low/NE− and genes associated with a NE phenotype are enriched in the sample H1_liver AR−/NE+.

Fig. 4. Intratumoral gene expression homogeneity and heterogeneity.

a Heatmap of 141 ROIs from 53 individual tumor cores grouped by 26 patients highlighting androgen receptor (AR), neuroendocrine (NE), cell cycle progression (CCP), FGFR/MEK, and RB1 gene signatures. Results are expressed as gene signature Z-scores and presented according to color scale. b Boxplot demonstrating interpatient, intrapatient, and intratissue correlation across 141 individual DSP ROIs from 53 tumors from 26 patients. c Volcano plot indicating intrapatient heterogeneity in sample 12-005K1. The green arrow highlights genes enriched in AR^low/NE− tumor core relative to the other two cores from the same tissue. d Data graphed as boxplots indicating differences in AR, NE, and CCP gene signature Z-scores across 141 ROIs from 26 patients included on the TMAs. Dotted lines separate each patient. The NE and CCP plots retain the patient ordering by the AR score as shown in the AR score plot. e Heatmap of six ROIs from two individual tumor cores (L3 and M2) from patient 15-096. Results are expressed as gene signature Z-scores and log₂ negative-normalized (NN) gene expression and presented according to the color scales in Fig. 2a. Exact intratumoral homogeneity was 40% based on the associated hypergeometric distribution for possible pairs of samples. Boxes in b and d represent the median and interquartile range (IQR) and the upper and lower whiskers extending to the values that are within 1.5 × IQR; data beyond the end of the whiskers are outliers and plotted as points.

Fig. 5. Intratumoral heterogeneity within full-tumor section 15-096M1.

a Hematoxylin and eosin (H&E) staining of 15-096M1 lymph node metastases. Distinct areas of morphology are demonstrated, cribriform well differentiated prostatic adenocarcinoma (lower left) and undifferentiated high-grade carcinoma (upper right). N = 1 tissue section for H&E staining. b Fluorescent labeling of 15-096M1 lymph node metastases. High PanCK staining is present in the lower left and low PanCK staining is present in the upper right. N = 1 tissue section for fluorescent labeling. c Individual tumor region of interest (ROIs) (200–500 µm) with varying levels of PanCK intensity and differential tumor morphology. N = 1 tissue section for fluorescent labeling and ROI selection. Expression plots of genes known to be associated with AR+/NE− (d), AR−/NE− (e), and AR−/NE+ (f) phenotypes of ROIs 7–12 from 15-096M1. Counts were Q3 normalized and scaled (Z-score) to enable plotting of all genes on the same axes. g–k Comparison of transcript levels of specific genes in ROIs 1–6 (N = 6) from the CK+ tumor region and ROIs 10–12 (N = 3) distant from the CK+ region. Counts were log₂ Q3 normalized. Significance was determined by two-sided Wilcoxon-rank tests (g–k: p = 0.024). Boxes represent the median and interquartile range (IQR) and the upper and lower whiskers extending to the values that are within 1.5 × IQR; data beyond the end of the whiskers are outliers and plotted as points.

Fig. 6. Alternative splicing of AR isoforms is present across metastases, as determined by DSP.

a Heatmap of 141 individual tumor ROIs grouped by 26 patients comparing bulk RNAseq AR-V7 expression to DSP. Results are expressed as gene signature Z-scores and log₂ negative-normalized (NN) gene expression and presented according to color scales. b DSP–RNA probe design of AR full length and AR-V7 variant isoform. c Fluorescent labeling and ROI selection of three cores from two sites of metastasis from patients 14-043 that demonstrated divergent AR-V7 expression determined by RNAseq and DSP. N = 1 TMA section for fluorescent labeling. d Heatmap demonstrating discordance in AR-V7 expression DSP AR-V7 and AR-V7 measured by bulk RNAseq and AR-V7 immunohistochemistry (IHC) in six tumor ROIs from two sites (CC3 and H3) within patient 14-043. Results are presented according to color scales in a.

Fig. 7. DSP describes immune cell microenvironments of distinct phenotypes of mCRPC.

a Heatmap of DSP immune signaling genes across 141 individual regions of interest (ROIs) from 26 patients. Results are expressed as gene signature Z-scores and log₂ mean-centered gene expression and presented according to color scales. b Fluorescently labeled patient core with matched hematoxylin and eosin (H&E) staining representing high and low levels of inflammatory infiltrate. High—17-081P2 is comprised of 70% tumor, 30% stroma, with 100 CD3+ leukocytes present, and 13-012M2 is comprised of 80% tumor, 20% stroma, with 40 CD3+ leukocytes present. Low—15-096M2 is comprised of 90% tumor, 10% stroma with zero CD3+ cells, and 13-104K2 is comprised of 90% tumor, 10% stroma with three CD3+ cells present. Immune cells counted based on CD3 immunohistochemical staining. N = 1 TMA section for fluorescent labeling and H&E staining. c DSP protein depicts overall low levels of intratumoral immune cells. Data are graphed as log₂ signal-to-noise ratio (SNR). d Consistently high expression of B7-H3 is present in the ARpos_NEneg phenotype when compared to other CRPC phenotypes. Expression is consistent across RNA and protein DSP in B7-H3, PD-L1, PD1, with slightly higher expression in protein DSP observed for CTLA4. Data are presented as log₂ normalized counts.