Quantification of spatial pharmacogene expression heterogeneity in breast tumors

Abstract

Background: Chemotherapeutic drug concentrations vary across different regions of tumors and this is thought to be involved in development of chemotherapy resistance. Insufficient drug delivery to some regions of the tumor may be due to spatial differences in expression of genes involved in the disposition, transport, and detoxification of drugs (pharmacogenes). Therefore, in this study, we analyzed the spatial expression of 286 pharmacogenes in six breast cancer tissues using the recently developed Visium spatial transcriptomics platform to (1) determine if these pharmacogenes are expressed heterogeneously across tumor tissue and (2) to determine which pharmacogenes have the most spatial expression heterogeneity.

Methods and results: The spatial transcriptomics technology sequences the transcriptome of 55 um diameter barcoded sections (spots) across a tissue sample. We analyzed spatial gene expression profiles of four biobank-sourced breast tumor samples in addition to two breast tumor sample datasets from 10× Genomics. We define heterogeneity as the interquartile range of read counts. Collectively, we identified 8887 spots in tumor regions, 3814 in stroma, 44 in lymphocytes, and 116 in normal regions based on pathologist annotation of the tissues. We showed statistically significant differences in expression of pharmacogenes in tumor regions compared to surrounding non-tumor regions. We also observed that the most heterogeneously expressed genes within tumor regions were involved in reactive oxygen species (ROS) handling and detoxification mechanisms. GPX4, GSTP1, MGST3, SOD1, CYP4Z1, CYB5R3, GSTK1, and NAT1 showed the most heterogeneous expression within tumor regions.

Conclusions: The heterogeneous expression of these pharmacogenes may have important implications for cancer therapy due to their ability to impact drug distribution and efficacy throughout the tumor. Our results suggest that chemoresistance caused by expression of GPX4, GSTP1, MGST3, and SOD1 may be intrinsic, not acquired, since the heterogeneity is not specific to chemotherapy-treated samples or cell type. Additionally, we identified candidate chemoresistance pharmacogenes that can be further tested through focused follow-up studies.

Keywords: chemotherapy; pharmacogene; resistance; spatial; transcriptomics.

FIGURE 1

Boxplot of all pharmacogenes, across all samples, in tumor regions only, with an interquartile range greater than zero. Each data point on the y‐axis is the number of unique‐UMI reads from a single barcoded spot. The genes on the X‐axis are sorted by interquartile range in descending order

FIGURE 2

Intensity maps showing subsets of pharmacogenes with an interquartile range greater than 0 in any of the groups on the x‐axis. (A) Interquartile range of unique‐UMI reads for each gene for each sample in tumor spots. (B) Interquartile range of UMI‐normalized unique‐UMI reads for each gene for each sample in tumor spots. (C) Interquartile range of unique‐UMI reads for each gene for each annotated region. (D) Interquartile range of UMI‐normalized unique‐UMI reads for each gene for each annotated region

FIGURE 3

Boxplots of UMI‐normalized unique‐UMI reads for subsets of pharmacogenes with the top 31 largest interquartile ranges. Groups are divided into tumor (tumor + DCIS + cellular tumor + desmoplastic tumor), normal, lymphocytes, and stroma. Y‐axis scale is arbitrary, as this figure is intended to show relative differences. DCIS, ductal carcinoma in situ

FIGURE 4

Log2 fold change comparisons for UMI‐normalized unique‐UMI reads, for tumor regions (per sample) versus combined non‐tumor regions across all samples. The X‐axis is in ascending order by average fold change per group. p‐values are based two‐sided Welch's t‐tests, are Bonferroni corrected, and are colored by intensity. Shapes denote which sample the fold change is based on. Values above and below the black bars represent divide‐by‐zero values of infinity or negative infinity; these values were kept because one of the groups in the comparison had an average of zero expression while the other group had non‐zero expression, indicating a potentially significant differential expression. (A) S‐plot showing the distribution for each gene‐sample combination. (B) Global comparison by pharmacogene. (C) Subset of 35 genes with the lowest p‐values. (D) Subset of genes involved in reactive oxygen species handling. (E) ABC transporter subset. (F) SLC transporters subset