Abstract
Drug resistance and relapse remain key challenges in pancreatic cancer. Here, we have used RNA-seq, ChIP-seq, and genome-wide CRISPR analysis to map molecular dependencies of pancreatic cancer stem cells, highly therapy-resistant cells that preferentially drive tumorigenesis and progression. This integrated genomic approach revealed an unexpected utilization of immuno-regulatory signals by pancreatic cancer epithelial cells. In particular, the nuclear hormone receptor RORγ, known to drive inflammation and T-cell differentiation, was upregulated during pancreatic cancer progression, and its genetic or pharmacologic inhibition led to a striking defect in pancreatic cancer growth, and a marked improvement in survival. Further, a large-scale retrospective analysis in patients revealed that RORγ expression may predict pancreatic cancer aggressiveness, as it positively correlated with advanced disease and metastasis. Collectively, these data identify an orthogonal co-option of immuno-regulatory signals by pancreatic cancer stem cells, suggesting that autoimmune drugs should be evaluated as novel treatment strategies for pancreatic cancer patients.
Graphical Abstract
In Brief
Pancreatic cancer stem cells co-opt immunoregulatory pathways, a vulnerability that could be exploited therapeutically by agents currently in trials for autoimmune diseases
Introduction
While cytotoxic agents remain the standard of care for most cancers, their use is often associated with initial efficacy, followed by disease progression. This is particularly true for pancreatic cancer, a highly aggressive disease, where current multidrug chemotherapy regimens result in tumor regression in 30% of patients, quickly followed by disease progression in the vast majority of cases (Conroy et al., 2011). This progression is largely due to the inability of chemotherapy to successfully eradicate all tumor cells, leaving behind subpopulations that can trigger tumor re-growth. Thus, identifying the cells that are preferentially drug resistant, and understanding their vulnerabilities, is critical to improving patient outcome and response to current therapies.
In previous work, several groups have focused on identifying the most tumorigenic populations within pancreatic cancer. Through this, subpopulations of cells marked by expression of CD24+/CD44+/ESA+ (Li et al., 2007), cMet (Li et al., 2011), CD133 (Hermann et al., 2007), Nestin (Kawamoto et al., 2009), ALDH (Rasheed et al., 2010), and more recently DCLK1(Bailey et al., 2014) and Musashi (Fox et al., 2016), have been shown to harbor stem cell characteristics, in being enriched for the capacity to drive tumorigenesis and recreate the heterogeneity of the original tumor (Reya et al., 2001). Importantly, these tumor propagating cells or cancer stem cells have been shown to be highly resistant to cytotoxic therapies, such as gemcitabine, consistent with the finding that cancer patients with a high cancer stem cell signature have poorer prognosis relative to those with a low stem cell signature (Grosse-Wilde et al., 2015). Although pancreatic cancer stem cells are epithelial in origin, these cells frequently express EMT-associated programs, which may in part explain their over-representation in circulation and propensity to seed metastatic sites (Fox et al., 2016; Hermann et al., 2007). Because these studies define stem cells as a population that present a particularly high risk for disease progression defining the molecular signals that sustain them remains an essential goal for achieving complete and durable responses.
Here we have used a combination of RNA-seq, ChIP-seq and genome-wide CRISPR screening to define the molecular framework that sustains the aggressive nature of pancreatic cancer stem cells. These studies identified a network of key nodes regulating pancreatic cancer stem cells, and revealed an unanticipated role for immuno-regulatory genes in their self-renewal and maintenance. Among these, the retinoic acid receptor-related orphan receptor gamma (RORγ) a nuclear hormone receptor known for its role in Th17 cell specification and regulation of inflammatory cytokine production (Ivanov et al., 2006), emerged as a key regulator of stem cells. RORγ expression increased with progression, and its blockade via genetic or pharmacologic approaches depleted the cancer stem cell pool and profoundly inhibited human and mouse tumor propagation, in part by suppressing a super-enhancer (SE) associated oncogenic network. Finally, sustained treatment with a RORγ inhibitor led to a significant improvement in autochthonous models of pancreatic cancer. Together, our studies offer a unique comprehensive map of pancreatic cancer stem cells and identify critical vulnerabilities that may be exploited to improve therapeutic targeting of aggressive, drug resistant pancreatic cells.
Results
Transcriptomic and epigenetic map of pancreatic cancer cells reveals a unique stem cell state
In previous work, we used the KPf/fC mouse model (Hingorani et al., 2003) of pancreatic ductal adenocarcinoma (PDAC) to show that a reporter mouse designed to mirror expression of the stem cell signal Musashi (Msi) could identify tumor cells that are preferentially drug resistant and can drive tumor re-growth (Fox et al., 2016). Consistent with this, Msi2+ tumor cells were 209-fold enriched in the ability to give rise to organoids in limiting dilution assays (Figure 1A, S1A–B) (Boj et al., 2015). Because Msi+ cells were enriched for tumor propagation and drug resistance -classically defined properties of cancer stem cells- we postulated that Msi reporters could be used as a tool to understand the molecular underpinnings of this aggressive subpopulation within pancreatic cancer.
Figure 1. Transcriptomic and epigenetic map of pancreatic cancer cells reveals a unique stem cell state.
(A) Tumor organoid formation from primary Msi2+ and Msi2- REM2-KPf/fC tumor cells. Representative images, scale=100um.
(B) RNA-seq and ChIP-seq of EpCAM+GFP+ and Epcam+GFP- REM2-KPf/fC tumor cells (n=3 RNA-seq, n=1 ChIP-seq).
(C) Principal component analysis of KPf/fC stem (purple) and non-stem (gray) cells.
(D) Transcripts enriched in stem (red, pink) and non-stem cells (dark blue, light blue). Pink, light blue, lfdr<0.3; red, dark blue, lfdr<0.1.
(E-J) GSEA cell states and corresponding heat-maps associated with development (E and F), metabolism (G and H), and cancer relapse (I and J). (E, G, and I) Red denotes overlapping gene signatures; blue denotes nonoverlapping gene signatures. (F, H, and J) Red, over-represented gene expression; blue, under-represented gene expression; shades denote fold change.
(K) Single-cell sequencing of KPR172H/+C tumors (left) and map of Msi2 expression in ETC and EMT clusters (right); CAF, cancer-associated fibroblasts (red); EMT, mesenchymal tumor cells (yellow/green); Endo, endothelial cells (green); ETC, epithelial tumor cells (blue); TAM, tumor-associated macrophages (magenta).
(L) Hockey stick plots of H3K27ac occupancy ranked by signal density. Stem cell SEs (left), or shared SEs (right) are demarcated by highest ranking and intensity signals.
(M-O) H3K27ac ChIP-seq reads across genes marked by stem cell SEs (M), shared SEs (N), or non-stem SEs (O).
To map the functional genomic landscape of the stem cell state, we utilized a combination of RNA-seq, ChIPseq and genome-wide CRISPR screening (Sanjana et al., 2014). Pancreatic cancer cells were isolated from Msi2-reporter (REM2) KPf/fC mice based on GFP and EpCAM expression and analyzed by RNA-seq (Figure 1B). Principal component analysis showed that KPf/fC reporter+ tumor cells were distinct from reporter- tumor cells at a global transcriptional level, and were defined by a unique set of programs in turn driven by the differential expression of over a thousand genes (Figure 1C–D). We focused on genes enriched in stem cells in order to understand the transcriptional programs that may functionally maintain the stem cell state. Gene Set Enrichment Analysis (GSEA) (Subramanian et al., 2005) was used to compare this PDAC stem cell transcriptomic signature with other cell signatures (Table S1). This revealed that the transcriptional state of PDAC stem cells mapped closely with other developmental and stem cell states, indicating molecular features aligned with their observed functional traits (Figure 1E–F). Additionally, the transcriptional signature of PDAC stem cells was inversely correlated with cell proliferation signatures (Figure S1C–D), consistent with our finding that stem cells become quiescent following chemotherapy (Figure S1E). Stem cells also harbored metabolic signatures associated with tumor aggressiveness including increased sulfur amino acid metabolism (Ryu et al., 2011), and enhanced glutathione synthesis, pathways that enable survival following radiation and chemotherapy (Lu et al., 2017) (Figure 1G–H). Finally, the stem cell transcriptome bore similarities to signatures from relapsed cancers of the breast, liver, and colon programs (Figure 1I–J); aligned with this, stem cells showed a significant overlap with mesenchymal cells in single cell RNA-seq analysis of pancreatic tumors (Figure 1K). These molecular properties may collectively underlie the ability of PDAC stem cells to survive chemotherapy and drive tumor recurrence.
Analysis of H3 lysine-27 acetylation (H3K27ac, Figures 1B, S1F), a histone mark associated with active enhancers (Hnisz et al., 2013), revealed that the differential gene expression programs in stem cells and nonstem cells were driven by changes at the chromatin level. Thus, genomic regions enriched for H3K27ac coincided with regions where gene expression was increased in each cell type (Figure S1G–J; stem cells: R2=0.28, p=7.1×10−14, non-stem cells R2=0.46, p=22×10−16). Because SEs have been proposed to be key drivers of cell identity (Hnisz et al., 2013; Whyte et al., 2013), we mapped shared and unique SEs in stem and non-stem cells (Figure 1L–O). This analysis revealed that SE associated H3K27ac marks were predominantly restricted to either stem cells or non-stem cells, with 65% of all SEs being unique to each population (364 unique SEs in stem cells/388 unique SEs in non-stem cells). In contrast, almost all promoter and conventional enhancerassociated H3K27ac marks were shared between stem and non-stem cells, with less than 5% being unique. Further, while SEs in the stem cell population were clearly demarcated by peaks with substantially greater relative enrichment than the same regions in non-stem cells (Figure 1M), the SEs found in non-stem cells showed a peak intensity that was only marginally greater than the corresponding regions in stem cells (Figure 1O). These data suggest that stem cells in pancreatic cancer have a more specialized SE landscape than nonstem cells, and raise the possibility that SE linked genes and their regulators may serve to control stem cell identity in pancreatic cancer. In support of this, key transcription factors and programs that underlie developmental and stem cell states, such as Tead4, Wnt7b and Msi2 (Figure 1L) and Foxp, Klf7, and Hmga1 (Table S2) were associated with SEs in KPf/fC stem cells.
Genome-wide CRISPR screen identifies core functional programs in pancreatic cancer
To define which of the programs uncovered by the transcriptional and epigenetic analyses represented true functional dependencies of stem cells, we carried out a genome-wide CRISPR screen. Thus, primary cell cultures enriched for stem cells (Figure S2A) were derived from REM-KPf/fC mice and transduced with the mouse GeCKO CRISPRv2 sgRNA library (Sanjana et al., 2014) (Figure 2A). The screen was multiplexed in order to identify genes required in conventional 2 dimensional cultures, as well as in 3 dimensional stem cell sphere cultures (Rovira et al., 2010) that selectively allow stem cell growth (Fox et al., 2016) (Figure 2A). The screens showed clear evidence of selection, with 807 genes depleted in 2D (Figure 2B, C) and an additional 178 in 3D stem cell cultures (Figure 2B, D). Importantly, the screens showed a loss of oncogenes and an enrichment of tumor suppressors in conventional cultures (Figures 2C, S2B), and a loss of stem cell signals and gain of negative regulators of stem signals in stem cell conditions (Figures 2D, S2C).
Figure 2. Genome-scale CRISPR screen identifies core stem cell programs in pancreatic cancer.
(A) Schematic of CRISPR screen.
(B) Number of guides in each replicate following lentiviral infection (gray bars), puromycin selection (red bars), and sphere formation (blue bars).
(C and D) Volcano plots of guides depleted in 2D (C) and 3D (D). Genes indicated on plots, p<0.005.
(E) Network propagation integrating transcriptomic, epigenetic and functional analysis of stem cells. Stemenriched genes by RNA-seq (log2FC>2) and depleted in 3D (FDR-adjusted p<0.5) were used to seed the network (triangles), then analyzed for protein-protein interactions. Each node represents a single gene; color denotes RNA-seq fold change; stem enriched, red; non-stem enriched, blue; not differentially expressed, gray. Labels shown for genes enriched in stem cells by RNA-seq and ChIP-seq (Up/Up) or enriched in nonstem cells by RNA-seq and ChIP-seq (Down/Down); RNA log2FC absolute value>2.0, ChIP-seq FDR<0.01. Seven core programs were defined by gene groups with high connectivity; annotated by GO analysis (FDR<0.05).
Computational integration of the transcriptomic and CRISPR-based functional genomic data was carried out using a network propagation method similar to one developed previously (Vanunu et al., 2010). First, the network was seeded with genes that were preferentially enriched in stem cells and also identified as essential for stem cell growth (Figure 2E). The genes most proximal to the seeds were then determined using the mouse STRING interactome (Szklarczyk et al., 2015) based on known and predicted protein-protein interactions using network propagation. Fold-change in RNA expression from the RNA-seq was overlaid onto the resulting subnetwork. The network was subsequently clustered into functional communities based on high interconnectivity between genes, and gene set over-representation analysis was performed on each community; this analysis identified seven subnetworks built around distinct biological pathways, thus providing a sytems level view of core programs that may be involved in driving pancreatic cancer growth. These programs identified stem and pluripotency pathways, developmental and proteasome signals, lipid metabolism/nuclear receptors, cell adhesion/cell-matrix/cell migration, and immuno-regulatory signaling as pathways integral to the stem cell state (Figures 2E, S2D).
Hijacked immunoregulatory programs as direct regulators of pancreatic cancer cells
Ultimately the power of such a map is the ability to identify and understand key new functional dependencies. Thus, we used the network map as a framework to select an integrated gene set based on the transcriptomic, epigenomic and the CRISPR analysis (Table S3). Selected genes were subsequently targeted via viral shRNA delivery into KPf/fC cells, and the impact on pancreatic cancer propagation assessed by sphere assays in vitro or tracking tumor growth in vivo. While many genes within the pluripotency/development core program were known to be important in pancreatic cancer (e.g. Wnt, Hedgehog and Hippo pathways), others such as Onecut3 and Tudor3, genes previously implicated in motor neuron development or in stress response, presented new opportunities for discovery and emerged as signals essential for pancreatic cancer stem cell growth (Figures 3A, S3A, Table S4). Further, novel metabolic factors such as Sptssb, a key contributor to sphingolipid metabolism (Zhao et al., 2015), and Lpin2, an enzyme involved in generation of pro-inflammatory very-low density lipoproteins (Dwyer et al., 2012), were found to be critical stem cell dependencies, implicating lipid metabolism as a key point of control in pancreatic cancer (Figure 3B, Table S4). This analysis also identified new gene families in pancreatic cancer: thus within the adhesion/cell-matrix core program (Figures 3C–J, S3B), several members of the multiple EGF repeat (MEGF) subfamily of orphan adhesion G protein coupled receptors (8 of 12) were preferentially expressed in stem cells (Figure 3E). Among this set, inhibition of Celsr1, Celsr2 (Figure S3C–D), and Pear1/Jedi (Figure S3E) triggered apoptosis, depleted Msi+ stem cells and potently blocked cancer propagation in vitro and in vivo (Figures 3G–J, S3F–J, Table S4). These pathways will likely be important to explore further, especially because GPCRs can frequently serve as effective drug targets.
Figure 3. Identification of novel pathway dependencies of pancreatic cancer stem cells.
(A-D) Genes from developmental processes (A), lipid metabolism (B), and cell adhesion, motility, and matrix components (C and D) were inhibited via shRNA in KPf/fC cells and sphere or flank tumor growth assessed. Sphere, n=3–6; flank transplant, n=4.
(E-I) Relative RNA expression of MEGF family and related (*Celsr1) genes in KPf/fC stem and non-stem cells (E). Red, over-represented; blue, under-represented; color denotes fold change from median values. Impact of inhibiting Celsr1, Celsr2, and Pear1 on KPf/fC sphere formation (F) and flank transplants (G-I). Sphere, n=3–6; flank transplant, n=4.
(J) Impact of shRNA-mediated inhibition of Pear1 in human FG cells on colony formation (n=3) and flank tumor propagation assessed (n=4).
(K-L) Single-cell sequencing of KPR172H/+C tumors (K) and tumor cells expressing IL10Rβ, IL34, and Csf1R (L). CAF, cancer-associated fibroblasts (red); EMT, mesenchymal tumor cells (yellow/green); Endo, endothelial cells (green); ETC, epithelial tumor cells (blue); TAM, tumor-associated macrophages (magenta)
(M) Impact of shRNA-mediated inhibition of IL-10rβ and Csf1R on sphere formation of KPf/fC cells, n=3–6.
(N) Impact of shRNA-mediated inhibition of IL-10rβ and Csf1R on stem content (Msi2-GFP+) of KPf/fC cells; in 3D culture assessed, n=3.
(O-P) Impact of shRNA-mediated inhibition of IL-10rβ (Ο) and Csf1R (P) on KPf/fC flank transplant growth, n=4.
(Q) Impact of shRNA-mediated inhibition of IL10Rβ in human FG cells on sphere formation, n=3, or flank transplant, n=4.
(R) Impact of shRNA-mediated inhibition of IL-10 and IL-34 on KPf/fC sphere formation, n=3.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
An unexpected discovery from this map was the identification of immune pathways/cytokine signaling as a core program. In line with this, retrospective analysis of the RNA-seq and ChIP-seq analysis revealed that multiple immunoregulatory cytokine receptors and their ligands were expressed in stem/non-stem tumor epithelial cells (Figure S3K). This was of particular interest because many genes associated with this program, such as IL-10, IL-34 and CSF1R have been previously studied in context of the tumor microenvironment, but have not been reported to be produced by, or to functionally impact, pancreatic epithelial cells directly. Single-cell RNA-seq analysis of KPR172H/+C tumor cells (Figure 1K, 3K, S3L) confirmed the presence of IL10Rβ, IL34 and CSF1R in epithelial tumor cells (Figure 3L), as well as in Msi2+ cancer stem cells (Figure S3M). Consistent with expression in stem cells, inhibition of IL10Rβ and CSF1R led to a marked loss of sphere forming capacity and reduced stem cells (Figure 3M, N, S3N–O) in vitro and impaired tumor growth and propagation in vivo (Figure 3O–Q, Figure S3P–Q). The activity of IL10Rβ and CSF1R may, at least in part, be ligand-dependent as their ligands were both expressed in epithelial cells (Figure S3R), and the impact of ligand and receptor inhibition mirrored each other (Figure 3R). Collectively, these findings demonstrate an orthogonal co-option of inflammatory mediators by pancreatic cancer stem cells and suggest that agents that modulate cytokine networks may directly impact pancreatic cancer propagation.
RORγ, a mediator of T cell fate, is a critical dependency in pancreatic cancer
To understand how the gene networks defined above are controlled, we focused on transcription factors because of their broad role in initiating programs key to cell fate and identity (Neph et al. 2012). Of the 53 transcription factors identified within the map, 12 were found to be enriched in stem cells by transcriptomic and epigenetic parameters (Figure S4A), and included several pro-tumorigenic pioneer factors such as Sox9 (Kopp et al., 2012) and Foxa2 (Bailey et al., 2016). Among transcription factors with no known role in pancreatic cancer (Arntl2, Nr1d1 and RORγ), only RORγ was actionable in the near term, with clinical-grade antagonists currently available (Table S5) (Gege, 2016). Motif enrichment analysis identified RORγ sites as preferentially enriched in chromatin regions uniquely open in stem cells (Figure S4B) and in open chromatin regions that corresponded with enriched gene expression in stem cells (Figure S4B). These findings were consistent with RORγ having a preferential role in controlling gene expression programs important for defining the stem cell state in pancreatic cancer.
RORγ was an unanticipated dependency as it is a nuclear hormone receptor that has been predominantly studied in Th17 cell differentiation (Ivanov et al., 2006) as well as in metabolism in context of the circadian rhythm (Cook et al., 2015); consistent with this, it mapped to both the hijacked cytokine signaling/immune subnetwork and the nuclear receptor/metabolism subnetwork (Figures 2E, S2D). While RORγ expression was low in normal murine pancreas (data not shown), it rose dramatically in KPf/fC tumors. Within epithelial tumor cells, RORγ expression was highly enriched in stem cells relative to non-stem cells (Figures 4A, S4C–D), mapping to individual EpCAM+Msi+ cells in single cell RNA-seq analysis (Figure S4E). RORγ was also expressed in KPR172H/+C tumor cells (not shown) suggesting it is active across models of pancreatic cancer. Importantly, RORγ expression in mouse models was predictive of expression in human pancreatic cancer: thus, while RORγ expression was low in the normal human pancreas and in pancreatitis, its expression increased significantly in epithelial tumor cells with disease progression (Figures 4B–C, S4F). Interestingly RORγ levels decreased with inhibition of IL1R signaling suggesting that the upstream regulators of RORγ in pancreatic cancer and in Th17 cells may be shared (Figure S4G). Functionally, shRNA-mediated knockdown (Figure S4H) confirmed the role of RORγ identified by the genetic CRISPR-based screen as it decreased stem cell sphere formation in both KPR172H/+C and KPf/fC cells (Figure 4D–E). At a cellular level, RORγ inhibition led to increased cell death (Figure S4I), decreased proliferation (Figure S4I) and an ultimate depletion of Msi+ stem cells (Figure 4F). Importantly, tumor cells lacking RORγ showed a striking defect in tumor initiation and propagation in vivo, with a 11-fold reduction in final tumor volume (Figures 4G, S4J). Finally, analysis of KPf/fC mice crossed to either RORγ null (Ivanov et al., 2006) or wild type controls revealed that targeted genetic deletion of RORγ can trigger an overall decrease in tumor burden; this ranged from reduced tumor weight or cellularity to the presence of more normal and benign PanIN lesions and reduced areas of adenocarcinoma in the pancreata (Figure 4H–I).
Figure 4. The immuno-regulatory gene RORγ is a critical dependency of pancreatic cancer.
(A) Rorc expression in stem and non-stem REM2-KPf/fC tumor cells; representative of three biological replicates.
(B) Representative images of RORγ expression in normal adjacent human pancreas (left), PanINs (middle), and PDAC (right). RORγ (green), E-Cadherin (red), Dapi (blue), scale=50um.
(C) Frequency of RORγ+ cells within E-Cadherin+ epithelial fraction in patient samples quantified by immunofluorescence; Normal adjacent, n=3; pancreatitis, n=8; PanIN 1, n=10; PanIN 2, n=6; PDAC, n=8.
(D-E) Impact of shRNA-mediated RORγ inhibition on 3D growth of KPR172H/+C (D) and KPf/fC (E) cells, n=3 per shRNA.
(F) Impact of shRNA-mediated RORγ inhibition on Msi2-GFP stem cell content of in KPf/fC cells in 3D culture (H) n=3.
(G) Impact of shRNA-mediated RORγ inhibition on flank tumor growth of KPf/fC cells, n=4.
(H-I) Reduced tumor burden in Rorc−/−KP f/f C mice. Age-matched WT KP f/f C and Rorc−/−KP f/f C set shown displayed reduced tumor cell number (H), and reduced adenocarcinoma content (I); low-grade PanIN indicated with red arrow, PDAC indicated with black arrow, scale=100um; n=3 mice from 8–10 wks of age; representative plots and images from matched mice are displayed.
(J and K) Relative RNA expression of stem cell programs (J) and pro-tumorigenic factors (K) in KPf/fC cells transduced with shCtrl or shRorc. Red, over-represented; blue, under-represented; color denotes fold change. (L) Venn diagram of genes downregulated with RORγ loss (q<0.05, purple), Stem specific SE-associated genes (green), and genes associated with H3K27ac peaks with RORγ consensus binding sites (orange).
(M) Number of RORγ, CBFB, Sox9 binding sites found in stem cell SEs relative to random genomic background of equivalent base-pair coverage (p<0.05).
(N) Relative RNA expression of SE-associated oncogenes in KPf/fC cells transduced with shCtrl or shRorc. Red, over-represented; blue, under-represented; color denotes fold change from median values.
(O) H3K27ac ChIP-seq reads for genes marked by stem cell SEs and downregulated in RORγ-depleted KPf/fC cells. Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
To define the transcriptional programs RORγ controls in pancreatic cancer cells, we used a combination of ChIP-seq and RNA-seq and found that RORγ knockdown led to extensive changes in transcriptional programs key to driving cancer growth: this included stem cell signals such as Wnt, BMP, and Fox (Figure 4J), and pro-tumorigenic signals such as Hmga2 (Figure 4K). Further, 28% of stem cell SE linked genes were downregulated in cells lacking RORγ (Figure 4L). Consistent with this, ChIP-seq analysis of active chromatin regions identified RORγ binding sites as disproportionately present in stem cell SEs compared to other transcription factors such as CBFB, or even the pioneer factor Sox9 (Figure 4M). Additional SE-linked stem cell genes regulated by RORγ included Msi2, Klf7 and Ehf (Figures 4N–O), potent oncogenic signals that can control cell fate. Mechanistically, loss of RORγ did not markedly impact the stem cell SE landscape in two independent KPf/fC derived lines (Figure S4K–M), suggesting that it may instead bind a preexisting landscape to preferentially impact transcriptional changes. These data collectively suggest that RORγ is an upstream regulator of a powerful SE-linked oncogenic network in pancreatic cancer stem cells.
The finding that RORγ is a key dependency in pancreatic cancer was particularly exciting, as multiple inhibitors have been developed to target this pathway in autoimmune disease (Huh and Littman, 2012). Pharmacologic blockade of RORγ using the inverse agonist SR2211 (Kumar et al., 2012) decreased sphere and organoid formation in both KPf/fC and KPR172H/+C cells (Figures 5A–D). To assess the impact of the inhibitor in vivo, SR2211 was delivered, either alone or in combination with gemcitabine, into immunocompetent KPf/fC derived tumor bearing mice (Figures 5E, S5A). SR2211 significantly reduced tumor growth as a single agent (Figures 5F–G); further, while gemcitabine alone had no impact on the stem cell burden, SR2211 triggered a 3-fold depletion in CD133+ and Msi+ cells, and an 11-fold depletion of CD133+ and 6 fold depletion of Msi2+ cells in combination with gemcitabine (Figure 5H, I). This suggests the exciting possibility that SR2211 can eradicate chemotherapy resistant cells (Figure 5H, I). Finally, to assess any impact on survival, we delivered the RORγ inhibitor into autochthonous, tumor-bearing KPf/fC mice; while none of the vehicle-treated mice were alive 25 days after the initiation of treatment, 75% of mice that received SR2211 were still alive at this point and 50% were alive even at 45 days after treatment initiation. SR2211 not only doubled median survival −18 days for vehicle-treated mice and 38.5 days for SR2211-treated mice- but also led to a 6-fold reduction in the risk of death (Figure 5J, HR=0.16). Hmga2, identified originally from the RNA-seq as a downstream target of RORγ, was downregulated in pancreatic epithelial cells following SR2211 delivery in vivo, suggesting effective target engagement at mid-point during treatment, although this was less apparent in end stage tumors and may explain why treated mice ultimately succumbed to disease (Figure S5B–C). Collectively, these data show that pancreatic cancer stem cells are profoundly dependent on RORγ and suggest that its inhibition may lead to a significant improvement in disease control. Further, the fact that its impact on tumor burden was amplified several fold when combined with gemcitabine suggests that it may synergize with chemotherapy to more effectively target tumors that remain refractory to therapy.
Figure 5. Pharmacologic targeting of RORγ impairs progression and improves survival in mouse models of pancreatic cancer.
(A and B) 3D growth of KPf/fC cells (A) and KPR172H/+C cells (B) in the presence of the SR2211 or vehicle (n=3).
(C and D) KPf/fC organoid formation in the presence of SR2211 or vehicle. Representative images (C) and quantification (D), scale=100um.
(E-I) Analysis of flank KPf/fC tumor-bearing mice treated with SR2211 or vehicle for 3 weeks. Strategy (E). Total live cells (F), total EpCAM+ tumor cells (G), total EpCAM+/CD133+ stem cells (H), and total EpCAM+/Msi2+ stem cells (I) (n=4 vehicle, n=2 vehicle+gemcitabine, n=4 SR2211, n=3 SR2211+gemcitabine).
(J) Survival of KPf/fC mice treated daily with vehicle (gray) or SR2211 (black) (p=0.051, Hazard ratio=.16, Median survival: vehicle=18 days, SR2211=38.5 days).
(K) Live imaging of REM2-KPf/fC mice treated with vehicle or SR2211 for 8 days (n=2). Msi2-reporter (green), VE-Cadherin (magenta), Hoecsht (blue); Msi2-reporter+ stem cells, gray box.
(L) Quantification of stem cell clusters from REM2-KPf/fC live imaging (n=2; 6–10 frames analyzed per mouse).
To visualize whether RORγ blockade impacts tumor progression by targeting stem cells, SR2211 was delivered in REM2-KPf/fC mice with late-stage autochthonous tumors and responses tracked via live imaging. In vehicletreated mice, large stem cell clusters could be readily identified throughout the tumor based on GFP expression driven by the Msi reporter (Figures 5K–L). SR2211 led to a marked depletion of the majority of large stem cell clusters within 1 week of treatment (Figures 5K–L), with no increased necrosis observed in surrounding tissues. This unique spatio-temporal analysis suggests that stem cell depletion is an early consequence of RORγ blockade and highlights the REM2-KPf/fC model as an effective platform to assess the impact of new agents on therapy resistant cells.
Since treatment with the inhibitor in immunocompetent mice or in patients in vivo could have an impact on both cancer cells and immune cells we tested the effect of SR2211 in the context of an immunocompromised environment. SR2211 significantly impacted growth of KPf/fC tumors in an immunodeficient background (Figure 6A–B) suggesting that inflammatory T cells were not necessary for its effect. Further, in chimeric mice where wild type tumors were transplanted into either wild type or RORγ null recipients, tumors grew equivalently (Figure 6C–D), suggesting that loss of RORγ in only the immune cells (such as Th17) and microenvironment, has no detectable impact on tumor growth. Finally, we delivered SR2211 into these chimeric mice to test if RORγ antagonists may influence tumor growth via Th17 cells, and found that the impact of SR2211 on tumor growth, cellularity, and stem cell content was equivalent in chimeric wild type and RORγ recipient mice (Figure 6E–L). These data collectively suggest that most of the observed effect of RORγ inhibition is tumor cell specific and not indirect through an environmental/Th17 dependence on RORγ (Figure 6E–L, S6A–B). Consistent with a primarily epithelial cell impact, we did not detect any significant impact of SR2211 on nonneoplastic cells such as CD45+, CD31+, MDSC, macrophage, dendritic, or T cells within the tumors at early time points (Figure S6C–M). These data do not preclude the possibility that RORγ inhibitors may act on both tumor cells and immune cells in human tumors if more inflammatory T cells were present.
Figure 6. RORγ is a direct dependency of pancreatic tumor epithelial cells.
(A-B) Analysis of flank KPf/fC tumor-bearing NSG mice treated with SR2211 or vehicle for 2 weeks. Strategy (A). Flank tumor growth following treatment with vehicle or SR2211 for 2 weeks (B). Fold change in tumor volume relative to volume at enrollment (n=4–6).
(C-D) Strategy (C). Growth of KPf/fC flank tumors in WT or RORγ−/− recipient mice (D) (n=3–4).
(E-L). Strategy (E). Flank tumor growth in WT recipients treated with vehicle or SR2211 for 2 weeks (F). Flank tumor growth in RORγ−/−- recipients treated with vehicle or SR2211 for 2 weeks (G). Tumor mass (H), total live cells (I), total EpCAM+ tumor cells (J), total EpCAM+/CD133+ stem cells (K), and total Th17 cells (L) in WT and RORγ−/− recipients (n=5–7).
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
To further explore the functional relevance of RORγ to human pancreatic cancer, RORγ was inhibited through both genetic and through pharmacologic means in human PDAC cells. CRISPR based disruption of RORγ led to a ~3 to 9-fold loss of colony formation in human FG PDAC cells (Figure 7A). To test if RORγ inhibition could block human tumor growth in vivo, we transplanted human PDAC cells into the flanks of immunocompromised mice and allowed tumors to become palpable before beginning treatment (Figure 7B). Compared to vehicle-treatment, SR2211 delivery was highly effective and tumor growth was essentially extinguished with a nearly 6-fold reduction in growth in mice receiving SR2211 (Figure 7C). Primary patient tumor cells were also remarkably sensitive to RORγ blockade, with a ~300 fold reduction in total organoid volume following SR2211 treatment (Figures 7D–F) and a severe reduction of in vivo tumor growth in primary patient-derived xenografts (Figure 7G). Mechanistically, RNA-seq and Gene Ontology (GO) analysis of human FG and KPf/fC cells identified a set of cytokines/growth factors as key common RORγ driven programs: thus Semaphorin 3c, its receptor Neuropilin2, Oncostatin M, and Angiopoietin, all highly pro-tumorigenic factors harboring RORγ binding motifs were shared targets of RORγ in mouse and human pancreatic cancer (Figure S7A–D). The dependence of human pancreatic tumors on RORγ function are exciting in light of the fact that genomic amplification of RORC occurs in ~12% of pancreatic cancer patients (Figure 7H). This raises the possibility that RORC status could serve as a biomarker for patients who may be particularly responsive to RORC inhibition.
Figure 7. RORγ is required for human pancreatic cancer growth and predicts advanced disease.
(A) Human FG colony formation after RORC CRISPR knockdown, n=5.
(B) Representative images of RORγ expression in human FG tumors, RORγ (green), E-Cadherin (red), Dapi (blue).
(C) Human FG tumor growth in mice treated with gemcitabine and either vehicle or SR2211 for 2.5 weeks. Tumor volume fold change is relative to volume at enrollment.
(D-F) Primary patient organoid growth in the presence of vehicle or SR2211. Representative images of organoids in Matrigel (D) (scale=100um), following recovery from Matrigel (E) (scale=50um), and quantification of organoid circumference (F, left) or volume (F, right).
(G) Growth of primary patient-derived xenografts treated with vehicle or SR2211 for 1.5 weeks. (n=4)
(H) RORC amplification in tumors of patients diagnosed with various malignancies.
(I-L) Representative TMAs of PDAC and PanINs illustrating scoring for negative, cytoplasmic, and cytoplasmic + nuclear RORγ staining (I). Correlation between RORγ staining and tumor stage (J), lymphatic invasion (K), and lymph node status (L).
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
Lastly, to determine whether expression of RORγ could serve as a prognostic for specific clinicopathologic features, we performed RORγ immunohistochemistry on tissue microarrays from a clinically annotated retrospective cohort of 116 PDAC patients (Table S6). For 69 patients, matched pancreatic intraepithelial neoplasia (PanIN) lesions were available. RORγ protein was detectable (cytoplasmic expression only/low or cytoplasmic and nuclear expression/high, Figure 7I) in 113 PDAC cases and 55 PanIN cases, respectively, and absent in 3 PDAC cases and 14 PanIN cases, respectively. Compared to cytoplasmic expression, nuclear RORγ expression in PDAC cases was significantly correlated with higher pathological tumor (pT) stages at diagnosis (Figures 7J). In addition, RORγ expression in PanIN lesions was positively correlated with lymphatic vessel invasion (L1, Figure 7K) and lymph node metastasis (pN1, pN2, Figure 7L) by the invasive carcinoma. These results indicate that RORγ expression in PanIN lesions and nuclear RORγ localization in invasive carcinoma could be useful markers to predict PDAC aggressiveness.
Discussion
It is an unfortunate truth that the most common outcome for pancreatic cancer patients following a response to cytotoxic therapy is not cure, but eventual disease progression and death driven by drug resistant stem cellenriched populations (Fox et al., 2016; Van den Broeck et al., 2013). The work we report here has allowed us to develop a comprehensive molecular map of the core dependencies of pancreatic cancer stem cells by integrating their epigenetic, transcriptomic and functional genomic landscape. This dataset thus provides a novel resource for understanding therapeutic resistance and relapse, and for discovering new vulnerabilities in pancreatic cancer. As an example, the MEGF family of orphan receptors represent a potentially actionable family of adhesion GPCRs, as this class of signaling receptors have been considered druggable in cancer and other diseases (Lappano and Maggiolini, 2011). Importantly, our epigenetic analyses revealed a significant relationship between SE-associated genes and functional dependencies in stem cell conditions; stem cell-unique SE-associated genes were more likely to drop out in the CRISPR screen in stem cell conditions compared to SE-associated genes in non-stem cells (Figure S7D). This provides additional evidence for the epigenetic and transcriptomic link to functional dependencies in cancer stem cells, and further supports previous findings that SE-linked genes may be more important for maintaining cell identity and more sensitive to perturbation (Whyte et al., 2013).
From the screens presented here, we identified an unexpected dependence of KPf/fC stem cells on inflammatory and immune mediators, such as the CSF1R/IL-34 axis and IL-10R signaling. While these have been previously thought to act primarily on immune cells in the microenvironment (Guillonneau et al., 2017; Wang et al., 2017), our data suggest that stem cells may have evolved to co-opt this cytokine-rich milieu, allowing them to resist effective immune-based elimination. These findings also suggest that agents targeting CSF1R, which are under investigation for pancreatic cancer (Sankhala KK, 2017), may act not only on the tumor microenvironment but also directly on pancreatic epithelial cells themselves. Our studies also raise the possibility that therapies designed to activate the immune system to attack tumors may have effects on tumor cells directly: just as we have learned chemotherapy can kill tumor cells but may also impair the immune system, therapies designed to activate the immune system such as IL-10 may also promote the growth of tumor cells. This dichotomy of action will need to be considered in order to better optimize immunomodulatory treatment strategies.
A major new discovery driven by the network map was the identification of RORγ as a key immuno-regulatory pathway hijacked in pancreatic cancer. This together with prior work implicating RORγ in prostate cancer models (Wang et al., 2016) suggests that this pathway may not be restricted to pancreatic cancer but may be more broadly utilized in other epithelial cancers. Interestingly, while cytokines such as IL17, IL21, IL22, and CSF2 are known targets of RORγ in Th17 cells, none of these were downregulated in Rorc-deficient pancreatic tumor cells. The fact that RORγ regulated potent oncogenes marked by SEs in stem cells, suggests it may be critical for defining the stem cell state in pancreatic cancer. The basis of this intriguing epithelial specific activity of RORγ will be an important area for future exploration. In addition, the network of genes impacted by RORγ inhibition included other immune-modulators such as CD47, raising the possibility that it may also mediate interaction with the surrounding niche and immune system cells. Finally, one particularly exciting aspect of this work is the identification of RORγ as a potential therapeutic target in pancreatic cancer. Given that inhibitors of RORγ are currently in Phase II trials for autoimmune diseases (Gege, 2016), our findings suggest that repositioning these agents as pancreatic cancer therapies warrants further investigation.
STAR METHODS
Contact for Reagent and Resource Sharing
Further information and requests for resources and reagents should be directed to and will be fulfilled by the Lead Contact, Tannishtha Reya (treya@ucsd.edu).
Experimental Model and Subject Details
Mice
REM2 (Msi2eGFP/+) reporter mice were generated as previously described (Fox et al., 2016); all of the reporter mice used in experiments were heterozygous for the Msi2 allele. The LSL-KrasG12D mouse, B6.129S4-Krastm4Tyj/J (Stock No: 008179), the p53flox/flox mouse, B6.129P2- Trp53tm1Brn/J (Stock No: 008462), and the RORγ -knockout mouse (Stock No: 007571), were purchased from The Jackson Laboratory. Dr. Chris Wright provided Ptf1a-Cre mice as previously described (Kawaguchi et al., 2002). LSL-R172H mutant p53, Trp53R172H mice were provided by Dr. Tyler Jacks as previously described (Olive et al., 2004) (JAX Stock No: 008183). The mice listed above are immunocompetent, with the exception of RORγ -knockout mice which are known to lack TH17 T-cells as described previously (Ivanov et al., 2006); these mice were maintained on antibiotic water (sulfamethoxazole and trimethoprim) when enrolled in flank transplantation and drug studies as outlined below. Immune compromised NOD/SCID (NOD.CB17-Prkdcscid/J, Stock No: 001303) and NSG (NOD.Cg-PrkdcscidIL2rgtm1Wji/SzJ, Stock No: 005557) mice purchased from The Jackson Laboratory. All mice were specific-pathogen free, and bred and maintained in the animal care facilities at the University of California San Diego. Animals had access to food and water ad libitum, and were housed in ventilated cages under controlled temperature and humidity with a 12 hour light-dark cycle. All animal experiments were performed according to protocols approved by the University of California San Diego Institutional Animal Care and Use Committee. No sexual dimorphism was noted in all mouse models. Therefore, males and females of each strain were equally used for experimental purposes and both sexes are represented in all data sets. All mice enrolled in experimental studies were treatment-naïve and not previously enrolled in any other experimental study.
Both REM2-KPf/fC and WT-KPf/fC mice (REM2; LSL-KrasG12D/+; Trp53f/f; Ptf1a-Cre and LSLKrasG12D/+, ; Trp53f/f; Ptf1a-Cre respectively) were used for isolation of tumor cells, establishment of primary mouse tumor cell and organoid lines, and autochthonous drug studies as described below. REM2-KPf/fC and KPf/fC mice were enrolled in drug studies between 8 to 11 weeks of age, and were used for tumor cell sorting and establishment of cell lines when they reached end-stage disease between 10 and 12 weeks of age. REM2KPf/fC mice were used for in vivo imaging studies between 9.5–10.5 weeks of age. KPR172HC (LSL-KrasG12D/+, ; Trp53R172h/+; Ptf1a-Cre) mice were used for cell sorting and establishment of tumor cell lines when they reached end-stage disease between 16–20 weeks of age. In some studies, KPf/fC-derived tumor cells were transplanted into the flanks of immunocompetent littermates between 5–8 weeks of age. Littermate recipients (WT or REM2LSL-KrasG12D/+, ; Trp53f/f or Trp53f/f mice) do not develop disease or express Cre. NOD/SCID and NSG mice were enrolled in flank transplantation studies between 5 to 8 weeks of age; KPf/fC derived cell lines and human FG cells were transplanted subcutaneously for tumor propagation studies in NOD/SCID recipients and patientderived xenografts and KPf/fC derived cell lines were transplanted subcutaneously in NSG recipients as described in detail below.
Human and mouse pancreatic cancer cell lines
Mouse primary pancreatic cancer cell lines and organoids were established from end-stage, treatment-naïve KPR172HC and WT- and REM2-KPf/fC mice as follows: tumors from endpoint mice (10–12 weeks of age for KPf/fC or 16–20 weeks of age for KPR172HC mice) were isolated and dissociated into single cell suspension as described below. Cells were then either plated in 3D sphere or organoid culture conditions detailed below, or plated in 2D in 1x DMEM containing 10% FBS, 1x pen/strep, and 1x non-essential amino acids. At the first passage in 2D, cells were collected and resuspended in HBSS (Gibco, Life Technologies) containing 2.5% FBS and 2 mM EDTA, then stained with FC block followed by 0.2 µg/106 cells anti-EpCAM APC (eBioscience). EpCAM+ tumor cells were sorted then re-plated for at least one additional passage. To evaluate any cellular contamination and validate the epithelial nature of these lines, cells were analyzed by flow cytometry again at the second passage for markers of blood cells (CD45-PeCy7, eBioscience), endothelial cells (CD31-PE, eBioscience), and fibroblasts (PDGFR-PacBlue, Biolegend). Cell lines were derived from both female and male KPR172HC and WT- and REM2-KPf/fC mice equivalently; both sexes are equally represented in the cell-based studies outlined below. Functional studies were performed using cell lines between passage 2 and passage 6. Human FG cells were originally derived from a PDAC metastasis and have been previously validated and described (Morgan et al., 1980). Patient-derived xenograft cells and organoids were derived from originallyconsented (now deceased) PDAC patients and use was approved by UCSD’s IRB; cells were de-identified and therefore no further information on patient status, treatment or otherwise, is available. FG cell lines were cultured in 2D conditions in 1x DMEM (Gibco, Life Technologies) containing 10% FBS, 1x pen/strep (Gibco, Life Technologies), and 1x non-essential amino acids (Gibco, Life Technologies). 3D in vitro culture conditions for all cells and organoids are detailed below.
Patient cohort for PDAC tissue microarray
The PDAC patient cohort and corresponding TMAs used for RORγ immunohistochemical staining and analysis have been reported previously (Wartenberg et al., 2018). Patient characteristics are detailed in Table S6. Briefly, a total of 4 TMAs with 0.6 mm core size was constructed: three TMAs for PDACs, with samples from the tumor center and invasive front (mean number of spots per patient: 10.5, range: 2–27) and one TMA for matching PanINs (mean number of spots per patient: 3.7, range: 1–6). Tumor samples from 116 patients (53 females and 63 males; mean age: 64.1 years, range: 34–84 years) with a diagnosis of PDAC were included. Matched PanIN samples were available for 69 patients. 99 of these patients received some form of chemotherapy; 14 received radiotherapy. No sexual dimorphism was observed in any of the parameters assessed, including overall survival (p=0.227), disease-free interval (p=0.3489) or RORγ expression in PDAC (p=0.9284) or PanINs (p=0.3579). The creation and use of the TMAs were reviewed and approved by the Ethics Committee at the University of Athens, Greece, and the University of Bern, Switzerland, and included written informed consent from the patients or their living relatives.
Method Details
In vitro and in vivo experimental strategies
Tissue dissociation, cell isolation, and FACS analysis
Mouse pancreatic tumors were washed in MEM (Gibco, Life Technologies) and cut into 1–2 mm pieces immediately following resection. Tumor pieces were collected into a 50 ml Falcon tube containing 10 ml Gey’s balanced salt solution (Sigma), 5 mg Collagenase P (Roche), 2 mg Pronase (Roche), and 0.2 µg DNAse I (Roche). Samples were incubated for 20 minutes at 37°C, then pipetted up and down 10 times and returned to 37°C. After 15 more minutes, samples were pipetted up and down 5 times, then passaged through a 100 µm nylon mesh (Corning). Red blood cells were lysed using RBC Lysis Buffer (eBioscience) and the remaining tumor cells were washed, then resuspended in HBSS (Gibco, Life Technologies) containing 2.5% FBS and 2 mM EDTA for staining, FACS analysis, and cell sorting. Analysis and cell sorting were carried out on a FACSAria III machine (Becton Dickinson), and data were analyzed with FlowJo software (Tree Star). For analysis of cell surface markers by flow cytometry, 5×105 cells were resuspended in HBSS containing 2.5% FBS and 2 mM EDTA, then stained with FC block followed by 0.5 µl of each antibody. For intracellular staining, cells were fixed and permeabilized using the BrdU flow cytometry kit (BD Biosciences); Annexin V apoptosis kit was used for analysis of apoptotic cells (eBioscience). The following rat antibodies were used: anti-mouse EpCAM-APC (eBioscience), anti-mouse CD133-PE (eBioscience), anti-mouse CD45-PE and PE/Cy7 (eBioscience), anti-mouse CD31-PE (BD Bioscience), anti-mouse Gr-1-FITC (eBioscience), antimouse F4/80-PE (Invitrogen), anti-mouse CD11b-APC (Affymetrix), anti-mouse CD11c-BV421 (Biolegend), anti-mouse CD4-FITC (eBioscience) and CD4-Pacific blue (Bioglegend), anti-mouse CD8-PE (eBioscience), anti-mouse IL-17-APC (Biolegend), anti-mouse BrdU-APC (BD Biosciences), and anti-mouse Annexin-V-APC (eBioscience). Propidium-iodide (Life Technologies) was used to stain for dead cells.
In vitro growth assays
We describe below the distinct growth assays used for pancreatic cancer cells. Colony formation is an assay in Matrigel (thus adherent/semi-adherent conditions), while tumorsphere formation is an assay in non-adherent conditions. We have found that cell types from different sources grow better in different conditions. For example, the murine KPR172H/+C and the human FG cell lines grow much better in Matrigel, while KPf/fC cell lines often grow well in non-adherent, sphere conditions (though they can also grow in Matrigel).
Pancreatic tumorsphere formation assay
Pancreatic tumorsphere formation assays were performed and modified from (Rovira et al., 2010). Briefly, lowpassage (<6 passages) WT or REM2-KPf/fC cell lines were infected with lentiviral particles containing shRNAs; positively infected (red) cells were sorted 72 hours after transduction. 100–300 infected cells were suspended in tumorsphere media: 100 µl DMEM F-12 (Gibco, Life Technologies) containing 1x B-27 supplement (Gibco, Life Technologies), 3% FBS, 100 µM Β-mercaptoethanol (Gibco, Life Technologies), 1x non-essential amino acids (Gibco, Life Technologies), 1x N2 supplement (Gibco, Life Technologies), 20 ng/ml EGF (Gibco, Life Technologies), 20 ng/ml bFGF2 (Gibco, Life Technologies), and 10 ng/ml ESGRO mLIF (Thermo Fisher). Cells in media were plated in 96-well ultra-low adhesion culture plates (Costar) and incubated at 37°C for 7 days. KPf/fC in vitro tumorsphere formation studies were conducted at a minimum of n=3 independent wells per cell line across two independent shRNA of n=3 wells; however, the majority of these experiments were additionally completed in >1 independently-derived cell lines n=3, at n=3 wells per shRNA. shRNA sequences and average knockdown efficiencies are available in Table S7.
Matrigel colony assay
For FG and KPR172H/+C cells, 300–500 cells were resuspended in 50 µl tumorsphere media as described below, then mixed with Matrigel (BD Biosciences, 354230) at a 1:1 ratio and plated in 96-well ultra-low adhesion culture plates (Costar). After incubation at 37°C for 5 min, 50 µl tumorsphere media was placed over the Matrigel layer. Colonies were counted 7 days later. For RORγ inhibitor studies, SR2211 or vehicle was added to cells in tumorsphere media, then mixed 1:1 with Matrigel and plated. SR2211 or vehicle was also added to the media that was placed over the solidified Matrigel layer. For FG colony formation, n=5 independent wells across 5 independent CRISPR sgRNA and two independent non-targeting gRNA. KPR172H/+C cells were plated at n=3 wells per shRNA from one cell line.
Organoid culture assays
Tumors from 10–12 week old end stage REM2-KPf/fC mice were harvested and dissociated into a single cell suspension as described above. Tumor cells were stained with FC block then 0.2 µg/106 cells anti-EpCAM APC (eBioscience). Msi2+/EpCAM+ (stem) and Msi2-/EpCAM+ (non-stem) cells were sorted, resuspended in 20 µl Matrigel (BD Biosciences, 354230). For limiting dilution assay, single cells were resuspended in matrigel at the indicated numbers from 20,000 to 10 cells/20uL and were plated as a dome in a pre-warmed 48 well plate. After incubation at 37°C for 5 min, domes were covered with 300 µl PancreaCult Organoid Growth Media (Stemcell Technologies). Organoids were imaged and quantified 6 days later. Limiting dilution analysis for stemness assessment was performed using web based- extreme limiting dilution analysis (ELDA) software (Hu and Smyth, 2009). Msi2+/EpCAM+ (stem) and Msi2-/EpCAM+ (non-stem) organoids were derived from n=3 independent mice and plated at the indicated cell numbers.
Organoids from REM2-KPf/fC were passaged at ~1:2 as previously described (Boj et al., 2015). Briefly, organoids were isolated using Cell Recovery Solution (Corning 354253), then dissociated using Accumax AM105), and plated in 20 µl matrigel (BD Biosciences, 354230) domes on a pre-warmed 48-well plate. After incubation at 37°C for 5 min, domes were covered with 300 µl PancreaCult Organoid Growth Media (Stemcell Technologies). SR2211 (Cayman Chemicals 11972) was resuspended in DMSO at 20 mg/ml, diluted 1:10 in 0.9% NaCl containing 0.2% acetic acid, and further diluted in PancreaCult Organoid Media (Stemcell Technologies) to the indicated dilutions. Organoids were grown in the presence of vehicle or SR2211 for 4 days, then imaged and quantified, n=3 independent wells plated per dose per treatment group.
Primary patient organoids were established and provided by Dr. Andrew Lowy. Briefly, patient-derived xenografts were digested for 1 hour at 37°C in RPMI containing 2.5% FBS, 5mg/ml Collagenase II, and 1.25mg/ml Dispase II, then passaged through a 70 µM mesh filter. Cells were plated at a density of 1.5 × 105 cells per 50 µl Matrigel. After domes were solidified, growth medium was added as follows: RPMI containing 50% Wnt3a conditioned media, 10% R-Spondin1-conditioned media, 2.5% FBS, 50 ng/ml EGF, 5 mg/ml Insulin, 12.5 ng/ml hydrocortisone, and 14 µM Rho Kinase Inhibitor. After establishment, organoids were passaged and maintained as previously described (Boj et al., 2015). Briefly, organoids were isolated using Cell Recovery Solution (Corning 354253), then dissociated into single cell suspensions with TrypLE Express (ThermoFisher 12604) supplemented with 25 µg/ml DNase I (Roche) and 14 µM Rho Kinase Inhibitor (Y27632, Sigma). Cells were split 1:2 into 20 µl domes plated on pre-warmed 48 well plates. Domes were incubated at 37°C for 5 min, then covered with human complete organoid feeding media (Boj et al., 2015) without Wnt3a-conditioned media. SR2211 was prepared as described above, added at the indicated doses, and refreshed every 3 days. Organoids were grown in the presence of vehicle or SR2211 for 7 days, then imaged and quantified, n=3 independent wells plated per dose per treatment group. All images were acquired on a Zeiss Axiovert 40 CFL. Organoids were counted and measured using ImageJ 1.51s software.
Flank tumor transplantation studies
For the flank transplantation studies outlined below, investigators blinded themselves when possible to the assigned treatment group of each tumor for analysis; mice were de-identified after completion of flow cytometry analysis. The number of tumors transplanted for each study is based on past experience with studies of this nature, where a group size of 10 is sufficient to determine if pancreatic cancer growth is significantly affected when a regulatory signal is perturbed (see Fox et al., 2016).
For shRNA-infected pancreatic tumor cell propagation in vivo, cells were infected with lentiviral particles containing shRNAs and positively infected (red) cells were sorted 72 hours after transduction. 1000 low passage, shRNA-infected KPf/fC, or 2×105 shRNA-infected FG cells were resuspended in 50 µl culture media, then mixed 1:1 with matrigel (BD Biosciences). Cells were injected subcutaneously into the left or right flank of 5–8 week-old NOD/SCID recipient mice. Subcutaneous tumor dimensions were measured with calipers 1–2x weekly for 6–8 weeks, and two independent transplant experiments were conducted for each shRNA at n=4 independent tumors per group.
For drug-treated KPf/fC flank tumors, 2×104 low passage REM2-KPf/fC tumor cells were resuspended in 50 µl culture media, then mixed 1:1 with matrigel (BD Biosciences). Cells were injected subcutaneously into the left or right flank of 5–8 week-old non-tumor bearing, immunocompetent littermates or NSG mice. Tumor growth was monitored twice weekly; when tumors reached 0.1–0.3 cm3, mice were randomly enrolled in treatment groups and were treated for 3 weeks as described below. After 3 weeks of therapy, tumors were removed, weighed, dissociated, and analyzed by flow cytometry. Tumor volume was calculated using the standard modified ellipsoid formula ½ (Length x Width2); n=2–4 tumors per treatment group in immunocompetent littermate recipients and n=4–6 tumors per treatment group in NSG recipients.
For chimeric transplantation studies, 2×104 low passage REM2-KPf/fC tumor cells were resuspended in 50 µl culture media, then mixed 1:1 with matrigel (BD Biosciences). Cells were injected subcutaneously into the left or right flank of 5–8 week-old RORγ-knockout or wild-type recipients; recipient mice were maintained on antibiotic water (sulfamethoxazole and trimethoprim). Tumor growth was monitored twice weekly; when tumors reached 0.1–0.3 cm3, mice were randomly enrolled in treatment groups and were treated for 3 weeks as described below. After 3 weeks of therapy, tumors were removed, weighed, dissociated, and analyzed by flow cytometry. Tumor volume was calculated using the standard modified ellipsoid formula ½ (Length x Width2); n=5–7 tumors per treatment group.
For drug-treated human pancreatic tumors 2×104 human pancreatic FG cancer cells or 2×106 patientderived xenograft cells were resuspended in 50 µl culture media, then mixed 1:1 with matrigel (BD Biosciences). Cells were injected subcutaneously into the left or right flank of 5–8 week-old NSG recipient mice. Mice were randomly enrolled in treatment groups and were treated for 3 weeks as described below. After 3 weeks of therapy, tumors were removed, weighed, and dissociated. Subcutaneous tumor dimensions were measured with calipers 1–2x weekly. Tumor volume was calculated using the standard modified ellipsoid formula ½ (Length x Width2); at minimum n=4 tumors per treatment group.
In vivo and in vitro drug therapy
The RORγ inverse agonists SR2211 (Cayman Chemicals, 11972, or Tocris, 4869) was resuspended in DMSO at 20 mg/ml or 50 mg/ml, respectively, then mixed 1:20 in 8% Tween80-PBS prior to use. Gemcitabine (Sigma, G6423) was resuspended in H2O at 20 mg/ml. For in vitro drug studies, low passage (<6 passage) WT- or REM2-KPf/fC cells, (<10 passage) KPR172H/+C cells, or FG cells were plated in non-adherent tumorsphere conditions or Matrigel colony conditions for 1 week in the presence of SR2211 or vehicle. For KPf/fC littermate, NSG mice, and RORγ-knockout mice bearing KPf/fC-derived flank tumors and for NSG mice bearing flank patient-derived xenograft tumors, mice were treated with either vehicle (PBS) or gemcitabine (25 mg/kg i.p., 1x weekly) alone or in combination with vehicle (5% DMSO, 8% Tween80-PBS) or SR2211 (10 mg/kg i.p., daily) for 3 weeks. RORγ-knockout mice and paired wild-type littermates were maintained on antibiotic water (sulfamethoxazole and trimethoprim). For NOD/SCID mice bearing flank FG tumors, mice were treated with either vehicle (5% DMSO in corn oil) or SR2211 (10 mg/kg i.p., daily) for 2.5 weeks. All flank tumors were measured 2x weekly and mice were sacrificed if tumors were >2cm3, in accordance with IACUC protocol. For KPf/fC autochthonous survival studies, 8 week old tumor-bearing KPf/fC mice were enrolled in either vehicle (10% DMSO, 0.9% NaCl with 0.2% acetic acid) or SR2211 (20 mg/kg i.p., daily) treatment groups, and treated until moribund, where n=4 separate mice per treatment group. For all drug studies, tumor-bearing mice were randomly assigned into drug treatment groups; treatment group size was determined based on previous studies (Fox et al., 2016).
Immunofluorescence staining
Pancreatic cancer tissue from KPf/fC mice was fixed in Z-fix (Anatech Ltd, Fisher Scientific) and paraffin embedded at the UCSD Histology and Immunohistochemistry Core at The Sanford Consortium for Regenerative Medicine according to standard protocols. 5 µm sections were obtained and deparaffinized in xylene. The human pancreas paraffin embedded tissue array was acquired from US Biomax, Inc (BIC14011a). For paraffin embedded mouse and human pancreas tissues, antigen retrieval was performed for 40 minutes in 95–100°C 1x Citrate Buffer, pH 6.0 (eBioscience). Sections were blocked in PBS containing 0.1% Triton X100 (Sigma- Aldrich), 10% Goat Serum (Fisher Scientific), and 5% bovine serum albumin (Invitrogen).
KPf/fC cells and human pancreatic cancer cell lines were suspended in DMEM (Gibco, Life Technologies) supplemented with 50% FBS and adhered to slides by centrifugation at 500 rpm. 24 hours later, cells were fixed with Z-fix (Anatech Ltd, Fisher Scientific), washed in PBS, and blocked with PBS containing 0.1% Triton X-100 (Sigma-Aldrich), 10% Goat serum (Fisher Scientific), and 5% bovine serum albumin (Invitrogen). All incubations with primary antibodies were carried out overnight at 4°C. Incubation with Alexafluor-conjugated secondary antibodies (Molecular Probes) was performed for 1 hour at room temperature. DAPI (Molecular Probes) was used to detect DNA and images were obtained with a Confocal Leica TCS SP5 II (Leica Microsystems). The following primary antibodies were used: chicken anti-GFP (Abcam, ab13970) 1:500, rabbit anti-RORγ (Thermo Fisher, PA5–23148) 1:500, mouse anti-E-Cadherin (BD Biosciences, 610181) 1:500, anti-Keratin (Abcam, ab8068) 1:15, anti-Hmga2 (Abcam. Ab52039) 1:100, anti-Celsr1 (EMD Millipore abt119) 1:1000, anti-Celsr2 (BosterBio A06880) 1:250.
Tumor imaging
9.5–10.5 week old REM2-KPf/fC mice were treated either vehicle or SR2211 (10 mg/kg i.p., daily) for 8 days. For imaging, mice were anesthetized by intraperitoneal injection of ketamine and xylazine (100/20 mg/kg). In order to visualize blood vessels and nuclei, mice were injected retro-orbitally with AlexaFluor 647 anti-mouse CD144 (VE-cadherin) antibody and Hoechst 33342 immediately following anesthesia induction. After 25 minutes, pancreatic tumors were removed and placed in HBSS containing 5% FBS and 2mM EDTA. 80–150 µm images in 1024 × 1024 format were acquired with an HCX APO L20x objective on an upright Leica SP5 confocal system using Leica LAS AF 1.8.2 software. GFP cluster sizes were measure using ImageJ 1.51s software. 2 mice per treatment group were analyzed in this study; 6–10 frames were analyzed per mouse.
Analysis of tissue microarrays
Immunohistochemistry (IHC) and staining analysis
TMAs were sectioned to 2.5 µm thickness. IHC staining was performed on a Leica BOND RX automated immunostainer using BOND primary antibody diluent and BOND Polymer Refine DAB Detection kit according to the manufacturer’s instructions (Leica Biosystems). Pre-treatment was performed using citrate buffer at 100°C for 30 min, and tissue was stained using rabbit anti-human RORγ(t) (polyclonal, #PA5–23148, Thermo Fisher Scientific) at a dilution of 1:4000. Stained slides were scanned using a Pannoramic P250 digital slide scanner (3DHistech). RORγ(t) staining of individual TMA spots was analyzed in an independent and randomized manner by two board-certified surgical pathologists (C.M.S and M.W.) using Scorenado, a custommade online digital TMA analysis tool. Interpretation of staining results was in accordance with the “reporting recommendations for tumor marker prognostic studies” (REMARK) guidelines. Equivocal and discordant cases were re-analyzed jointly to reach a consensus. RORγ(t) staining in tumor cells was classified microscopically as 0 (absence of any cytoplasmic or nuclear staining), 1+ (cytoplasmic staining only), and 2+ (cytoplasmic and nuclear staining). For patients in whom multiple different scores were reported, only the highest score was used for further analysis. Spots/patients with no interpretable tissue (less than 10 intact, unequivocally identifiable tumor cells) or other artifacts were excluded.
Statistical analysis of TMA data
Descriptive statistics were performed for patients’ characteristics. Frequencies, means, and range values are given. Association of RORγ(t) expression with categorical variables was performed using the Chi-square or Fisher’s Exact test, where appropriate, while correlation with continuous values was tested using the nonparametric Kruskal-Wallis or Wilcoxon test. Univariate survival time differences were analyzed using the Kaplan-Meier method and log-rank test. All p-values were two-sided and considered significant if <0.05.
shRNA lentiviral constructs and production
Short hairpin RNA (shRNA) constructs were designed and cloned into pLV-hU6-mPGK-red vector by Biosettia. The target sequences are listed in Supplementary Table S7. Virus was produced in 293T cells transfected with 4 µg shRNA constructs along with 2 µg pRSV/REV, 2 µg pMDLg/pRRE, and 2 µg pHCMVG constructs (Dull et al., 1998; Sena-Esteves et al., 2004). Viral supernatants were collected for two days then concentrated by ultracentrifugation at 20,000 rpm for 2 hours at 4°C. Knockdown efficiency for the shRNA constructs used in this study varied from 45–95% (Table S7).
RT-qPCR analysis
RNA was isolated using RNeasy Micro and Mini kits (Qiagen) and converted to cDNA using Superscript III (Invitrogen). Quantitative real-time PCR was performed using an iCycler (BioRad) by mixing cDNAs, iQ SYBR Green Supermix (BioRad) and gene specific primers. Primer sequences are available in Table S7. All real time data was normalized to B2M or Gapdh.
Genome-wide profiling and bioinformatic analysis
Primary Msi2+ and Msi2- KPf/fC RNA-seq, data analysis, and visualization
Stem and non-stem tumor cell isolation followed by RNA-sequencing
Tumors from three independent 10–12 week old REM2-KPf/fC mice were harvested and dissociated into a single cell suspension as described above. Tumor cells were stained with FC block then 0.2 µg/106 cells anti-EpCAM APC (eBioscience). 70,00–100,00 Msi2+/EpCAM+ (stem) and Msi2-/EpCAM+ (non-stem) cells were sorted and total RNA was isolated using RNeasy Micro kit (Qiagen). Total RNA was assessed for quality using an Agilent Tapestation, and all samples had RIN 7.9. RNA libraries were generated from 65 ng of RNA using Illumina’s TruSeq Stranded mRNA Sample Prep Kit following manufacturer’s instructions, modifying the shear time to 5 minutes. RNA libraries were multiplexed and sequenced with 50 basepair (bp) single end reads (SR50) to a depth of approximately 30 million reads per sample on an Illumina HiSeq2500 using V4 sequencing chemistry.
RNA-seq analysis
RNA-seq fastq files were processed into transcript-level summaries using kallisto (Bray et al., 2016), an ultrafast pseudo-alignment algorithm with expectation maximization. Transcript-level summaries were processed into gene-level summaries by adding all transcript counts from the same gene. Gene counts were normalized across samples using DESeq normalization (Anders and Huber 2010) and the gene list was filtered based on mean abundance, which left 13,787 genes for further analysis. Differential expression was assessed with an R package limma (Ritchie et al., 2015) applied to log2-transformed counts. Statistical significance of each test was expressed in terms of local false discovery rate lfdr (Efron and Tibshirani, 2002) using the limma function eBayes (Lönnstedt, I., and Speed, T. 2002). lfdr, also called posterior error probability, is the probability that a particular gene is not differentially expressed, given the data.
Cell state analysis
For cell state analysis, Gene Set Enrichment Analysis (GSEA) (Subramanian et al., 2005) was performed with the Bioconductor GSVA (Hänzelmann et al., 2013) and the Bioconductor GSVAdata c2BroadSets gene set collection, which is the C2 collection of canonical gene sets from MsigDB3.0 (Subramanian et al., 2005). Briefly, GSEA evaluates a ranked gene expression data-set against previously defined gene sets. GSEA was performed with the following parameters: mx.diff=TRUE, verbose=TRUE, parallel.sz=1, min.sz=5, max.sz=500, rnaseq=F.
Primary Msi2+ and Msi2- KPf/fC ChIP-seq for histone H3K27ac
Stem and non-stem tumor cell isolation followed by H3K27ac ChIP-sequencing
70,000 Msi2+/EpCAM+ (stem) and Msi2-/EpCAM+ (non-stem) cells were freshly isolated from a single mouse as described above. ChIP was performed as described previously (Deshpande et al., 2014); cells were pelleted by centrifugation and crosslinked with 1% formalin in culture medium using the protocol described previously (Deshpande et al., 2014). Fixed cells were then lysed in SDS buffer and sonicated on a Covaris S2 ultrasonicator. The following settings were used: Duty factor: 20%, Intensity: 4 and 200 Cycles/burst, Duration: 60 seconds for a total of 10 cycles to shear chromatin with an average fragment size of 200–400 bp. ChIP for H3K27Acetyl was performed using the antibody ab4729 (Abcam, Cambridge, UK) specific to the H3K27Ac modification. Library preparation of eluted chromatin immunoprecipitated DNA fragments was performed using the NEBNext Ultra II DNA library prep kit (E7645S and E7600S- NEB) for Illumina as per the manufacturer’s protocol. Library prepped DNA was then subjected to single-end, 75-nucleotide reads sequencing on the Illumina NexSeq500 sequencer at a sequencing depth of 20 million reads per sample.
H3K27ac signal quantification from ChIP-seq data
Pre-processed H3K27ac ChIP sequencing data was aligned to the UCSC mm10 mouse genome using the Bowtie2 aligner (version 2.1.0 (Langmead and Salzberg, 2012), removing reads with quality scores of <15. Non-unique and duplicate reads were removed using samtools (version 0.1.16, Li et al., 2009) and Picard tools (version 1.98), respectively. Replicates were then combined using BEDTools (version 2.17.0). Absolute H3K27ac occupancy in stem cells and non-stem cells was determined using the SICER-df algorithm without an input control (version 1.1; (Zang et al., 2009), using a redundancy threshold of 1, a window size of 200bp, a fragment size of 150, an effective genome fraction of 0.75, a gap size of 200bp and an E-value of 1000. Relative H3K27ac occupancy in stem cells vs non-stem cells was determined as above, with the exception that the SICER-df-rb algorithm was used.
Determining the overlap between peaks and genomic features
Genomic coordinates for features such as coding genes in the mouse mm10 build were obtained from the Ensembl 84 build (Ensembl BioMart). The observed vs expected number of overlapping features and bases between the experimental peaks and these genomic features (datasets A and B) was then determined computationally using a custom python script, as described in (Cole et al., 2017). Briefly, the number of base pairs within each region of A that overlapped with each region of B was computed. An expected background level of expected overlap was determined using permutation tests to randomly generate >1000 sets of regions with equivalent lengths and chromosomal distributions to dataset B, ensuring that only sequenced genomic regions were considered. The overlaps between the random datasets and experimental datasets were then determined, and p values and fold changes were estimated by comparing the overlap occurring by chance (expected) with that observed empirically (observed). This same process was used to determine the observed vs expected overlap of different experimental datasets.
RNA-Seq/ChIP-Seq correlation
Overlap between gene expression and H3K27ac modification
Genes that were up- or down-regulated in stem cells were determined using the Cuffdiff algorithm, and H3K27ac peaks that were enriched or disfavoured in stem cells were determined using the SICER-df-rb algorithm. The H3K27ac peaks were then annotated at the gene level using the ‘ChippeakAnno’ (Zhu et al., 2010) and ‘org.Mm.eg.db’ packages in R, and genes with peaks that were either exclusively up-regulated or exclusively down-regulated (termed ‘unique up’ or ‘unique down’) were isolated. The correlation between upregulated gene expression and up-regulated H3K27ac occupancy, or down-regulated gene expression and down-regulated H3K27ac occupancy, was then determined using the Spearman method in R.
Creation of composite plots
Composite plots showing RNA expression and H3K27ac signal across the length of the gene were created. Up- and down-regulated RNA peaks were determined using the FPKM output values from Tophat2 (Kim et al., 2013), and up- and down-regulated H3K27ac peaks were determined using the SICER algorithm. Peaks were annotated with nearest gene information, and their location relative to the TSS was calculated. Data were then pooled into bins covering gene length intervals of 5%. Overlapping up/up and down/down sets, containing either up- or down-regulated RNA and H3K27ac, respectively, were created, and the stem and non-stem peaks within these sets were plotted in Excel.
Super-enhancer identification
Enhancers in stem and non-stem cells were defined as regions with H3K27ac occupancy, as described in Hnisz et al. 2013. Peaks were obtained using the SICER-df algorithm before being indexed and converted to .gff format. H3K27ac Bowtie2 alignments for stem and non-stem cells were used to rank enhancers by signal density. Super-enhancers were then defined using the ROSE algorithm, with a stitching distance of 12.5kb and a TSS exclusion zone of 2.5kb. The resulting super-enhancers for stem or non-stem cells were then annotated at the gene level using the R packages ‘ChippeakAnno’(Zhu et al., 2010) and ‘org.Mm.eg.db’, and overlapping peaks between the two sets were determined using ‘ChippeakAnno’. Super-enhancers that are unique to stem or non-stem cells were annotated to known biological pathways using the Gene Ontology (GO) overrepresentation analysis functionality of the tool WebGestalt (Wang et al., 2017).
Genome-wide CRISPR screen
CRISPR library amplification and viral preparation
The mouse GeCKO CRISPRv2 knockout pooled library (Sanjana et al., 2014) was acquired from Addgene (catalog# 1000000052) as two half-libraries (A and B). Each library was amplified according to the Zhang lab library amplification protocol (Sanjana et al., 2014) and plasmid DNA was purified using NucleoBond Xtra Maxi DNA purification kit (Macherey-Nagel). For lentiviral production, 24 x T225 flasks were plated with 21×106 293T each in 1x DMEM containing 10% FBS. 24 hours later, cells were transfected with pooled GeCKOv2 library and viral constructs. Briefly, media was removed and replaced with 12.5 ml warm OptiMEM (Gibco). Per plate, 200 µl PLUS reagent (Life Technologies), 10 µg library A, and 10 µg library B was mixed in 4 ml OptiMEM along with 10 µg pRSV/REV (Addgene), 10 µg pMDLg/pRRE (Addgene), and 10 µg pHCMVG (Addgene) constructs. Separately, 200 µl Lipofectamine (Life Technologies) was mixed with 4 ml OptiMEM. After 5 minutes, the plasmid mix was combined with Lipofectamine and left to incubate at room temperature for 20 minutes, then added dropwise to each flask. Transfection media was removed 22 hours later and replaced with DMEM containing 10% FBS, 5 mM MgCl2, 1 U/ml DNase (Thermo Scientific), and 20mM HEPES pH 7.4. Viral supernatants were collected at 24 and 48 hours, passaged through 0.45 µm filter (corning), and concentrated by ultracentrifugation at 20,000 rpm for 2 hours at 4°C. Viral particles were resuspended in DMEM containing 10% FBS, 5 mM MgCl2, and 20 mM HEPES pH 7.4, and stored at −80°C.
CRISPR screen in primary KPf/fC cells
3 independent primary REM2-KPf/fC cell lines were established as described above and maintained in DMEM containing 10% FBS, 1x non-essential amino acids, and 1x pen/strep. At passage 3, each cell line was tested for puromycin sensitivity and GeCKOv2 lentiviral titer was determined. At passage 5, 1.6×108 cells from each cell line were transduced with GeCKOv2 lentivirus at an MOI of 0.3. 48 hours after transduction, 1×108 cells were harvested for sequencing (“T0”) and 1.6×108 were re-plated in the presence of puromycin according to previously tested puromycin sensitivity. Cells were passaged every 3–4 days for 3 weeks; at every passage, 5×107 cells were re-plated to maintain library coverage. At 2 weeks post-transduction, cell lines were tested for sphere forming capacity. At 3 weeks, 3×107 cells were harvested for sequencing (“2D; cell essential genes”), and 2.6×107 cells were plated in sphere conditions as described above (“3D; stem cell essential genes”). After 1 week in sphere conditions, tumorspheres were harvested for sequencing.
Analysis of the 2D data sets revealed that while some genes were required for growth in 2D, other genes that were not (detectably) required for growth in 2D were still required for growth in 3D (for example, Rorc Sox4, Foxo1, Wnt1 and ROBO3). These findings suggested that growth in 3D is dependent on a distinct or additional set of pathways. Since only stem cells give rise to 3D spheres, targets within the 3D datasets were prioritized for subsequent analyses. Of the genes that significantly dropped out in 3D, some also dropped out in 2D either significantly or as a trend.
DNA isolation, library preparation, and sequencing
Cells pellets were stored at −20°C until DNA isolation using Qiagen Blood and Cell Culture DNA Midi Kit (13343). Briefly, per 1.5×107 cells, cell pellets were resuspended in 2 ml cold PBS, then mixed with 2 ml cold buffer C1 and 6 ml cold H2O, and incubated on ice for 10 minutes. Samples were pelleted 1300 x g for 15 minutes at 4°C, then resuspended in 1 ml cold buffer C1 with 3 ml cold H2O, and centrifuged again. Pellets were then resuspended in 5 ml buffer G2 and treated with 100 µl RNAse A (Qiagen 1007885) for 2 minutes at room temperature followed by 95 µl Proteinase K for 1 hour at 50°C. DNA was extracted using Genomic-tip 100/G columns, eluted in 50°C buffer QF, and spooled into 300 µl TE buffer pH 8.0. Genomic DNA was stored at 4°C. For sequencing, gRNAs were first amplified from total genomic DNA isolated from each replicate at T0, 2D, and 3D (PCR1). Per 50 µl reaction, 4 µg gDNA was mixed with 25 µl KAPA HiFi HotStart ReadyMIX (KAPA Biosystems), 1 µM reverse primer1, and 1 µM forward primer1 mix (including staggers). Primer sequences are available upon request. After amplification (98°C 20 seconds, 66°C 20 seconds, 72°C 30 seconds, × 22 cycles), 50 µl of PCR1 products were cleaned up using QIAquick PCR Purification Kit (Qiagen). The resulting ~200bp products were then barcoded with Illumina Adaptors by PCR2. 5 µl of each cleaned PCR1 product was mixed with 25 µl KAPA HiFi HotStart ReadyMIX (KAPA Biostystems), 10 µl H2O, 1 µM reverse primer2, and 1 µM forward primer2. After amplification (98°C 20 seconds, 72°C 45 seconds, × 8 cycles), PCR2 products were gel purified, and eluted in 30 µl buffer EB. Final concentrations of the desired products were determined and equimolar amounts from each sample was pooled for Next Generation Sequencing.
Processing of the CRISPR screen data
Sequence read quality was assessed using fastqc (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Prior to alignment, 5’ and 3’ adapters flanking the sgRNA sequences were trimmed off using cutadapt v1.11 (Martin, 2011) with the 5’-adapter TCTTGTGGAAAGGACGAAACACCG and the 3’ adapter GTTTTAGAGCTAGAAATAGCAAGTT, which came from the cloning protocols of the respective libraries deposited on Addgene (https://www.addgene.org/pooled-library/). Error tolerance for adapter identification was set to 0.25, and minimal required read length after trimming was set to 10 bp. Trimmed reads were aligned to the GeCKO mouse library using Bowtie2 (Langmead and Salzberg, 2012) in the --local mode with a seed length of 11, an allowed seed mismatch of 1 and the interval function set to ‘S,1,0.75’. After completion, alignments were classified as either unique, failed, tolerated or ambiguous based on the primary (‘AS’) and secondary (‘XS’) alignment scores reported by Bowtie2. Reads with the primary alignment score not exceeding the secondary score by at least 5 points were discarded as ambiguous matches. Read counts were normalized by using the “size-factor” method as described in (Li et al., 2014). All of this was done using implementations in the PinAPL-Py webtool (Spahn et al., 2017), with detailed code available at https://github.com/LewisLabUCSD/PinAPL-Py.
gRNA growth and decay analysis
We used a parametric method in which the cell population with damaged gene i grows as Ni (t) = Ni(0)e(α0+δi)t, where α0 is the growth rate of unmodified cells and δi is the change of the growth rate due to the gene deletion. Since the aliquot extracted at each time point is roughly the same and represents only a fraction of the entire population, the observed sgRNA counts ni do not correspond to Ni directly. The correspondence is only relative: if we define ci ≡ ni/∑knk as the compositional fraction of sgRNA species , the correspondence is ci = Ni/∑kNk. As a result, the exponential can only be determined up to a multiplicative constant, e–δti = A∙ci(0)/ cit . The constant is determined from the assumption that a gene deletion typically does not affect the growth rate. Mathematically, | 1 = A med[ci(0)/ci(t)]. We define the statistic that measures the effect of gene deletion as xi ≡ e−δit and calculate it for every gene i from
Since we were interested in genes essential for growth, we performed a single-tailed test for xi . We collected the three values of xi , one from each biological replicate, into a vector xi . A statistically significant effect would have all three values large (>1) and consistent. If xi were to denote position of a point in a threedimensional space, we would be interested in points that lie close to the body diagonal and far away from the origin. A suitable statistic is , where is the unit vector in the direction of the body diagonal and · denotes scalar product. A q-value (false discovery rate) for each gene was estimated as the number of s-statistics not smaller than si expected in the null model divided by the observed number of s-statistics not smaller than si in the data. The null model was simulated numerically by permuting gene labels in xi for every experimental replicate, independently of each other, repeated 103 times.
STRING Interactome Network Analysis
The results from the CRISPR 3D experiment were integrated with the RNA-seq results using a network approach. We identified likely CRISPR-essential genes by filtering to include genes which had a falsediscovery rate corrected p-value of less than 0.5, resulting in 94 genes. We chose a relaxed filter here because the following filtering steps would help eliminate false positives, and our network analysis method would help to amplify weak signals. These genes were further filtered in two ways: first, we included only genes which were expressed in the RNA-seq data (this resulted in 57 genes), and second, we further restricted by genes which had enriched expression in stem cells by >2 log fold change in the RNA-seq (this resulted in 10 genes). These results were used to seed the network neighborhood exploration. We used the STRING mouse interactome (Szklarczyk et al., 2015) as our background network, including only high confidence interactions (edge weight > 700). The STRING interactome contains known and predicted functional protein-protein interactions. The interactions are assembled from a variety of sources, including genomic context predictions, high throughput lab experiments, and co-expression databases. Interaction confidence is a weighted combination of all lines of evidence, with higher quality experiments contributing more. The high confidence STRING interactome contains 13,863 genes, and 411,296 edges. Because not all genes are found in the interactome, our seed gene sets were further filtered when integrated with the network. This resulted in 39 CRISPR-essential, RNA-expressed seed genes, and 5 CRISPR-essential, RNA differentially-expressed seed genes. After integrating the seed genes with the background interactome, we employed a network propagation algorithm to explore the network neighborhood around these seed genes. Network propagation is a powerful method for amplifying weak signals by taking advantage of the fact that genes related to the same phenotype tend to interact. We implemented the network propagation method developed in (Vanunu et al., 2010), which simulates how heat would diffuse, with loss, through the network by traversing the edges, starting from an initially hot set of ‘seed’ nodes. At each step, one unit of heat is added to the seed nodes, and is then spread to the neighbor nodes. A constant fraction of heat is then removed from each node, so that heat is conserved in the system. After a number of iterations, the heat on the nodes converges to a stable value. This final heat vector is a proxy for how close each node is to the seed set. For example, if a node was between two initially hot nodes, it would have an extremely high final heat value, and if a node was quite far from the initially hot seed nodes, it would have a very low final heat value. This process is described by the following as in (Vanunu et al., 2010):
Where Ft is the heat vector at time t, Y is the initial value of the heat vector, W’ is the normalized adjacency matrix, and α ∈ (0,1) represents the fraction of total heat which is dissipated at every timestep. We examined the results of the subnetwork composed of the 500 genes nearest to the seed genes after network propagation. This is referred to as the ‘hot subnetwork’. In order to identify pathways and biological mechanisms related to the seed genes, we applied a clustering algorithm to the hot subnetwork, which partitioned the network into groups of genes which are highly interconnected within the group, and sparsely connected to genes in other groups. We used a modularity maximization algorithm for clustering (Blondel et al., 2008), which has proven effective in detecting modules, or clusters, in protein-protein interaction networks. These clusters were annotated to known biological pathways using the over-representation analysis functionality of the tool WebGestalt (Wang et al., 2017). We used the 500 genes in the hot subnetwork as the background reference gene set. To display the networks, we used a spring-embedded layout, which is modified by cluster membership (along with some manual adjustment to ensure non-overlapping labels) (Figure 2E). Genes belonging to each cluster were laid out radially along a circle, to emphasize the within cluster and between cluster connections. VisJS2jupyter (Rosenthal et al., 2018) was used for network propagation and visualization. Node color is mapped to the RNA-seq log fold change, with down-regulated genes displayed in blue, upregulated genes displayed in red, and genes with small fold changes displayed in gray. Labels are shown for genes which have a log fold change with absolute value greater than 3.0. Seed genes are shown as triangles with white outlines, while all other genes in the hot subnetwork are circles. The clusters have been annotated by selecting representative pathways from the enrichment analysis.
KPR172HC single cell analysis
Freshly harvested tumors from two independent KPR172hC mice were subjected to mechanical and enzymatic dissociation using a Miltenyi gentleMACS Tissue Dissociator to obtain single cells. The 10X Genomics Chromium Single Cell Solution was employed for capture, amplification and labeling of mRNA from single cells and for scRNA-Seq library preparation. Sequencing of libraries was performed on a Illumina HiSeq 2500 system. Sequencing data was input into the Cell Ranger analysis pipeline to align reads and generate gene-cell expression matrices. Finally, Custom R packages were used to perform gene-expression analyses and cell clustering projected using the t-SNE (t-Distributed Stochastic Neighbor Embedding) clustering algorithm. scRNA-seq datasets from the two independent KPR127hC tumor tissues generated on 10xGenomics platform were merged and utilized to explore and validate the molecular signatures of the tumor cells under dynamic development. The tumor cells that were used to illustrate the signal of Il10rb, Il34 and Csf1r etc. were characterized from the heterogeneous cellular constituents using SuperCT method developed by Dr. Wei Lin and confirmed by the Seurat FindClusters with the enriched signal of Epcam, Krt19 and Prom1 etc (Xie et al., 2018). The tSNE layout of the tumor cells was calculated by Seurat pipeline using the single-cell digital expression profiles.
KPf/fC single cell analysis
Three age-matched KPf/fC pancreatic tumors were collected and freshly dissociated, as described above. Tumor cells were stained with rat anti-mouse CD45-PE/Cy7 (eBioscience), rat anti-mouse CD31-PE (eBioscience), and rat anti-mouse PDGFRα-PacBlue (eBioscience) and tumor cells negative for these three markers were sorted for analysis. Individual cells were isolated, barcoded, and libraries were constructed using the 10x genomics platform using the Chromium Single Cell 3’ GEM library and gel bead kit v2 per manufacturer’s protocol. Libraries were sequenced on an Illumina HiSeq4000. The Cell Ranger software was used for alignment, filtering and barcode and UMI counting. The Seurat R package was used for further secondary analysis using default settings for unsupervised clustering and cell type discovery.
shRorc vs. shCtrl KPf/fC RNA-seq
Primary WT-KPf/fC cell lines were established as described above. WT-KPf/fC cells derived from an individual low passage cell line (<6 passage) were plated and transduced in triplicate with lentiviral particles containing shCtrl or shRorc. Positively infected (red) cells were sorted 5 days after transduction. Total RNA was isolated using the RNeasy Micro Plus kit (Qiagen). RNA libraries were generated from 200 ng of RNA using Illumina’s TruSeq Stranded mRNA Sample Prep Kit (Illumina) following manufacturer’s instructions. Libraries were pooled and single end sequenced (1X75) on the Illumina NextSeq 500 using the High output V2 kit (Illumina Inc., San Diego CA).
Read data was processed in BaseSpace (basespace.illumina.com). Reads were aligned to Mus musculus genome (mm10) using STAR aligner (https://code.google.com/p/rna-star/) with default settings. Differential transcript expression was determined using the Cufflinks Cuffdiff package (Trapnell et al., 2012) (https://github.com/cole-trapnell-lab/cufflinks). Differential expression data was then filtered to represent only significantly differentially expressed genes (q value<0.05). This list was used for pathway analysis and heatmaps of specific significantly differentially regulated pathways.
shRorc vs. shCtrl KPf/fC ChIP-seq for histone H3K27ac
Primary WT-KPf/fC cell lines were established as described above. Low passage (<6 passages) WT-KPf/fC cells from two independent cell lines were plated and transduced in triplicate with lentiviral particles containing shCtrl or shRorc. Positively infected (red) cells were sorted 5 days after transduction. ChIP-seq for histone H3K27-ac, signal quantification, and determination of the overlap between peaks and genomic features was conducted as described above.
Super-enhancers in control and shRorc-treated KPf/fC cell lines as well as Musashi stem cells were determined from H3K27ac ChIPseq data using the ROSE algorithm (http://younglab.wi.mit.edu/super_enhancer_code.html). The Musashi stem cell super-enhancer peaks were then further refined to include only those unique to the stem cell state (defined as present in stem cells but not nonstem cells) and/or those with RORγ binding sites within the peaks. Peak sequences were extracted using the ‘getSeq’ function from the ‘BSGenome.MMusculus.UCSC.mm10’ R package. RORγ binding sites were then mapped using the matrix RORG_MOUSE.H10MO.C.pcm (HOCOMOCO database) as a reference, along with the ‘matchPWM’ function in R at 90% stringency. Baseline peaks were then defined for each KPf/fC cell line as those overlapping each of the four Musashi stem cell peaklists with each KPC control SE list, giving eight in total. The R packages ‘GenomicRanges’ and ‘ChIPpeakAnno’ were used to assess peak overlap with a minimum overlap of 1bp used. To estimate the proportion of super-enhancers that are closed on RORC knockdown, divergence between each baseline condition and the corresponding KPf/fC shRorc super-enhancer list was assessed by quantifying the peak overlap and then expressing this as a proportion of the baseline list (‘shared%’). The proportion of unique peaks in each condition was then calculated as 100%-shared% and plotted.
sgRORC vs sgNT human RNA-seq
Human FG cells were plated and transduced in triplicate with lentiviral particles containing Cas9 and nontargeting guide RNA or guide RNA against Rorc. Positively infected (green) cells were sorted 5 days after transduction. Total RNA was isolated using the RNeasy Micro Plus kit (Qiagen). RNA libraries were generated from 200 ng of RNA using Illumina’s TruSeq Stranded mRNA Sample Prep Kit (illumina) following manufacturer’s instructions. Libraries were pooled and single end sequenced (1X75) on the Illumina NextSeq 500 using the High output V2 kit (Illumina Inc., San Diego CA).
Comparative RNA-seq and cell state analysis
RORC knockdown and control RNA-seq fastq files in mouse KPf/fC and human FG cells were processed into transcript-level summaries using kallisto (Bray et al., 2016). Transcript-level summaries were processed into gene-level summaries and differential gene expression was performed using sleuth with the Wald test (Pimentel et al., 2017). GSEA was performed as detailed above (Subramanian et al., 2005). Gene ontology analysis was performed using Metascape using a custom analysis with GO biological processes and default settings with genes with a FDR < 5% and a beta value >0.5.
cBioportal
RORC genomic amplification data from cancer patients was collected from the Memorial Sloan Kettering Cancer Center cBioPortal for Cancer Genomics (http://www.cbioportal.org).
Quantification and Statistical Analysis
Statistical analyses were carried out using GraphPad Prism software version 7.0d (GraphPad Software Inc.). Sample sizes for in vivo drug studies were determined based on the variability of pancreatic tumor models used. For flank transplant and autochthonous drug studies, tumor bearing animals within each group were randomly assigned to treatment groups. Treatment sizes were determined based on previous studies (Fox et al., 2016). Data are shown as the mean ± SEM. Two-tailed unpaired Student’s t-tests with Welch’s correction or One-way analysis of variance (ANOVA) for multiple comparisons when appropriate were used to determine statistical significance (*P<0.05, **P<0.01, ***P<0.001, ****P<0.0001).
The level of replication for each in vitro and in vivo study is noted in the figure legends for each figure and described in detail in the Method Details section above. However to summarize briefly, in vitro tumorsphere or colony formation studies were conducted with n=3 independent wells per cell line across two independent shRNA of n=3 wells; however, the majority of these experiments were additionally completed in >1 independently derived cell line, n=3 wells per shRNA. For limiting dilution assays, organoids were derived from 3 independent mice; drug-treated mouse and human organoids were plated at n=3 wells per dose per treatment condition. Flank shRNA studies were conducted twice independently, with n=4 tumors per group in each experiment. Flank drug studies were conducted at n=2–7 tumors per treatment group; autochthonous KPf/fC survival studies were conducted with a minimum of 4 mice enrolled in each treatment group. Live imaging studies were carried out with two mice per treatment group.
Statistical considerations and bioinformatic analysis of large data-sets generated are explained in great detail above. In brief, primary KPf/fC RNA-seq was performed using Msi2+ and Msi2- cells sorted independently from three different end-stage KPf/fC mice. Primary KPf/fC ChIP-seq was performed using Msi2+ and Msi2- cells sorted from an individual end-stage KPf/fC mouse. The genome-wide CRISPR screen was conducted using three biologically independent cell lines (derived from three different KPf/fC tumors). Singlecell analysis of tumors represents merged data from ~10,000 cells across two KPR172HC and three KPf/fC mice. RNA-seq for shRorc and shCtrl KPf/fC cells was conducted in triplicate, while ChIP-seq was conducted in single replicates from two biologically independent KPf/fC cell lines.
Supplementary Material
(A) Tumor organoid formation from primary isolated Musashi2+ (REM2+) and Musashi2- (REM2-) KPf/fC tumor cells. Number of cells plated is indicated above representative images, scale=200um.
(B) Limiting dilution frequency (left) calculated for Msi2+ (black) and Msi2- (red) organoid formation. Table (right) indicates cell doses tested in biological replicates.
(C-D) Gene set enrichment analysis (GSEA) of stem and non-stem gene signatures. Cell states (C), and corresponding heat-maps (D) of selected genes related to cell cycle. (C) Red denotes overlapping gene signatures; blue denotes non-overlapping gene signatures. (D) Red, over-represented gene expression; blue, under-represented gene expression; shades denote fold change from median values.
(E) Frequency of proliferating (Ki67+) Msi2+ (left) and Msi2- (right) tumor cells in untreated 10–12 week old REM2-KPf/fC mice (n=3), or treated with gemcitabine for 72 hours (n=1) or 6 days (n=1) prior to analysis; 200 mg/kg gemcitabine i.p. was delivered every 72 hours.
(F) Overlap of H3K27ac peaks and genomic features. For each genomic feature, frequency of H3K27ac peaks in stem cells (blue) and non-stem cells (gray) are represented as ratio of observed peak distribution/expected random genomic distribution.
(G and H) Concordance of H3K27ac peaks with RNA expression in stem cells (G; p=7.1×10−14) and non-stem cells (H; p<22×10−16).
(I and J) Ratio of observed/expected overlap in gene expression and H3K27ac enrichment comparing stem and non-stem cells. Down/Up, gene expression enriched in non-stem/H3K27ac enriched in stem; Up/Down, gene expression enriched in stem/H3K27ac enriched in non-stem; Down/Down, both gene expression and H3K27ac enriched in non-stem; Up/Up, both gene expression and H3K27ac enriched in stem.
Output of super-enhancer analysis on Msi2+ and Mis2- H3K27-acetyl ChIP-seq; all identified super-enhancers unique to stem-cells, unique to non-stem cells, and shared between stem and non-stem cells are listed in separated tabs.
shRNA target sequences and qPCR primer sequences used for each gene are listed here as separate tabs. Average knockdown efficiency is also listed for each shRNA construct per gene for in vitro experiments.
Table S3. Selected genes from stem cell networks. Related to Figure 2.
Table S4. Selected novel genes in pancreatic cancer. Related to Figure 3.
Table S5. Clinical and compound tool antagonists. Related to Figures 3 and 4.
Table S6. PDAC patients' characteristics (n=116). Related to Figure 7.
Complete output of gene set enrichment analysis conducted on stem and non-stem cell RNA-seq data-sets; gene signatures are include along with fold enrichment in stem cells, significance, and specific genes included in each signature.
(A) Establishment of three independent REM2-KPf/fC cell lines from end-stage REM2-KPf/fC mice for genomewide CRISPR-screen analysis. Stem cell content of freshly-dissociated REM2-KPf/fC tumors (A, left), and after puromycin selection in standard growth conditions (A, right).
(B and C) Volcano plots of guides enriched in 2D (B, tumor suppressors) and 3D (C, negative regulators of stem cells). Genes indicated on plots, p<0.005.
(D) Network propagation analysis integrating transcriptomic, epigenetic and functional analysis of stem cells. Genes enriched in stem cells by RNA-seq (stem/non-stem log2 fold-change >2) and depleted in 3D stem cell growth conditions (FDR<0.5) were used to seed the network (triangles), then analyzed for known and predicted protein-protein interactions and restricted to genes enriched in stem cells by RNA-seq (stem/non-stem log2 foldchange >2). Each node represents a single gene; node color is mapped to the RNA-seq fold change; stem cell enriched genes in red. Seven core programs were defined by groups of genes with high interconnectivity; each core program is annotated by Gene Ontology analysis (FDR<0.05).
(A and B) Sphere forming capacity of KPf/fC cells following shRNA knockdown. Selected genes involved in stem and developmental processes (A) or cell adhesion, cell motility, and matrix components (B).
(C-D) Immunofluorescence analysis of Celsr1 (C) and Celsr2 (D) in EpCAM+ stem (CD133+) and non-stem (CD133-) primary tumor cells isolated from KPf/fC mice. Three frames were analyzed per slide, and the frequency of Celsr1-high or Celsr2-high cells determined, scale=25um.
(E) KPf/fC cells were infected with shRNA against Pear1and protein knockdown efficiency determined five days post-transduction by western blot.
(F-H) Independent replicates for impact of shRNA inhibition of target genes on tumor growth in vivo. Celsr1 (F), Celsr2 (G), and Pear1 (H) were inhibited via shRNA delivery in KPf/fC cells, and impact on tumor growth assessed by tracking flank transplants in vivo, n=4 per condition.
(I) Pear1 was inhibited via shRNA in REM-KPf/fC cells in sphere culture and impact on Msi+ stem cell content assessed by FACS, n=3 per condition, p = 0.0629.
(J) Pear1 was inhibited via shRNA in KPf/fC cells and impact on apoptosis in sphere culture as marked by Annexin-V assessed by FACS, n=3 per condition.
(K) Heat map of relative RNA expression of cytokines and related receptors in KPf/fC stem and non-stem cells (left) and average RNA-seq TPM values in Msi2- and Msi2+ cells (right). Red, over-represented; blue, underrepresented; color denotes fold change from median values.
(L) Single cell RNA Sequencing maps of KPR172H/+C tumors. Tumor cells defined by expression of EpCAM (far left), Krt19 (left center), Cdh1 (right center), and Cdh2 (far right).
(M) Left, KPR172H/+C tumor single-cell sequencing map of cells expressing Msi2 within the EpCAM+ tumor cell fraction. Right, KPR172H/+C tumor single-cell sequencing map of cells expressing IL10Rβ, IL34, and CSF1R within the EpCAM+Msi2+ stem cell fraction.
(N-O) Independent replicates for impact of shRNA inhibition of target genes on tumor growth in vivo. IL10Rβ, (N) and CSF1R (O) were inhibited via shRNA delivery in KPf/fC cells, and impact on tumor growth assessed by tracking flank transplants in vivo, n=4 per condition.
(P) Cytokine receptors IL10Rβ and CSF1R were inhibited by shRNA delivery in KPf/fC cells and plated in sphere culture for one week. Increased apoptosis in KPf/fC cells with shIL10Rb (p<.05) and shCSF1R (trend). Frequency of apoptotic cells determined by Annexin-V staining and FACS analysis, n=3 per condition.
(Q) Representative FACS plots for stem content analysis in vitro. IL-10rβ and Csf1R were inhibited via shRNA delivery in KPf/fC cells, and impact on stem content (Msi2-GFP+ cells) in sphere culture assessed by FACS, n=3 per condition.
(R) ELISA based quantification (Quantikine, R&D Systems) of IL-10, IL-34, and CSF-1 in media (left) and KPf/fC cell lystate (right). Cytokines were quantified in fresh sphere culture media, KPf/fC stem and non-stem cell conditioned media, and KPf/fC epithelial cell lysate. Conditioned media was generated by culturing sorted CD133- or CD133+ KPf/fC cells in sphere media for 48 hours; media was filtered and assayed immediately. Cell lysate was collected in RIPA buffer and assayed at 2 mg/mL for ELISA. n=3 per condition.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01 by Student’s t-test or One-way ANOVA.
(A) Heat map of transcription factors in KPf/fC stem and non-stem identified as possible pancreatic cancer stem cell dependencies within the network map (see Figure 2E). Red, over-represented; blue, under-represented; color denotes fold change from median values.
(B) Distribution of RORγ consensus binding sites in genomic regions associated with H3K27ac. Down/Down, both gene expression and H3K27ac enriched in non-stem cells; Up/Up, both gene expression and H3K27ac enriched in stem cells.
(C) Biological replicates showing qPCR analysis of RORγ expression in primary KPf/fC stem and non-stem tumor cells isolated from REM2-KPf/fC mice.
(D) Immunofluorescence analysis of RORγ in primary KPf/fC EpCAM+ CD133+ and CD133- tumor cells. Three frames were analyzed per slide, and the frequency of RORγ -high cells determined.
(E) KPf/fC tumor single-cell sequencing map of cells expressing RORγ within the EpCAM+Msi2+ cell fraction (n=3 mice represented).
(F) RORγ expression within E-Cadherin- stromal cells in patient samples.
(G) IL1R1 was inhibited by CRISPR-mediated deletion in KPf/fC cells, and impact on Rorc expression assessed by qPCR. Two distinct guide RNAs (sgIL1r1–1 and sgIL1r1–2) were used to knockout IL1R1; expression was quantified by qPCR and is shown relative to control (non-targeting guide RNA), n=3 per condition.
(H) Knockdown efficiency of RORγ in KPf/fC cells infected with Rorc shRNA determined five days post-transduction. Relative expression in western blots quantified relative to tubulin loading control.
(I) Impact of shRNA-mediated RORγ inhibition on apoptosis and proliferation of in KPf/fC cells in 3D culture n=3.
(J) Independent replicate of shRNA Rorc impact on KPf/fC tumor propagation as assessed by tracking flank transplants in vivo, n=4 per condition.
(K-M) SE analysis of shRorc KPf/fC cells. KPf/fC cells were infected with shRorc, and used for H3K27ac ChIPseq and SE analysis, schematic (K). H3K27ac peaks were analyzed to assess SE overlap in shCtrl and shRorc samples (L). SEs lost in shRorc samples were crossed to stem-enriched and stem-unique SEs identified in primary Msi2-GFP+ KPf/fC tumors cells, and further restricted to SEs containing RORγ binding motifs (M). Majority of super-enhancer landscape remained unchanged with RORγ loss, and landscape changes that did occur were not enriched in SEs with RORγ binding sites. ChIP-seq analysis was conducted in two independent KPf/fC cell lines.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
(A) Size of flank KPf/fC tumors in immunocompetent mice prior to enrollment into RORγ targeted therapy. Group 1, vehicle; group 2, SR2211; group 3, vehicle + gemcitabine; group 4, SR2211 + gemcitabine.
(B) Target engagement following acute RORγ inhibition in vivo. 9.5 wk tumor-bearing KPf/fC mice were treated with vehicle or SR2211 for two weeks (midpoint), after which tumors were isolated, fixed, and analyzed for target engagement of Hmga2 in epithelial cells by immunofluorescence. Representative images (left) and quantification (right) of Hmga2+ Keratin+ epithelial cells in vehicle or SR2211 treated tumors. Four frames were analyzed per mouse, n=2–4 mice per condition, Hmga2 (red), Keratin (green), scale=25um.
(C) Target engagement in endpoint tumors following continuous RORγ inhibition in vivo. 8 wk tumor-bearing KPf/fC mice were treated till endpoint with either vehicle or SR2211, after which tumors were isolated, fixed, and analyzed for target engagement of Hmga2 in epithelial cells by immunofluorescence. Representative images (left) and quantification (right) of Hmga2+ Keratin+ epithelial cells in vehicle or SR2211 treated tumors. Four frames were analyzed per mouse, n=2–4 mice per condition, Hmga2 (red), Keratin (green).
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA. Grubb’s test (p=0.1). was used to remove an outlier from the midpoint SR2211 treated group, scale=25um.
(A-B) Analysis of T cell subsets in KPf/fC tumors transplanted into wild-type or Rorc-knockout recipient mice (vehicle treated groups shown). Absolute cell numbers of the following populations were evaluated: (A) CD45+/CD3+/CD8+ or CD8+ T cells, (B) CD45+/CD3+/CD4+ or CD4+ T cells.
(C-L) FACS analysis of non-neoplastic cell populations in autochthonous tumors from KPf/fC mice treated with vehicle or SR2211 for 1 week. Schematic (C). Absolute cell numbers of the following populations were evaluated: CD45+ cells (D), CD11b+/F480+ cells (macrophage) (E), CD11b+/Gr-1+ cells (MDSC) (F), CD11c+ cells (dendritic) (G), CD45+/CD3+ T cells (H), CD3+/CD8+ T cells (I), CD3+/CD4+ T cells (J), CD4+/IL-17+ Th17 cells (K), CD31+ cells (endothelial) (L). (n=3 per condition).
(M) In vivo imaging of tumor vasculature of KPf/fC mice treated with vehicle or SR2211, vasculature is marked by in vivo delivery of anti-VE-Cadherin (magenta), scale=75um. Data represented as mean +/− S.E.M. * p<0.05 by Student’s t-test or One-way ANOVA.
(A-D) Gene ontology and gene set enrichment analysis of RNA-seq in human and mouse pancreatic cancer cells to identify common genes and pathways regulated by RORγ. Gene ontology analysis of KPf/fC RNA-seq showing genes downregulated with shRorc were enriched for cytokine-mediated signaling pathway GO term (A). Differentially expressed genes in KPf/fC within cytokine-mediated signaling pathway (B) were crossed with differentially expressed genes identified by RNA-seq analysis of human pancreatic cancer cells (FG) where RORC was knocked out using CRISPR. Gene set enrichment analysis of mouse and human RNA-seq shows common cytokine gene sets regulated by Rorc across species (D).
(E) Analysis of CRISPR guide depletion in stem cell conditions for SE associated genes expressed in stem or non-stem cells.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
Highlights.
Map of PDAC dependencies using RNA-seq, ChIP-seq and genome-wide CRISPR screening
Expression and direct utilization of cytokine/immune signals in PDAC stem cells
Nuclear hormone receptor RORγ regulates mouse and human pancreatic cancer
Pharmacologic blockade of RORγ reduces tumor burden and improves survival
Acknowledgments
We are grateful to Olivier Harismendy, Prashant Mali and Kristen Jepsen for help with the CRISPR screen design, Dan Littman for comments on the manuscript, Michael Karin and Ron Evans for scientific advice, Christopher Wright for providing the Ptf1a-Cre mice, and Armin Ahmadi and Kendall Chambers for technical support. N.K.L. received support from T32 GM007752 and an NRSA F31 CA206416; L.P.F received support from T32 GM007752. The project was partially supported by the NIH Grant UL1TR001442, as well as by CA169281 to H.H. and D.V.H., a grant from the National Foundation for Cancer Research to D.V.H., R35 GM119850 and NNF10CC1016517 to N.E.L., CA155620 to A.M.L, and R35 CA197699 and R01 CA186043 to T.R. This work was also supported by a SU2C–CRUK–Lustgarten Foundation Pancreatic Cancer Dream Team Research Grant (SU2C-AACR-DT-20–16) to D.V.H, H.H., A.M.L., and T.R, and a SU2C-Lustgarten Foundation pancreatic cancer collective grant to A.M.L. and T.R. The content is solely the responsibility of the authors and does not necessarily represent the official views of the NIH.
Footnotes
Declaration of Interests
This work includes data filed as a provisional application SD2019–243. TR is a consultant for Orphagen Pharmaceuticals, and DVH serves on the Scientific Advisory Board of Five Prime and is a consultant for Celgene.
Data and Software Availability
The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request. Singe cell, Genome-wide CRISPR screen, H3K27ac ChIP, and RNA sequencing data have been deposited at NCBI GEO:
Primary Msi2+ and Msi2− KPf/fC RNA-seq https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE114906
Primary Msi2+ and Msi2− KPf/fC ChIP-seq for histone H3K27ac https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE113712
Genome-wide CRISPR screen https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE114914
shRorc vs shControl KPf/fC ChIP-seq for histone H3K27ac https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE126536
shRorc vs. shCtrl KPf/fC RNA-seq https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE126538
sgRORC vs sgNT human RNA-seq https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE126537
KPf/fC single cell analysis https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE126539
KPR172HC single cell analysis https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE126388
Code availability
Custom code developed for CRISPR screen analysis and network propagation were deposited to github.com and can be accessed at https://github.com/ucsd-ccbb/crispr_network_analysis.
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
References
- Anders S and Huber W (2010). Differential expression analysis for sequence count data. Genome biology, 11(10), p.R106. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Bailey JM, Alsina J, Rasheed ZA, McAllister FM, Fu YY, Plentz R, Zhang H, Pasricha PJ, Bardeesy N, Matsui W, et al. (2014). DCLK1 marks a morphologically distinct subpopulation of cells with stem cell properties in preinvasive pancreatic cancer. Gastroenterology 146, 245–256. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Bailey P, Chang DK, Nones K, Johns AL, Patch AM, Gingras MC, Miller DK, Christ AN, Bruxner TJ, Quinn MC, et al. (2016). Genomic analyses identify molecular subtypes of pancreatic cancer. Nature 531, 47–52. [DOI] [PubMed] [Google Scholar]
- Blondel VD, J.L., G., R., L., and E., L. (2008). Fast unfolding of communities in large networks. In Journal of Statistical Mechanics (P10008), pp. arXiv:0803.0476.
- Boj SF, Hwang CI, Baker LA, Chio II, Engle DD, Corbo V, Jager M, Ponz-Sarvise M, Tiriac H, Spector MS, et al. (2015). Organoid models of human and mouse ductal pancreatic cancer. Cell 160, 324–338. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Bray NL, Pimentel H, Melsted P, and Pachter L (2016). Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol 34, 525–527. [DOI] [PubMed] [Google Scholar]
- Cole JJ, Robertson NA, Rather MI, Thomson JP, McBryan T, Sproul D, Wang T, Brock C, Clark W, Ideker T, et al. (2017). Diverse interventions that extend mouse lifespan suppress shared age-associated epigenetic changes at critical gene regulatory regions. Genome Biol 18, 58. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Conroy T, Desseigne F, Ychou M, Bouche O, Guimbaud R, Becouarn Y, Adenis A, Raoul JL, Gourgou-Bourgade S, de la Fouchardiere C, et al. (2011). FOLFIRINOX versus gemcitabine for metastatic pancreatic cancer. N Engl J Med 364, 1817–1825. [DOI] [PubMed] [Google Scholar]
- Cook DN, Kang HS, and Jetten AM (2015). Retinoic Acid-Related Orphan Receptors (RORs): Regulatory Functions in Immunity, Development, Circadian Rhythm, and Metabolism. Nucl Receptor Res 2. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Deshpande AJ, Deshpande A, Sinha AU, Chen L, Chang J, Cihan A, Fazio M, Chen CW, Zhu N, Koche R, et al. (2014). AF10 regulates progressive H3K79 methylation and HOX gene expression in diverse AML subtypes. Cancer Cell 26, 896–908. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M and Gingeras TR (2013). STAR: ultrafast universal RNA-seq aligner. Bioinformatics, 29(1), pp.15–21. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Dull T, Zufferey R, Kelly M, Mandel RJ, Nguyen M, Trono D and Naldini L (1998). A thirdgeneration lentivirus vector with a conditional packaging system. Journal of virology, 72(11), pp.8463–8471. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Dwyer JR, Donkor J, Zhang P, Csaki LS, Vergnes L, Lee JM, Dewald J, Brindley DN, Atti E, Tetradis S, et al. (2012). Mouse lipin-1 and lipin-2 cooperate to maintain glycerolipid homeostasis in liver and aging cerebellum. Proc Natl Acad Sci U S A 109, E2486–2495. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Efron B, and Tibshirani R (2002). Empirical bayes methods and false discovery rates for microarrays. Genet Epidemiol 23, 70–86. [DOI] [PubMed] [Google Scholar]
- Fox RG, Lytle NK, Jaquish DV, Park FD, Ito T, Bajaj J, Koechlein CS, Zimdahl B, Yano M, Kopp J, et al. (2016). Image-based detection and targeting of therapy resistance in pancreatic adenocarcinoma. Nature 534, 407–411. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Gege C (2016). Retinoid-related orphan receptor gamma t (RORgammat) inhibitors from Vitae Pharmaceuticals (WO2015116904) and structure proposal for their Phase I candidate VTP-43742. Expert Opin Ther Pat 26, 737–744. [DOI] [PubMed] [Google Scholar]
- Grosse-Wilde A, Fouquier d’Herouel A, McIntosh E, Ertaylan G, Skupin A, Kuestner RE, del Sol A, Walters KA, and Huang S (2015). Stemness of the hybrid Epithelial/Mesenchymal State in Breast Cancer and Its Association with Poor Survival. PLoS One 10, e0126522. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Guillonneau C, Bezie S, and Anegon I (2017). Immunoregulatory properties of the cytokine IL-34. Cell Mol Life Sci 74, 2569–2586. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Hänzelmann S, Castelo R and Guinney J (2013). GSVA: gene set variation analysis for microarray and RNA-seq data. BMC bioinformatics, 14(1), p.7. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Hermann PC, Huber SL, Herrler T, Aicher A, Ellwart JW, Guba M, Bruns CJ, and Heeschen C (2007). Distinct populations of cancer stem cells determine tumor growth and metastatic activity in human pancreatic cancer. Cell Stem Cell 1, 313–323. [DOI] [PubMed] [Google Scholar]
- Hingorani SR, Petricoin EF, Maitra A, Rajapakse V, King C, Jacobetz MA, Ross S, Conrads TP, Veenstra TD, Hitt BA, et al. (2003). Preinvasive and invasive ductal pancreatic cancer and its early detection in the mouse. Cancer Cell 4, 437–450. [DOI] [PubMed] [Google Scholar]
- Hnisz D, Abraham BJ, Lee TI, Lau A, Saint-Andre V, Sigova AA, Hoke HA, and Young RA (2013). Super-enhancers in the control of cell identity and disease. Cell 155, 934–947. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Hu Y, and Smyth GK (2009). ELDA: extreme limiting dilution analysis for comparing depleted and enriched populations in stem cell and other assays. J Immunol Methods 347, 70–78. [DOI] [PubMed] [Google Scholar]
- Huh JR, and Littman DR (2012). Small molecule inhibitors of RORgammat: targeting Th17 cells and other applications. Eur J Immunol 42, 2232–2237. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Ivanov II, McKenzie BS, Zhou L, Tadokoro CE, Lepelley A, Lafaille JJ, Cua DJ, and Littman DR (2006). The orphan nuclear receptor RORgammat directs the differentiation program of proinflammatory IL17+ T helper cells. Cell 126, 1121–1133. [DOI] [PubMed] [Google Scholar]
- Kawaguchi Y, Cooper B, Gannon M, Ray M, MacDonald RJ, and Wright CV (2002). The role of the transcriptional regulator Ptf1a in converting intestinal to pancreatic progenitors. Nat Genet 32, 128–134. [DOI] [PubMed] [Google Scholar]
- Kawamoto M, Ishiwata T, Cho K, Uchida E, Korc M, Naito Z, and Tajiri T (2009). Nestin expression correlates with nerve and retroperitoneal tissue invasion in pancreatic cancer. Hum Pathol 40, 189–198. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R and Salzberg SL (2013). TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome biology, 14(4), p.R36. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Kopp JL, von Figura G, Mayes E, Liu FF, Dubois CL, Morris J.P.t., Pan FC, Akiyama H, Wright CV, Jensen K, et al. (2012). Identification of Sox9-dependent acinar-to-ductal reprogramming as the principal mechanism for initiation of pancreatic ductal adenocarcinoma. Cancer Cell 22, 737–750. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Kumar N, Lyda B, Chang MR, Lauer JL, Solt LA, Burris TP, Kamenecka TM, and Griffin PR (2012). Identification of SR2211: a potent synthetic RORgamma-selective modulator. ACS Chem Biol 7, 672677. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Langmead B, and Salzberg SL (2012). Fast gapped-read alignment with Bowtie 2. Nat Methods 9, 357–359. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Lappano R, and Maggiolini M (2011). G protein-coupled receptors: novel targets for drug discovery in cancer. Nat Rev Drug Discov 10, 47–60. [DOI] [PubMed] [Google Scholar]
- Li C, Heidt DG, Dalerba P, Burant CF, Zhang L, Adsay V, Wicha M, Clarke MF, and Simeone DM (2007). Identification of pancreatic cancer stem cells. Cancer Res 67, 1030–1037. [DOI] [PubMed] [Google Scholar]
- Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G and Durbin R (2009). The sequence alignment/map format and SAMtools. Bioinformatics, 25(16), pp.2078–2079. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Li C, Wu JJ, Hynes M, Dosch J, Sarkar B, Welling TH, Pasca di Magliano M, and Simeone DM (2011). c-Met is a marker of pancreatic cancer stem cells and therapeutic target. Gastroenterology 141, 22182227 e2215. [DOI] [PubMed] [Google Scholar]
- Li W, Xu H, Xiao T, Cong L, Love MI, Zhang F, Irizarry RA, Liu JS, Brown M, and Liu XS (2014). MAGeCK enables robust identification of essential genes from genome-scale CRISPR/Cas9 knockout screens. Genome Biol 15, 554. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Lönnstedt I, and Speed T (2002). Replicated microarray data. Statistica Sinica 12, 31–46. [Google Scholar]
- Lu H, Chen I, Shimoda LA, Park Y, Zhang C, Tran L, Zhang H, and Semenza GL (2017). Chemotherapy-Induced Ca(2+) Release Stimulates Breast Cancer Stem Cell Enrichment. Cell Rep 18, 19461957. [DOI] [PubMed] [Google Scholar]
- Macosko EZ, Basu A, Satija R, Nemesh J, Shekhar K, Goldman M, Tirsoh I, Bialas AR, Kamitaki N, Martersteck EM, et al. (2015). Highly parallel genome-wide expression profiling of individual cells using nanoliter droplets. Cell, 161(5), 1202–1214. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Martin M (2011). Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet. journal, 17(1), pp.10–12. [Google Scholar]
- Morgan RT, Woods LK, Moore GE, Quinn LA, McGavran L and Gordon SG (1980). Human cell line (COLO 357) of metastatic pancreatic adenocarcinoma. International journal of cancer, 25(5), pp.591–598. [DOI] [PubMed] [Google Scholar]
- Neph S, Stergachis AB, Reynolds A, Sandstrom R, Borenstein E, and Stamatoyannopoulos JA (2012). Circuitry and dynamics of human transcription factor regulatory networks. Cell 150, 1274–1286. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Olive KP, Tuveson DA, Ruhe ZC, Yin B, Willis NA, Bronson RT, Crowley D and Jacks T (2004). Mutant p53 gain of function in two mouse models of Li-Fraumeni syndrome. Cell, 119(6), pp.847–860 [DOI] [PubMed] [Google Scholar]
- Pimentel H, Bray NL, Puente S, Melsted P and Pachter L (2017). Differential analysis of RNA-seq incorporating quantification uncertainty. Nature methods, 14(7), p.687. [DOI] [PubMed] [Google Scholar]
- Rasheed ZA, Yang J, Wang Q, Kowalski J, Freed I, Murter C, Hong SM, Koorstra JB, Rajeshkumar NV, He X, et al. (2010). Prognostic significance of tumorigenic cells with mesenchymal features in pancreatic adenocarcinoma. J Natl Cancer Inst 102, 340–351. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Reya T, Morrison SJ, Clarke MF, and Weissman IL (2001). Stem cells, cancer, and cancer stem cells. Nature 414, 105–111. [DOI] [PubMed] [Google Scholar]
- Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, and Smyth GK (2015). limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res 43, e47. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Rosenthal SB, Len J, Webster M, Gary A, Birmingham A, and Fisch KM (2018). Interactive network visualization in Jupyter notebooks: visJS2jupyter. Bioinformatics 34, 126–128. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Rovira M, Scott SG, Liss AS, Jensen J, Thayer SP, and Leach SD (2010). Isolation and characterization of centroacinar/terminal ductal progenitor cells in adult mouse pancreas. Proc Natl Acad Sci U S A 107, 75–80. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Ryu CS, Kwak HC, Lee KS, Kang KW, Oh SJ, Lee KH, Kim HM, Ma JY, and Kim SK (2011). Sulfur amino acid metabolism in doxorubicin-resistant breast cancer cells. Toxicol Appl Pharmacol 255, 94–102. [DOI] [PubMed] [Google Scholar]
- Sanjana NE, Shalem O, and Zhang F (2014). Improved vectors and genome-wide libraries for CRISPR screening. Nat Methods 11, 783–784. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Sankhala KK, B.J.-Y., Ganjoo KN, Italiano A, Hassan AB, Kim TM, Ravi V, Cassier PA, Rutkowski P, Sankar N, Qazi I, et al. (2017). A phase I/II dose escalation and expansion study of cabiralizumab (cabira; FPA-008), an anti-CSF1R antibody, in tenosynovial giant cell tumor (TGCT, diffuse pigmented villonodular synovitis D-PVNS). J Clin Oncol 35, 11078. [Google Scholar]
- Sena-Esteves M, Tebbets JC, Steffens S, Crombleholme T and Flake AW (2004). Optimized large-scale production of high titer lentivirus vector pseudotypes. Journal of virological methods, 122(2), pp.131–139. [DOI] [PubMed] [Google Scholar]
- Spahn PN, Bath T, Weiss RJ, Kim J, Esko JD, Lewis NE, and Harismendy O (2017). PinAPL-Py: A comprehensive web-application for the analysis of CRISPR/Cas9 screens. Sci Rep 7, 15854. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, Paulovich A, Pomeroy SL, Golub TR, Lander ES, et al. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A 102, 15545–15550. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Szklarczyk D, Franceschini A, Wyder S, Forslund K, Heller D, Huerta-Cepas J, Simonovic M, Roth A, Santos A, Tsafou KP, et al. (2015). STRING v10: protein-protein interaction networks, integrated over the tree of life. Nucleic Acids Res 43, D447–452. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL and Pachter L (2012). Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nature protocols, 7(3), p.562. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Tuveson DA, Shaw AT, Willis NA, Silver DP, Jackson EL, Chang S, Mercer KL, Grochow R, Hock H, Crowley D, et al. (2004). Endogenous oncogenic K-ras(G12D) stimulates proliferation and widespread neoplastic and developmental defects. Cancer Cell 5, 375–387. [DOI] [PubMed] [Google Scholar]
- Van den Broeck A, Vankelecom H, Van Delm W, Gremeaux L, Wouters J, Allemeersch J, Govaere O, Roskams T, and Topal B (2013). Human pancreatic cancer contains a side population expressing cancer stem cell-associated and prognostic genes. PLoS One 8, e73968. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Van Laarhoven T, and Marchiori E (2012). Robust community detection methods with resolution parameter for complex detection in protein protein interaction networks. Pattern Recognition in Bioinformatics, 1–13.
- Vanunu O, Magger O, Ruppin E, Shlomi T, and Sharan R (2010). Associating genes and protein complexes with disease via network propagation. PLoS Comput Biol 6, e1000641. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Wang J, Vasaikar S, Shi Z, Greer M, and Zhang B (2017). WebGestalt 2017: a more comprehensive, powerful, flexible and interactive gene set enrichment analysis toolkit. Nucleic Acids Res 45, W130–W137. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Wang J, Zou JX, Xue X, Cai D, Zhang Y, Duan Z, Xiang Q, Yang JC, Louie MC, Borowsky AD, et al. (2016). ROR-gamma drives androgen receptor expression and represents a therapeutic target in castration-resistant prostate cancer. Nat Med 22, 488–496. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Wang X, Wong K, Ouyang W, and Rutz S (2017). Targeting IL-10 Family Cytokines for the Treatment of Human Diseases. Cold Spring Harb Perspect Biol [DOI] [PMC free article] [PubMed]
- Wartenberg M, Cibin S, Zlobec I, Vassella E, Eppenberger-Castori S, Terracciano L, Eichmann MD, Worni M, Gloor B, Perren A and Karamitopoulou E (2018). Integrated genomic and immunophenotypic classification of pancreatic cancer reveals three distinct subtypes with prognostic/predictive significance. Clinical cancer research, 24(18), pp.4444–4454. [DOI] [PubMed] [Google Scholar]
- Whyte WA, Orlando DA, Hnisz D, Abraham BJ, Lin CY, Kagey MH, Rahl PB, Lee TI, and Young RA (2013). Master transcription factors and mediator establish super-enhancers at key cell identity genes. Cell 153, 307–319. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Xie P, Gao M, Wang C, Noel P, Yang C, Von Hoff D, Han H, Zhang M and Lin W (2018). SuperCT: A supervised-learning-framework to enhance the characterization of single-cell transcriptomic profiles. bioRxiv, p.416719. [DOI] [PMC free article] [PubMed]
- Zang C, Schones DE, Zeng C, Cui K, Zhao K, and Peng W (2009). A clustering approach for identification of enriched domains from histone modification ChIP-Seq data. Bioinformatics 25, 1952–1958. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Zhao L, Spassieva S, Gable K, Gupta SD, Shi LY, Wang J, Bielawski J, Hicks WL, Krebs MP, Naggert J, et al. (2015). Elevation of 20-carbon long chain bases due to a mutation in serine palmitoyltransferase small subunit b results in neurodegeneration. Proc Natl Acad Sci U S A 112, 12962–12967 J. [DOI] [PMC free article] [PubMed] [Google Scholar]
- Zhu LJ, Gazin C, Lawson ND, Pagès H, Lin SM, Lapointe DS and Green MR (2010). ChIPpeakAnno: a Bioconductor package to annotate ChIP-seq and ChIP-chip data. BMC bioinformatics, 11(1), p.237. [DOI] [PMC free article] [PubMed] [Google Scholar]
Associated Data
This section collects any data citations, data availability statements, or supplementary materials included in this article.
Supplementary Materials
(A) Tumor organoid formation from primary isolated Musashi2+ (REM2+) and Musashi2- (REM2-) KPf/fC tumor cells. Number of cells plated is indicated above representative images, scale=200um.
(B) Limiting dilution frequency (left) calculated for Msi2+ (black) and Msi2- (red) organoid formation. Table (right) indicates cell doses tested in biological replicates.
(C-D) Gene set enrichment analysis (GSEA) of stem and non-stem gene signatures. Cell states (C), and corresponding heat-maps (D) of selected genes related to cell cycle. (C) Red denotes overlapping gene signatures; blue denotes non-overlapping gene signatures. (D) Red, over-represented gene expression; blue, under-represented gene expression; shades denote fold change from median values.
(E) Frequency of proliferating (Ki67+) Msi2+ (left) and Msi2- (right) tumor cells in untreated 10–12 week old REM2-KPf/fC mice (n=3), or treated with gemcitabine for 72 hours (n=1) or 6 days (n=1) prior to analysis; 200 mg/kg gemcitabine i.p. was delivered every 72 hours.
(F) Overlap of H3K27ac peaks and genomic features. For each genomic feature, frequency of H3K27ac peaks in stem cells (blue) and non-stem cells (gray) are represented as ratio of observed peak distribution/expected random genomic distribution.
(G and H) Concordance of H3K27ac peaks with RNA expression in stem cells (G; p=7.1×10−14) and non-stem cells (H; p<22×10−16).
(I and J) Ratio of observed/expected overlap in gene expression and H3K27ac enrichment comparing stem and non-stem cells. Down/Up, gene expression enriched in non-stem/H3K27ac enriched in stem; Up/Down, gene expression enriched in stem/H3K27ac enriched in non-stem; Down/Down, both gene expression and H3K27ac enriched in non-stem; Up/Up, both gene expression and H3K27ac enriched in stem.
Output of super-enhancer analysis on Msi2+ and Mis2- H3K27-acetyl ChIP-seq; all identified super-enhancers unique to stem-cells, unique to non-stem cells, and shared between stem and non-stem cells are listed in separated tabs.
shRNA target sequences and qPCR primer sequences used for each gene are listed here as separate tabs. Average knockdown efficiency is also listed for each shRNA construct per gene for in vitro experiments.
Table S3. Selected genes from stem cell networks. Related to Figure 2.
Table S4. Selected novel genes in pancreatic cancer. Related to Figure 3.
Table S5. Clinical and compound tool antagonists. Related to Figures 3 and 4.
Table S6. PDAC patients' characteristics (n=116). Related to Figure 7.
Complete output of gene set enrichment analysis conducted on stem and non-stem cell RNA-seq data-sets; gene signatures are include along with fold enrichment in stem cells, significance, and specific genes included in each signature.
(A) Establishment of three independent REM2-KPf/fC cell lines from end-stage REM2-KPf/fC mice for genomewide CRISPR-screen analysis. Stem cell content of freshly-dissociated REM2-KPf/fC tumors (A, left), and after puromycin selection in standard growth conditions (A, right).
(B and C) Volcano plots of guides enriched in 2D (B, tumor suppressors) and 3D (C, negative regulators of stem cells). Genes indicated on plots, p<0.005.
(D) Network propagation analysis integrating transcriptomic, epigenetic and functional analysis of stem cells. Genes enriched in stem cells by RNA-seq (stem/non-stem log2 fold-change >2) and depleted in 3D stem cell growth conditions (FDR<0.5) were used to seed the network (triangles), then analyzed for known and predicted protein-protein interactions and restricted to genes enriched in stem cells by RNA-seq (stem/non-stem log2 foldchange >2). Each node represents a single gene; node color is mapped to the RNA-seq fold change; stem cell enriched genes in red. Seven core programs were defined by groups of genes with high interconnectivity; each core program is annotated by Gene Ontology analysis (FDR<0.05).
(A and B) Sphere forming capacity of KPf/fC cells following shRNA knockdown. Selected genes involved in stem and developmental processes (A) or cell adhesion, cell motility, and matrix components (B).
(C-D) Immunofluorescence analysis of Celsr1 (C) and Celsr2 (D) in EpCAM+ stem (CD133+) and non-stem (CD133-) primary tumor cells isolated from KPf/fC mice. Three frames were analyzed per slide, and the frequency of Celsr1-high or Celsr2-high cells determined, scale=25um.
(E) KPf/fC cells were infected with shRNA against Pear1and protein knockdown efficiency determined five days post-transduction by western blot.
(F-H) Independent replicates for impact of shRNA inhibition of target genes on tumor growth in vivo. Celsr1 (F), Celsr2 (G), and Pear1 (H) were inhibited via shRNA delivery in KPf/fC cells, and impact on tumor growth assessed by tracking flank transplants in vivo, n=4 per condition.
(I) Pear1 was inhibited via shRNA in REM-KPf/fC cells in sphere culture and impact on Msi+ stem cell content assessed by FACS, n=3 per condition, p = 0.0629.
(J) Pear1 was inhibited via shRNA in KPf/fC cells and impact on apoptosis in sphere culture as marked by Annexin-V assessed by FACS, n=3 per condition.
(K) Heat map of relative RNA expression of cytokines and related receptors in KPf/fC stem and non-stem cells (left) and average RNA-seq TPM values in Msi2- and Msi2+ cells (right). Red, over-represented; blue, underrepresented; color denotes fold change from median values.
(L) Single cell RNA Sequencing maps of KPR172H/+C tumors. Tumor cells defined by expression of EpCAM (far left), Krt19 (left center), Cdh1 (right center), and Cdh2 (far right).
(M) Left, KPR172H/+C tumor single-cell sequencing map of cells expressing Msi2 within the EpCAM+ tumor cell fraction. Right, KPR172H/+C tumor single-cell sequencing map of cells expressing IL10Rβ, IL34, and CSF1R within the EpCAM+Msi2+ stem cell fraction.
(N-O) Independent replicates for impact of shRNA inhibition of target genes on tumor growth in vivo. IL10Rβ, (N) and CSF1R (O) were inhibited via shRNA delivery in KPf/fC cells, and impact on tumor growth assessed by tracking flank transplants in vivo, n=4 per condition.
(P) Cytokine receptors IL10Rβ and CSF1R were inhibited by shRNA delivery in KPf/fC cells and plated in sphere culture for one week. Increased apoptosis in KPf/fC cells with shIL10Rb (p<.05) and shCSF1R (trend). Frequency of apoptotic cells determined by Annexin-V staining and FACS analysis, n=3 per condition.
(Q) Representative FACS plots for stem content analysis in vitro. IL-10rβ and Csf1R were inhibited via shRNA delivery in KPf/fC cells, and impact on stem content (Msi2-GFP+ cells) in sphere culture assessed by FACS, n=3 per condition.
(R) ELISA based quantification (Quantikine, R&D Systems) of IL-10, IL-34, and CSF-1 in media (left) and KPf/fC cell lystate (right). Cytokines were quantified in fresh sphere culture media, KPf/fC stem and non-stem cell conditioned media, and KPf/fC epithelial cell lysate. Conditioned media was generated by culturing sorted CD133- or CD133+ KPf/fC cells in sphere media for 48 hours; media was filtered and assayed immediately. Cell lysate was collected in RIPA buffer and assayed at 2 mg/mL for ELISA. n=3 per condition.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01 by Student’s t-test or One-way ANOVA.
(A) Heat map of transcription factors in KPf/fC stem and non-stem identified as possible pancreatic cancer stem cell dependencies within the network map (see Figure 2E). Red, over-represented; blue, under-represented; color denotes fold change from median values.
(B) Distribution of RORγ consensus binding sites in genomic regions associated with H3K27ac. Down/Down, both gene expression and H3K27ac enriched in non-stem cells; Up/Up, both gene expression and H3K27ac enriched in stem cells.
(C) Biological replicates showing qPCR analysis of RORγ expression in primary KPf/fC stem and non-stem tumor cells isolated from REM2-KPf/fC mice.
(D) Immunofluorescence analysis of RORγ in primary KPf/fC EpCAM+ CD133+ and CD133- tumor cells. Three frames were analyzed per slide, and the frequency of RORγ -high cells determined.
(E) KPf/fC tumor single-cell sequencing map of cells expressing RORγ within the EpCAM+Msi2+ cell fraction (n=3 mice represented).
(F) RORγ expression within E-Cadherin- stromal cells in patient samples.
(G) IL1R1 was inhibited by CRISPR-mediated deletion in KPf/fC cells, and impact on Rorc expression assessed by qPCR. Two distinct guide RNAs (sgIL1r1–1 and sgIL1r1–2) were used to knockout IL1R1; expression was quantified by qPCR and is shown relative to control (non-targeting guide RNA), n=3 per condition.
(H) Knockdown efficiency of RORγ in KPf/fC cells infected with Rorc shRNA determined five days post-transduction. Relative expression in western blots quantified relative to tubulin loading control.
(I) Impact of shRNA-mediated RORγ inhibition on apoptosis and proliferation of in KPf/fC cells in 3D culture n=3.
(J) Independent replicate of shRNA Rorc impact on KPf/fC tumor propagation as assessed by tracking flank transplants in vivo, n=4 per condition.
(K-M) SE analysis of shRorc KPf/fC cells. KPf/fC cells were infected with shRorc, and used for H3K27ac ChIPseq and SE analysis, schematic (K). H3K27ac peaks were analyzed to assess SE overlap in shCtrl and shRorc samples (L). SEs lost in shRorc samples were crossed to stem-enriched and stem-unique SEs identified in primary Msi2-GFP+ KPf/fC tumors cells, and further restricted to SEs containing RORγ binding motifs (M). Majority of super-enhancer landscape remained unchanged with RORγ loss, and landscape changes that did occur were not enriched in SEs with RORγ binding sites. ChIP-seq analysis was conducted in two independent KPf/fC cell lines.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.
(A) Size of flank KPf/fC tumors in immunocompetent mice prior to enrollment into RORγ targeted therapy. Group 1, vehicle; group 2, SR2211; group 3, vehicle + gemcitabine; group 4, SR2211 + gemcitabine.
(B) Target engagement following acute RORγ inhibition in vivo. 9.5 wk tumor-bearing KPf/fC mice were treated with vehicle or SR2211 for two weeks (midpoint), after which tumors were isolated, fixed, and analyzed for target engagement of Hmga2 in epithelial cells by immunofluorescence. Representative images (left) and quantification (right) of Hmga2+ Keratin+ epithelial cells in vehicle or SR2211 treated tumors. Four frames were analyzed per mouse, n=2–4 mice per condition, Hmga2 (red), Keratin (green), scale=25um.
(C) Target engagement in endpoint tumors following continuous RORγ inhibition in vivo. 8 wk tumor-bearing KPf/fC mice were treated till endpoint with either vehicle or SR2211, after which tumors were isolated, fixed, and analyzed for target engagement of Hmga2 in epithelial cells by immunofluorescence. Representative images (left) and quantification (right) of Hmga2+ Keratin+ epithelial cells in vehicle or SR2211 treated tumors. Four frames were analyzed per mouse, n=2–4 mice per condition, Hmga2 (red), Keratin (green).
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA. Grubb’s test (p=0.1). was used to remove an outlier from the midpoint SR2211 treated group, scale=25um.
(A-B) Analysis of T cell subsets in KPf/fC tumors transplanted into wild-type or Rorc-knockout recipient mice (vehicle treated groups shown). Absolute cell numbers of the following populations were evaluated: (A) CD45+/CD3+/CD8+ or CD8+ T cells, (B) CD45+/CD3+/CD4+ or CD4+ T cells.
(C-L) FACS analysis of non-neoplastic cell populations in autochthonous tumors from KPf/fC mice treated with vehicle or SR2211 for 1 week. Schematic (C). Absolute cell numbers of the following populations were evaluated: CD45+ cells (D), CD11b+/F480+ cells (macrophage) (E), CD11b+/Gr-1+ cells (MDSC) (F), CD11c+ cells (dendritic) (G), CD45+/CD3+ T cells (H), CD3+/CD8+ T cells (I), CD3+/CD4+ T cells (J), CD4+/IL-17+ Th17 cells (K), CD31+ cells (endothelial) (L). (n=3 per condition).
(M) In vivo imaging of tumor vasculature of KPf/fC mice treated with vehicle or SR2211, vasculature is marked by in vivo delivery of anti-VE-Cadherin (magenta), scale=75um. Data represented as mean +/− S.E.M. * p<0.05 by Student’s t-test or One-way ANOVA.
(A-D) Gene ontology and gene set enrichment analysis of RNA-seq in human and mouse pancreatic cancer cells to identify common genes and pathways regulated by RORγ. Gene ontology analysis of KPf/fC RNA-seq showing genes downregulated with shRorc were enriched for cytokine-mediated signaling pathway GO term (A). Differentially expressed genes in KPf/fC within cytokine-mediated signaling pathway (B) were crossed with differentially expressed genes identified by RNA-seq analysis of human pancreatic cancer cells (FG) where RORC was knocked out using CRISPR. Gene set enrichment analysis of mouse and human RNA-seq shows common cytokine gene sets regulated by Rorc across species (D).
(E) Analysis of CRISPR guide depletion in stem cell conditions for SE associated genes expressed in stem or non-stem cells.
Data represented as mean +/− S.E.M. * p<0.05, ** p<0.01, *** p<0.001 by Student’s t-test or One-way ANOVA.