Multiplexed screening reveals how cancer-specific alternative polyadenylation shapes tumor growth in vivo

Abstract

Alternative polyadenylation (APA) is strikingly dysregulated in many cancers. Although global APA dysregulation is frequently associated with poor prognosis, the importance of most individual APA events is controversial simply because few have been functionally studied. Here, we address this gap by developing a CRISPR-Cas9-based screen to manipulate endogenous polyadenylation and systematically quantify how APA events contribute to tumor growth in vivo. Our screen reveals individual APA events that control mouse melanoma growth in an immunocompetent host, with concordant associations in clinical human cancer. For example, forced Atg7 3' UTR lengthening in mouse melanoma suppresses ATG7 protein levels, slows tumor growth, and improves host survival; similarly, in clinical human melanoma, a long ATG7 3' UTR is associated with significantly prolonged patient survival. Overall, our study provides an easily adaptable means to functionally dissect APA in physiological systems and directly quantifies the contributions of recurrent APA events to tumorigenic phenotypes.

Fig. 1. Global 3′ UTR length correlates with clinical outcomes across tumor types.

A TCGA cancer subtypes analyzed and schematic of quantification of gene-level 3′ UTR lengths for all genes per sample. We followed standard practice to compute the median 3′ UTR length per sample, normalized such that a value closer to 0 indicates, on average, globally shorter 3′ UTRs and a value closer to 1 indicates globally longer 3′ UTRs. Values are based on the percent of distal poly(A) site usage index computed by the DaPars algorithm. B Distribution of median 3′ UTR per sample for the 424 tumor samples in The Cancer Genome Atlas cutaneous melanoma cohort. Samples were binned into terciles corresponding to short (blue), medium (gray), and long (yellow) global 3′ UTR samples. The median of all samples is marked with a dashed black line. C–E Individual read coverage plots (read coverage reflective of TPM values) for three exemplar genes that exhibit significant differences in 3′ UTR length between samples with globally shorter (blue) and longer (yellow) UTRs as indicated in (B) for EGLN1, EIF3A and SRSF11 in cutaneous melanoma (SKCM). Poly(A) signal sequences (PAS) sites indicated as diamonds. F–H Violin plots of gene level 3′ UTR length for EGLN1, EIF3A and SRSF11 from (C–E) comparing short and long median 3′ UTR stratified cutaneous melanoma samples. A value closer to 0 indicates higher use of the proximal poly(A) site and a value closer to 1 indicates higher use of the distal poly(A) site (p values from two-sided Wilcoxon rank-sum test). Exact P values are P = 1.178188⁻³³, 2.047001⁻⁴⁸ and 2.354266⁻⁴⁰ for F, G and H, respectively. I Kaplan–Meier analysis comparing overall survival of TCGA cutaneous melanoma samples binned as short or long median 3′ UTR samples (p values from a two-sided logrank test). J Volcano plot of the log₂(hazard ratio) calculated from univariate cox regression models comparing overall survival of short versus long UTR stratified samples plotted against the −log₁₀(p value). P value from a Cox proportional hazard model. Cancer subtypes with a p value < 0.05 are indicated in red. Vertical dashed line indicates a log₂(hazard ratio) of 0, and the horizontal dashed line indicated a P value of 0.05.

Fig. 2. Identification of differentially polyadenylated RNAs in a murine model of melanoma.

A Description of two syngeneic model cell lines, B16-F10 cells (red) and Melan-A cells (gray) both derived from a C57BL/6 background melanocyte origin. B Differences in information provided by RNA-seq and Poly(A)-seq methods, and examples of read coverage plots reflective of the reads generated by each sequencing approach across the terminal exon including the 3′ UTR of an example gene. C Read coverage plots of RNA-seq and Poly(A)-seq completed for Melan-A cells (gray) B16-F10 cells (red) of the Sap30l terminal exon and 3′ UTR with the stop codon and the annotated poly(A) signal sequences (PAS). D BAM coverage plot of RNA-seq and Poly(A)-seq completed for Melan-A cells (gray) and B16-F10 cells (red) of the Higd1a terminal exon and 3′ UTR with the stop codon and the annotated poly(A) signal sequences (PAS). E Box plot demonstrating relative 3′ UTR length calculated as the log₂(distal reads / proximal reads) per sample for the gene Higd1a and Sap30l. Data reflects six Poly(A)-seq runs per cell line. P values from two-sided Wilcoxon rank-sum test. F Volcano plot of all differentially polyadenylated transcripts between B16-F10 cells and Melan-A cells quantified using the APAlyzer pipeline. Data reflects six Poly(A)-seq runs per cell line, significantly altered events were determined using a two-sided Student’s t test. Significantly shortened 3′ UTRs are indicated in blue and significantly lengthened 3′ UTRs are indicated in yellow. Vertical dashed line indicates a difference of 0, and the horizontal dashed line indicated a P value of 0.05. G Scatter plot of all genes identified as differentially polyadenylated from Poly(A)-seq data, comparing gene-level 3′ UTR length differences and gene expression differences between B16-F10 and Melan-A cells. Blue indicates genes that are significantly shortened in B16-F10 cells and display a significant difference in expression levels, yellow indicates genes that are significantly lengthened and display a significant difference in expression levels, and gray indicates the gene shows no significant change in expression between the two cell lines. The horizontal dashed lines indicate a log₂(fold-change) in gene expression of 1 and −1.

Fig. 3. A CRISPR-Cas9 paired-guide RNA strategy can force targeted distal poly(A) site usage.

A Schematic of CRISPR-Cas9 paired-guide RNA (pgRNA) approach to excise proximal poly(A) sites and force use of distal poly(A) sites. B Diagram of Sap30l terminal exon with two distinct 3′ UTRs which vary depending on use of a proximal or distal poly(A) signal. Schematic of two distinct sets of proximal poly(A) KO (pKO) pgRNAs designed to delete the proximal poly(A) signal in Sap30l, annotated as Sap30l pKO1 and pKO2. C Genotyping PCR of polyclonal B16-F10 Cas9-expressing cells treated with a non-targeting control (NTC) or one of two distinct Sap30l pKO pgRNAs, with the calculated percent of signal showing DNA excision (representative gel from n = 3 biological replicates). D Read coverage plots of Poly(A)-seq completed for Cas9-expressing B16-F10 cells treated with either a control pgRNA (gray), Sap30l pKO1 (green), or Sap30l pKO1 (blue). Schematic illustrates the Sap30l terminal exon and 3′ UTR with the annotated poly(A) signal sequences (PAS). E Box plots of the percent distal isoform (pdi) usage per cell line quantified from the nested RT-PCR of the Sap30l 3′ UTR for each monoclonal cell line grouped by genotype (WT, heterozygous or homozygous deletion of the proximal poly(A) signal). Data reflects six, seven, and five distinct monoclonal lines with wild-type, heterozygous, or homozygous deletion of the proximal poly(A) signal sequence, respectively. P value from a two-sided Wilcoxon rank-sum test. F Box plots of Sap30l mRNA abundance per cell line measured by q-RT-PCR for each monoclonal cell line grouped by genotype (WT, heterozygous or homozygous deletion of the proximal poly(A) signal). Data reflects six, seven, and five distinct monoclonal lines with wild-type, heterozygous, or homozygous deletion of the proximal poly(A) signal sequence, respectively. P value from a two-sided Wilcoxon rank-sum test. G Scatter plot and Pearson correlation of Sap30l percent distal isoform usage versus Sap30l mRNA abundance. Gray error band represents 95% confidence interval. P value from a Pearson correlation.

Fig. 4. A high-throughput functional CRISPR-Cas9 screen reveals APA events that influence melanoma growth.

A Schematic of CRISPR-Cas9 pgRNA library design (8–10 pgRNAs per target), cloning and paired in vitro and in vivo screen format, and timing. B Mean log₂ (fold-change) of all pgRNAs targeting a specific gene or proximal poly(A) site normalized to negative control pgRNAs targeting proximal poly(A) sites of unexpressed genes for paired in vitro and in vivo screens in Cas9-expressing B16-F10 cells (n = 8 biological replicates). Black outlines, significantly enriched (green) or depleted (purple) proximal poly(A) targets. No outlines, positive (light green) and negative (light purple) control genes. Dashed lines indicate a log₂(fold-change) of 0.5 and −0.5, respectively for in vitro (vertical) and in vivo (horizontal) data. C Log fold-changes associated with the indicated targets. Each point, reflects the mean value of a single pgRNA across n = 8 biological replicates. Dashed line indicates a normalized fold-change of 1. D As (C), but in vivo. Dashed line indicates a normalized fold-change of 1. E Read coverage plots of Poly(A)-seq data for Cas9-expressing B16-F10 cells treated with either a control pgRNA (gray) or Atg7 pKO1 (orange). Schematic shows the Atg7 terminal exon and 3′ UTR with the annotated poly(A) signal sequences (PAS). F Western blot of lysates from Cas9-expressing B16-F10 cells treated with a control pgRNA or Atg7 pKO pgRNA. ATG7 protein level normalized to alpha-tubulin control. G In vitro cell growth of Cas9-expressing B16-F10 cells treated with a control pgRNA or Atg7 pKO pgRNA measured by CellTiter-Glo. Measurement is average of three experimental replicates +/- S.E.M. Significance denoted as *p < 0.05, **p < 0.01 or ***p < 0.001 using a two-sided Student’s t test (exact p values 0.023, 0.00011, 9 × 10⁻⁸, and 0.00096). H As (E), but for Egln1 and Egln1 pKO1. I Western blot of lysates from Cas9-expressing B16-F10 cells treated with a control pgRNA or Egln1 pKO pgRNA. Egln1 protein level normalized to alpha-tubulin control. J In vitro cell growth of Cas9-expressing B16-F10 cells treated with a control pgRNA or Egln1 pKO pgRNA measured by CellTiter-Glo. Measurement is the average of three experimental replicates +/- S.E.M. Significance denoted as *p < 0.05, **p < 0.01 or ***p < 0.001 using a two-sided Student’s t test (exact p values 0.0005, 0.0011, 0.0079, and 0.059).

Fig. 5. Atg7 alternative poly(A) site selection alters melanoma cell growth in vitro and in vivo.

A Nested RT-PCR of the Atg7 3′ UTR in Cas9-expressing B16-F10 cell lines treated with 1 µg/mL Actinomycin D to inhibit transcription. Expected Atg7 3′ UTR size when utilizing the proximal poly(A) site (short) or the distal poly(A) site (long) are indicated. Representative gel from n = 2 biological replicates. B Quantification of the percent distal isoform from (A). Gel intensities quantified using FIJI. C Survival data from a cohort of C57BL/6 mice injected with Cas9-expressing B16-F10 cells treated with a control pgRNA targeting the poly(A) site of an unexpressed gene (control) or an Atg7 pKO pgRNA (n = 6 mice, 12 tumors, per condition). D Left, representative immunohistochemistry images of control or Atg7 pKO pgRNA tumor sections stained for Ki67. Right, images with nuclei classified using HALO as negative (blue), weak (yellow), moderate (orange) or strong (red) staining. E Stacked bar plot quantifying data from (D). n = 4 images from 4 distinct tumors per genotype. For each image, the entire slide is processed, only excluding areas if they are easily discernible as non-tumor tissue. P value calculated with a two-sided binomial proportion test. F Representative flow cytometry histograms of Cas9-expressing B16-F10 cells treated with a control or Atg7 pKO pgRNA cells stained with propidium iodide and then analyzed using Dean-Jett-Fox classification for cell cycle stage from FlowJo v10. G Rank order plot of loss-of-function (LOF) mutations normalized per kilobase per patient detected in genomic DNA sequencing from TCGA cutaneous melanoma cohort (n = 424 patients). Vertical dashed line, point at which all remaining ranked genes have no detected LOF mutations. H Kaplan–Meier analysis of progression-free survival in TCGA cutaneous melanoma cohort for patients binned into terciles based on ATG7 3′ UTR length. P value from a two-sided logrank test.