Abstract

The degree of metastatic disease varies widely among patients with cancer and affects clinical outcomes. However, the biological and functional differences that drive the extent of metastasis are poorly understood. We analyzed primary tumors and paired metastases using a multifluorescent lineage-labeled mouse model of pancreatic ductal adenocarcinoma (PDAC)—a tumor type in which most patients present with metastases. Genomic and transcriptomic analysis revealed an association between metastatic burden and gene amplification or transcriptional upregulation of MYC and its downstream targets. Functional experiments showed that MYC promotes metastasis by recruiting tumor-associated macrophages, leading to greater bloodstream intravasation. Consistent with these findings, metastatic progression in human PDAC was associated with activation of MYC signaling pathways and enrichment for MYC amplifications specifically in metastatic patients. Collectively, these results implicate MYC activity as a major determinant of metastatic burden in advanced PDAC.

Significance: Here, we investigate metastatic variation seen clinically in patients with PDAC and murine PDAC tumors and identify MYC as a major driver of this heterogeneity.

Introduction

Tumor heterogeneity, most commonly studied in a primary disease setting, is a critical driver of phenotypic diversity, culminating in metastatic, lethal cancers (1–5). In most cancers, prognosis and therapeutic decisions are defined by the presence or absence of metastasis. However, tumor heterogeneity is increasingly being questioned at the level of metastatic disease, with recent studies in several cancer types suggesting that metastasis is not a binary phenotype but rather a disease spectrum ranging from oligometastatic (limited) to polymetastatic (widespread) disesase (6–8). Heterogeneity in the manifestation of metastatic disease can guide decisions on the use of local–regional versus systemic therapies with emerging evidence of its importance in clinical outcome (9–11). Despite its clinical significance, the mechanisms that underlie this spectrum of metastatic states remain unclear and largely understudied.

Pancreatic ductal adenocarcinoma (PDAC) represents a disease entity well suited for the study of metastasis, as most PDACs present with metastatic disease that is associated with a dismal prognosis (12). Genomic studies have comprehensively cataloged core mutations responsible for primary tumor development in PDAC (e.g., KRAS, TRP53, CDKN2A, and SMAD4), paving the path for genomic investigations of metastatic disease and the identification of metastasis-promoting alterations. Indeed, recent sequencing studies and functional analysis in model systems have associated genomic amplification in mutant KRAS alleles with progression from the nonmetastatic (stage III) to metastatic (stage IV) disease state (13). However, genetic factors mediating metastasic heterogeneity in patients and, importantly, the downstream cellular mechanisms remain largely undefined (14–18). Furthermore, it is unclear whether metastasis-associated alterations perturb the tumor microenvironment, whose influence on metastatic behavior is well documented (19–31). Therefore, understanding the interplay between genetic alterations that influence metastatic behavior and the tumor biology that promotes it—via cell-autonomous and/or non–cell-autonomous mechanisms—is crucial for understanding metastasis as a distinct disease state and critical for the development of more effective treatments.

One barrier to understanding metastatic heterogeneity has been a paucity of model systems that capture this natural variation and allow for direct assessments of paired primary tumors and metastases in vivo. This has limited the ability to define factors intrinsic to primary tumors that influence the extent of metastatic spread. We previously developed an autochthonous model of PDAC—the KPCX model—that employs multiplexed fluorescence-based labeling to track the simultaneous development of multiple primary tumor cell lineages and follow them as they metastasize (32). Importantly, this technique facilitates confirmation of lineage relationships in vivo, such that primary tumor clones with substantial metastatic potential can be distinguished from those having poor metastatic potential. Here, we show that this system recapitulates the variation in metastatic burden found in human PDAC, and we use it to dissect molecular and cellular features contributing to metastatic heterogeneity.

ResultsMetastatic Burden Is Variable in Human and Murine PDAC

Although most patients with PDAC have metastases (principally liver and lung), the number of metastases is highly variable from patient to patient (33, 34). Importantly, data regarding metastases have largely been obtained at autopsy and thus are confounded by varying treatment histories and reseeding due to end-stage disease (4). Thus, we first sought to characterize the burden of metastases in treatment-naive patients. To this end, we performed a retrospective analysis of initial computed tomography (CT) scans from 55 patients newly diagnosed with metastatic (stage IV) PDAC (Fig. 1A). The total number of lesions in the lung and liver was counted by examining both coronal and sagittal planes for both organs and binned into groups of 10, revealing a wide distribution of metastatic burden (Fig. 1B). K-means clustering identified two metastatic subgroups: a MetLow subgroup (≤10 metastases, 25/55) and a MetHigh subgroup (>10 metastases, 30/55; Fig. 1B; Supplementary Fig. S1A). Primary tumor size, age, sex, and race were not correlated with differences in metastatic burden (Fig. 1C; Supplementary Fig. S1B). However, having a greater number of metastases was associated with worse overall survival (Fig. 1D). Thus, even among patients with stage IV PDAC, metastatic burden is variable and correlates with clinical outcome.

Figure 1.

Advanced pancreatic tumors exhibit intertumoral differences in their propensity for metastasis. A, CT imaging of human PDAC liver metastasis demonstrating heterogeneity in metastatic burden in stage IV disease. Arrowheads indicate solitary metastasis in the top image and selected metastases in the bottom. B, Density plot and histogram showing the distribution of total (liver and lung) metastases enumerated from CT scans of human stage IV PDAC at the time of diagnosis (n = 55). Values above each histogram bar represent the number of patients in each group. The vertical dotted line (red) represents the cutoff between MetLow tumors [≤10 metastases (mets)] and MetHigh tumors (>10 mets) determined by k-means clustering. C, Quantification of tumor area (based on tumor dimensions from largest cross-sectional plane on imaging) comparing MetLow and MetHigh cases from the cohort in B. D, Overall survival analysis of the cohort in B. E, Top, schematic view of the KPCXY model, showing multiple primary tumors distinguishable by color arising in the pancreas with matched metastases in the liver. Bottom, representative fluorescent stereomicroscopic images showing a YFP+ tumor adjoining a CFP+ tumor in the pancreas (left) and liver metastases derived from the CFP+ tumor in the same animal (right). F, Density plot and histogram showing the distribution of total (liver and lung) metastases enumerated at autopsy of KPCXY mice. Values above each histogram bar represent the number of tumors giving rise to the indicated number of metastases, based on color (n = 85 tumors from 30 KPCXY mice). The vertical dotted line (red) represents the cutoff between MetLow tumors (≤10 mets, n = 58) and MetHigh tumors (>10 mets, n = 28) determined by k-means clustering. G, Quantification of tumor area comparing MetLow and MetHigh tumors from the cohort in F. Statistical analysis by Student unpaired t test with P values indicated (ns, not significant). Box-and-whisker plots in C and G indicate mean and interquartile range. Scale bar for E, 1 mm.

We hypothesized that the differences in metastatic burden seen in human PDAC may also be present in autochthonous murine models. To test this, we used the KPCXY model—in which Cre-mediated recombination triggers expression of mutant KrasG12D and deletion of one allele of Trp53 in the pancreatic epithelium along with YFP and confetti (X) lineage tracers (Fig. 1E; Methods)—to measure metastatic heterogeneity in a cohort of tumor-bearing mice. By exploiting the multicolor features of the KPCX model, we previously showed that these mice harbor (on average) two to five independent primary tumor clones; importantly, the clonal marking of tumors with different fluorophores makes it possible to infer the lineages of primary tumors with different metastatic potential (32). In our earlier work with this model, we noted that in most tumor-bearing animals, even those with multiple primary tumors, liver and lung metastases were driven by a single tumor clone (Fig. 1E; Supplementary Fig. S1C). This suggested that tumor cell–intrinsic factors strongly influence the metastatic behavior of a tumor, even within a single animal.

To quantify differences in metastatic burden, we examined a panel of mice with at least two uniquely labeled fluorescent tumors in which most metastases could be attributed to a specific tumor on the basis of color (Fig. 1E; Supplementary Fig. S1C). A total of 85 primary tumors from 30 mice were examined, and gross metastases to the liver and lung arising from each tumor were then quantified by stereomicroscopy (Methods). Murine PDACs exhibited a wide distribution of metastatic burden, with a pattern resembling that of the human disease (Fig. 1F). Similarly, K-means clustering grouped murine samples into a low-metastasis subgroup (≤10 metastases, 58/85) and a high-metastasis subgroup (>10 metastases, 27/85), which we similarly refer to as MetLow and MetHigh, respectively (Fig. 1F; Supplementary Fig. S1D). As with the human disease, neither primary tumor size nor tumor cell proliferation correlated with metastatic burden (Fig. 1G; Supplementary Fig. S1E). Thus, the KPCXY model recapitulates the intertumoral metastatic heterogeneity seen in human PDAC and provides a unique experimental model for comparing highly metastatic and poorly metastatic tumor clones.

Individual Tumor Lineages in KPCXY Mice Correspond to Clones with Distinct Somatic Copy-Number Profiles

Although primary KPCXY tumors were easily distinguishable based on the expression of a distinct fluorophore, each tumor could have arisen via the clonal expansion of a single cell or through fusion of multiple tumors that happened to share the same color. Somatic copy-number alterations (SCNA) have been shown to provide an unambiguous picture of genomic heterogeneity and lineage relationships between primary tumors and matching metastases in human disease (35). Consequently, we performed copy-number analysis via genome sequencing on a set of 20 primary tumors, including multiregional sampling on a subset of the tumors in which sufficient tissue was available (nine tumors with two to four regions sampled per tumor; Fig. 2A; Supplementary Table S1). Tumors bearing different colors exhibited unique DNA copy-number profiles, indicating that they arose independently (Fig. 2B; Supplementary Fig. S2A; ref. 36). By contrast, multiregional sampling of monochromatic tumors revealed shared copy-number alterations, indicating that all subregions within a given tumor (defined by color) shared a common ancestral lineage (Fig. 2C; Supplementary Fig. S2B). In addition, subregion-specific alterations were also observed, suggesting that subclonal heterogeneity is also present in each tumor (Fig. 2C; Supplementary Fig. S2B). These results suggest that the monochromatic tumors observed in KPCXY mice are clonal in origin and continue to undergo subclonal evolution during tumor progression.

Figure 2.

SCNA analysis confirms fluorescence-based lineage relationships and reveals genetic heterogeneity in paired primary pancreatic tumors and liver metastases. A, Schematic representation of KPCXY pancreatic tumor and matching liver metastases with multiregion sampling for copy-number sequence analysis. B, Representative genome-wide copy-number profiles of MetHigh (CFP+ fluorescence) and MetLow (YFP+ fluorescence) tumors from mouse 832 (m832) as depicted in Fig. 1E. Gray shading denotes alterations that are unique to the MetHigh (CFP+) tumor. The y-axis illustrates normalized read count values (low ratio), which are directly proportional to genome copy number at a given chromosomal location. The copy-number profiles are centered around a mean of 1 with gains and deletions called for segments with values higher and lower than the mean, respectively (Methods). C, Representative genome-wide copy-number profiles of three subsampled tissue regions of the MetHigh (CFP+) primary tumor from m832. Gray shading denotes alterations that are found heterogeneously from multiregion sequencing of the primary tumor. D, Genome-wide heat map with hierarchical clustering based on copy-number alterations of matched primary and metastatic samples profiled from m832. E, Representative genome-wide copy-number profiles of fluorescently matched primary and metastatic tissue from two profiled mice (m832, left; m836, right), illustrating the shared clonal genetic lineage. F, Zoomed-in chromosomal views of copy-number alterations with distinguishing breakpoint patterns supporting shared genetic lineage. Panels are ordered as in E.

To ascertain the lineage relationships between primary tumors and metastases, we compared DNA copy-number profiles between liver metastases and primary tumors within a given mouse. This revealed that primary tumors and metastases of the same color shared common DNA copy-number profiles across the dataset, confirming on a genetic basis the fluorescence-based lineage relationships (Fig. 2D–F; Supplementary Fig. S2C). As most lung metastases were microscopic and difficult to isolate by dissection, they were not included in the molecular analysis. Together, these results indicate that the lineage history of metastases can be inferred by color and genomic analysis, allowing primary tumors with high versus low metastatic potential to be unambiguously classified.

Genomic and Transcriptional Analyses Identify Myc as a Potential Driver of Metastatic Phenotypes

We next sought to examine the molecular differences that distinguish primary tumors with high versus low metastatic potential. We began by examining large-scale (mega-base level as well as chromosome-wide) SCNAs in 20 MetHigh and MetLow primary tumor samples. This analysis revealed largely similar genome-wide copy-number patterns between MetHigh and MetLow primary tumors, with key PDAC-associated genes, such as loss of heterozygosity at Cdkn2a/b and Trp53 as well as chromosomal gain of Kras occurring at similar frequencies (Supplementary Fig. S3A). Thus, KPCXY tumors exhibit frequent copy-number alterations in canonical PDAC genes, but these alterations do not account for the variation in metastatic behavior between MetHigh and MetLow tumors.

We next asked whether other factors (genomic and/or transcriptional) may be acting to enhance metastasis in the MetHigh group. Focal amplifications in driver oncogenes—Cdk6 and Yap in breast cancer and mutant Kras in PDAC—have been linked to the acquisition of metastatic competence (13, 14, 37, 38). Consistent with prior studies, we observed focal amplicons at genomic regions encoding Cdk6, Yap, and Kras in our tumors (Fig. 3A; Supplementary Fig. S3A; refs. 13, 14, 37–39). However, in contrast to these amplifications, which occurred at equal frequencies in MetHigh and MetLow tumors, focal high-amplitude amplifications in Myc were found in 42.8% (3/7) of MetHigh tumors compared with 7.6% (1/13) of MetLow tumors (Fig. 3B). Thus, Myc amplifications are enriched in MetHigh tumors. In all cases, these amplifications were maintained in paired metastases (Supplementary Fig. S3B). In addition, RNA sequencing (RNA-seq) analysis demonstrated significantly higher levels of Myc transcripts in MetHigh tumors and metastases compared with MetLow tumors (Fig. 3C); overall, Myc was the third-most significantly upregulated gene in MetHigh tumors compared with MetLow tumors (Fig. 3D). Gene set enrichment analysis (GSEA) of the differentially expressed genes between MetHigh and MetLow tumors identified MYC and E2F signatures as highly enriched, along with other signatures that have been implicated in PDAC metastasis, including unfolded protein response, oxidative phosphorylation, and hypoxia (Fig. 3E; Supplementary Table S2; refs. 31, 40, 41). Moreover, MetaCore transcription factor enrichment analysis identified MYC as the transcription factor most significantly associated with genes overexpressed in MetHigh tumors (Fig. 3F), and Ingenuity Pathway Analysis placed Myc at the center of the interactome generated by these differentially expressed genes (Supplementary Fig. S4). Collectively, these results demonstrate a strong association between a tumor's metastatic behavior and the abundance and/or activity of MYC at the genomic and transcriptional levels.

Figure 3.

The MetHigh phenotype is associated with focal, high-amplitude Myc amplifications and elevated expression. A, Schematic representation of focal amplifications identified in profiled primary tumors. Vertical gray line denotes the location of amplicon and likely driver gene. B, Zoomed-in schematic representation of three identified Myc amplicons in MetHigh tumors illustrating the focal and high-amplitude nature of the event (left). Each event (amplicon) is illustrated by a different colored segment line. The shared amplified region between the different amplicons is denoted by the chromosomal cytoband top of panel and illustrated in a UCSC Genome Browser view (right) with RefSeq genes, including Myc, illustrated. C, Box-and-whisker plot showing Myc mRNA levels in MetHigh tumors (n = 7) and paired metastases (n = 34) compared with MetLow tumors (n = 13). FPKM, fragments per kilobase of exon per million. D, Volcano plot illustrating genes meeting cutoffs for differential expression [log fold change (logFC) >1, Padj. < 0.05] between MetHigh and MetLow tumors (n = 20 tumors used in the comparison). Genes upregulated in MetHigh tumors are highlighted in green, and genes upregulated in MetLow tumors are highlighted in red. E, Top 10 hallmark gene sets identified as enriched in MetHigh tumors compared with MetLow tumors using all differentially expressed genes (DEG; Padj. < 0.05). F, Top five transcription factor (TF) binding sites enriched in DEGs in MetHigh tumors compared with MetLow tumors (Padj. < 0.05) identified by Metacore prediction software. G, Heat map showing unsupervised clustering of DEGs (logFC >1, Padj. < 0.05) between MetHigh and MetLow tumors (n = 20) and their association with PDAC transcriptional subtypes previously reported by Collisson and colleagues (42), Moffitt and colleagues (15), and Bailey and colleagues (16). ADEX, aberrantly differentiated endocrine exocrine; QM-PDA, quasi-mesenchymal-pancreatic ductal adenocarcinoma. H, Kaplan–Meier analysis showing overall survival of patients with PDAC in the TCGA cohort stratified into those with a MetHigh signature (red line) versus those with a MetLow signature (green line). Signature based on DEGs with absolute logFC >0.58 and Padj. < 0.05 (736 up- and 1,036 downregulated genes). Statistical analysis in C was performed by Wilcoxon test (*, P = 3.9 × 10−4; **, P = 5.3 × 10−5). Box and whiskers represent median mRNA expression and interquartile range. Statistical analysis in H was performed by log-rank test.

Human PDAC can be grouped into two main transcriptomic subtypes—a well-differentiated classical/exocrine-like/progenitor (classical) subtype and a poorly differentiated squamous/quasi-mesenchymal/basal (basal-like) subtype (15, 16, 18, 42). We found that MetHigh tumors were associated with basal-like PDACs, in line with their more aggressive behavior (Fig. 3G). Likewise, applying murine MetLow and MetHigh signatures (see Methods) to human data from The Cancer Genome Atlas (TCGA) predicted a worse survival—indicative of disease recurrence—for patients with a MetHigh signature (Fig. 3H). These data indicate that murine MetHigh tumors correspond to the more aggressive subtypes of human PDAC.

A Panel of Cell Lines that Preserve the MetLow and MetHigh Phenotypes

To understand the mechanisms underlying these different metastatic properties, we generated a panel of cell lines from six MetHigh tumors and five MetLow tumors. Consistent with the parental in vivo tumors, Myc gene expression and MYC protein levels were higher in the MetHigh lines compared with the MetLow lines (Fig. 4A and B). SCNA analysis in these cells lines found that they retained most of the genomic alterations found in the matched primary samples, including Myc amplifications (Supplementary Fig. S5A and S5B). Furthermore, Myc amplifications were not found in any of the cell lines whose tumors were originally characterized as non–Myc amplified, indicating that in vitro culture does not select for this specific copy-number alteration. Importantly, elevations in MYC mRNA and protein were observed in both the Myc-amplified and nonamplified MetHigh lines, suggesting that elevated MYC expression is a stable phenotype of these cells in culture.

Figure 4.

MYC regulates metastasis by enhancing tumor cell intravasation. A, Bar graph showing Myc mRNA levels in cell lines derived from MetHigh and MetLow tumors, normalized to Gapdh (n = 6 MetHigh and n = 5 MetLow cell lines). B, Western blot showing corresponding MYC protein levels in cell lines derived from MetHigh and MetLow tumors shown in A. C, Representative fluorescent images of primary tumors and associated liver and lung metastases following orthotopic transplantation of the cell lines in A and B into NOD.SCID mice. The bar graph shows the total number of metastases (liver and lung) counted following orthotopic transplantation of five MetLow cell lines or five MetHigh cell lines (pooled data from n = 49 mice in total). D, Representative fluorescent images of primary tumors, liver metastases, and lung metastases following orthotopic transplantation of MetLow cell lines that were stably transduced with either a Myc_OE or an empty vector (EV) construct. The bar graph shows the total number of metastases (liver and lung) counted following orthotopic transplantations of Myc_OE or EV cells. Data were pooled from four independent MetLow lines transduced with either the Myc_OE or EV construct transplanted into 12 NOD.SCID (for the Myc_OE cells) or 10 NOD.SCID mice (for the EV cells). E, Quantification of CTCs in arterial blood derived from the orthotopic tumors depicted in C (n = 27 mice examined) and D (n = 12 mice examined). F, Representative fluorescent images of lung metastases following tail vein injection of cell lines derived from the MetLow and MetHigh primary tumor clones. The bar graph shows the total number of lung metastases counted following tail vein injection of five MetLow cell lines or five MetHigh cell lines (pooled data from n = 36 mice in total). Statistical analysis by Student unpaired t test with significance indicated (*, P = 0.0152; **, P = 0.013; ***, P = 0.0008; ns, not significant). Error bars indicate SEM (C–F). Scale bar, 1 mm (C, D, and F).

To investigate the metastatic properties of the MetHigh and MetLow lines in vivo, we performed orthotopic implantation of five MetHigh and five MetLow lines into the pancreas of NOD.SCID mice and examined distant organs for evidence of metastasis. Although the weights of MetHigh and MetLow tumors were not significantly different (Supplementary Fig. S6A), MetHigh tumors gave rise to 28-fold more liver and lung metastases compared with MetLow tumors (Fig. 4C). Consistent with the cell line expression differences, the orthotopic MetHigh tumors expressed higher levels of Myc compared with MetLow tumors (Supplementary Fig. S6B). To further confirm that differences in Myc expression were sufficient to drive the metastatic phenotype, we introduced a Myc overexpression (Myc_OE) construct into four MetLow lines and generated orthotopic tumors (Supplementary Fig. S6C). Myc overexpression led to a dramatic (22-fold) increase in liver and lung metastases (Fig. 4D). Thus, cell lines derived from spontaneously generated MetHigh and MetLow tumors retain their metastatic phenotypes upon implantation.

MYC Promotes Tumor Cell Intravasation through the Recruitment of Tumor-Associated Macrophages

To form distant metastases, cancer cells must navigate a series of events collectively referred to as the “metastatic cascade.” These events include (i) intravasation into the bloodstream or lymphatics, (ii) survival in the circulation, (iii) extravasation from the vessel, and (iv) growth and survival at the distant site (43). To determine the step(s) at which MYC was exerting its prometastatic effects, we began by measuring the number of circulating tumor cells (CTC) in orthotopically implanted MetHigh and MetLow tumors and in MetLow tumors engineered to overexpress Myc. Remarkably, CTCs arising from MetHigh and Myc_OE tumors were 38-fold and 17-fold more abundant than those arising from MetLow tumors (Fig. 4E), far greater than the approximately twofold increase in tumor weight resulting from Myc overexpression (Supplementary Fig. S6D). Next, we performed a tail vein metastasis assay, which bypasses the invasion step by introducing tumor cells directly into the bloodstream, and measured lung metastases. Surprisingly, in contrast to the orthotopic tumor experiment, there was no difference in the number of metastases between MetHigh and MetLow lines (Fig. 4F). Moreover, Myc overexpression had no effect on tumor cell survival in the circulation (Supplementary Fig. S6E–S6G). Taken together, these data suggest that MetHigh tumors achieve a higher metastatic rate principally by promoting cancer cell invasion into the circulation, which can be driven by increased Myc expression.

Beyond activation of tumor cell intrinsic programs, MYC can also affect tumor phenotypes by altering the tumor immune microenvironment (TiME; refs. 44–46). Thus, we sought to determine if differences in MYC levels between MetHigh and MetLow tumors were associated with distinct TiMEs. To this end, we examined the immune composition of parental primary tumors by staining for markers of immune cells previously implicated in metastasis of PDAC and other cancers. Although MetHigh and MetLow tumors had a similar degree of neutrophil infiltration, MetHigh tumors had lower numbers of CD3+ T cells but were highly enriched for F4/80+ macrophages (Fig. 5A). Thus, compared with MetLow tumors, the TiME of MetHigh tumors contains an increased number of tumor-associated macrophages (TAM).