A Visually Apparent and Quantifiable CT Imaging Feature Identifies Biophysical Subtypes of Pancreatic Ductal Adenocarcinoma

This article is featured in Highlights of This Issue, p. 5785

Pancreatic ductal adenocarcinoma (PDAC) is associated with early distant metastasis (DM) and resistance to chemotherapy and radiation. However, this disease has substantial heterogeneity in biophysical features and in outcomes for patients (1–5). The identification of a biomarker that distinguishes patients with more aggressive disease from those with less aggressive disease would enable more rational therapeutic choices.

One such biomarker in PDAC may involve the stroma, elements of which have been linked with suppression of metastatic spread of the cancer cells (6–8). Notably, the stroma changes in response to therapy, and the nature of those changes may be associated with more aggressive cancer biology (9). The native and dynamic roles of the stroma remain areas of intense scientific investigation. Clarifying the role of stroma is thought to have important clinical implications that could affect patient survival.

We previously showed that the stroma may represent a physical barrier to drug delivery in PDAC (3). We further demonstrated that mass transport properties derived from CT scans correlated with the delivery of, response to, and outcome after gemcitabine-based therapies (2, 3). These findings suggested that properties derived from CT scans could serve as biophysical markers of PDAC, although they could not provide insight a priori into the natural history of the disease or underlying physical and biological mechanisms related to patient outcomes.

Here, we hypothesized that morphologic characteristics of human PDAC tumors would indicate the aggressiveness of the disease, as with other cancers (10–17). We developed a method to quantify the morphology of PDAC, classifying patients based on this imaging-based feature and correlating that classification with pathologic and clinical findings. We investigated the biophysical mechanisms related to our observations by using a mathematical theory of macroscopic PDAC tumor growth (18) that considered the stroma as a determinant of cancer cell proliferation and migration rates.

All studies were conducted in accordance with U.S. Common Rule and with Institutional Review Board approval (PA14–0646), including waiver of informed consent. We reviewed patients in Supplementary Tables S1–S4. Consort flow diagrams are shown in Supplementary Fig. S1. Only patients with pancreatic protocol CT scans prior to treatment were included.

We conducted studies using data from a group of patients who enrolled on a previously described intraoperative gemcitabine infusion clinical trial (3). Briefly, patients with untreated, resectable pancreatic cancer were eligible (NCT01276613, Supplementary Table S1).

We included patients who underwent resection of treatment-naïve, localized PDAC between 2001 and 2013 (Supplementary Table S2).

Patients who enrolled on two phase II trials (19) of neoadjuvant therapy for resectable PDAC are described in Supplementary Table S3.

We included patients who presented between 2012 and 2014 with treatment-naïve stage IV PDAC, with metastatic disease either confirmed pathologically or radiographically (Supplementary Table S4).

The pancreatic protocol CT scan is a diagnostic test for pancreatic cancer, where iodine-based contrast is injected intravenously at a fixed rate (20). The test usually consists of a precontrast phase, an arterial phase (35–40 seconds after starting contrast infusion) and a portal-venous phase (65–70 seconds after starting contrast infusion). We describe the use of each phase of the pancreatic protocol CT for our measurements below.

The preoperative CT images (2.5-mm slice thickness) were imported into Velocity AI (Varian Medical Systems) for image registration of the arterial phase and portal venous phase of the pancreatic protocol CT scan for each patient. After performing deformable registration of the different phases of the scan, the interface between the PDAC tumor and the surrounding pancreas parenchyma was characterized by volumetrically contouring both the tumor and pancreas at the border, while avoiding the pancreatic duct, the surrounding fat space of the pancreas, and metal artifacts. The mean value of the Hounsfield unit (HU) distribution within each contour was compared, providing a difference in HU called “delta” (Fig. 1A and B).

The details of our original analysis method were described previously (2, 3). We analyzed the same registered CT scans as described for the delta measurement.

The histologic sections of pancreatic cancer tissue stained with hematoxylin and eosin (H&E; test set: 12 cases, see Supplementary Table S1 for patients’ characteristics; validation set: 17 cases, part of cohort in Supplementary Table S2) were scanned using a Vectra slide scanner (PerkinElmer). We used identical exposure times for all slides. The whole slide was scanned at 4×, then 100 high-power images within the tumor were scanned at 20×. All scanned images were reviewed (Vectra Review 1.0.5, PerkinElmer). Images were excluded if they were out of focus or contained artifacts. To establish the spectral library of H&E, a representative digitized image was imported to Nuance (PerkinElmer). The spectral profiles of each of these single stains were used to unmix images that contain H&E to quantify each stain. Then, we established two different image analysis methods using inForm (PerkinElmer).

inForm software was used to segment the nuclei of cells in the slide and subsequently identify proportions of stroma cells and lymphocytes. The software provides pattern recognition machine learning algorithms as well as cell phenotyping algorithms, which are trained to segment subcellular compartments and identify specific types of cells, respectively. A pathologist (DC) provided annotation of 1,950 cells [cancer cells: 563, stroma cells: 151, lymphocytes: 310, and normal pancreas cells (islet cells: 309, acinar cells: 382, ductal cells: 235)]. We established two different cell type identification training methods. First, we used four cell type classifiers (four types, 1: cancer cells; 2: stroma cells; 3: lymphocytes; 4: normal pancreas cells). This analysis system showed high accuracy of identification of stroma cells (97.2%) and lymphocytes (98.7%). However, the accurate identification of cancer cells and normal pancreas cells (e.g., acinar cells, ductal cells, islet cells) was limited, especially when normal pancreas cells were located around cancer cells. Cancer cells were correctly identified at 81.1%, and normal pancreas cells were identified correctly at 57.4%. For comparison, we applied three cell type classifiers (three types, 1: cancer/normal pancreas cells, 2: stroma cells, 3: lymphocytes). Both methods (three or four cell type classifiers) identified the same total number of cells on each of the slides that we analyzed (range, 22,202 to 90,976 cells/slide) while maintaining a high degree of accuracy in identifying stroma cells and lymphocytes. Therefore, we used the algorithm with three classifiers to count stroma cells and lymphocytes in each slide and used the total number of cells to calculate the proportions of stroma cells and lymphocytes in each tumor.

Specific images from each slide, which contained only cancer cells, stroma cells, and/or lymphocytes, were selected by a pathologist (D. Li). Then, the cell phenotyping algorithm with three classifiers, as described above, was applied to identify nuclei of pancreatic cancer cells (range, 1828 to 37966 cells/slide). The nucleus axis ratio was calculated as minor axis divided by major axis.

In addition to validation by a pathologist, we performed quality control measures (hematoxylin intensity of the nucleus, nucleus area, and major/minor axes of nucleus) to ensure that the cells were being identified correctly and that the measurements of the cell features were in line with expectation. ImageJ 1.48v (NIH, Bethesda, MD) was used to measure the major and minor axes of cells.

MDA-PATC50, 66, 69, 102, and 118 cell lines were derived from our pancreatic xenograft program (21–23). The Institutional Review Board and the Institutional Animal Care and Use Committee at The University of Texas MD Anderson Cancer Center approved these studies. Details of cell isolation, fingerprinting, culture, and Western blotting are in Supplementary Methods.

Exome DNA was captured using the Agilent 50M capture kit and tumor/normal sample pairs were sequenced with 100 base pair paired-end reads on an Illumina HiSeq 4000 sequencer at a mean depth of 150× for normal samples and 300× for tumor samples. Reads were aligned with the Burrows-Wheeler Alignment software (BWA) to the reference human genome version hg19 (24). We processed aligned reads from these tumor/normal pairs using MuTect to call point mutations and small insertions and deletions (25).

We use a multicomponent mixture modeling framework to simulate solid tumor progression (10, 12, 26–30). Tumor and host regions are described as a saturated medium, comprised of a mixture of viable and dead cells. Volume fractions, φ_V and φ_D, describe the relative amounts of these components, respectively. The governing equations consist of mass and momentum balance equations for each component. The mass balance equations are:

where u_i is the velocity, J_i is a flux, and S_i is a source term that accounts for proliferation and death. The source terms are:

where Λ_P is the mitosis rate, λ_N is the death rate associated with cell necrosis, H is the Heaviside function, σ is the concentration of cell substrates, which models the combined effects of oxygen, glucose, and growth factors, σ_N is the concentration of substrates needed for cell viability, and λ_L is the clearance rate of dead cells. Appealing to a generalized Fick's law (27), the flux |{{\bf{J}}_i} = - M{\phi _i}\nabla \mu $|⁠, where |\mu = f'( {{\phi _T}} ) - {\varepsilon ^2}{\nabla ^2}{\phi _T}$|⁠, with |{\phi _T} = {\phi _V} + {\phi _D}$|⁠, is a chemical potential. We take |f( {{\phi _T}} ) = {\frac{1}{4}}\phi _T^2{( {1 - {\phi _T}} )^2}$| to be a double well potential, which reflects the tendency of the tumor cells to adhere to one another. The parameter ϵ models the range of cell–cell interactions and introduces a finite thickness (proportional to ϵ) of the tumor–host interface.

We assume that the mixture is tightly packed so that viable and dead cells move with the same velocity, which is taken to be a generalized Darcy law (27):

where p is the solid pressure, Λ_M is a chemotaxis migration coefficient, σ is the concentration of life-sustaining cell substrates, which models the combined effects of oxygen, glucose, and growth factors, and |{F_a} = {\frac{\gamma }{\varepsilon }}\mu \nabla {\phi _T}$| is a cell–cell adhesion velocity where γ is a measure of cell–cell adhesion (27).

Assuming that the sum of the tumor and host volume fractions is constant, an equation for the pressure can be derived by summing the equations for the volume fractions (27):

Finally, the cell substrates are assumed to diffuse throughout the tissue and be supplied by an underlying vasculature:

where we have assumed that substrates diffuse sufficiently rapidly so as to enable a quasi-steady approximation, D is the diffusion coefficient, λ is the uptake rate, S_n is the nutrient source, p_H is the rate of substrate delivery, σ* is the effective level of substrates delivered and φ_H is the host volume fraction. Details of the mathematical simulations, and parameter values are in the Supplementary Methods (10, 12, 14, 26–30). Here, we only discuss the estimation of the proliferation and migration parameters responsible for the morphologic stability of the tumors.

JMP Pro 11.2.0 (SAS) was used to perform all statistical analyses. Mann–Whitney (Wilcoxon) test was used to compare distributions between groups. Survival curves were constructed using the Kaplan–Meier method. The Cox proportional-hazards model was used for univariate and multivariate survival analyses. Overall survival was defined from time of diagnosis to time of death or last follow-up for all cohorts. Distant metastasis-free survival was defined from time of diagnosis to time of first distant metastasis or last follow-up for all cohorts. To build multivariate Cox proportional hazards models, we started with the potential covariates which demonstrated a P<0.2 in the univariate models, and backward model selection technique was used. C-index was also calculated for the multivariate models. We considered a P value less than 0.05 as significant.

We previously demonstrated (10, 11, 13, 28), using a biophysical theory of tumor morphologic stability (12), that the features of the tumor at the microscopic scale were described by the competition of mechanisms that oppose infiltrative growth (termed “relaxation mechanisms,” e.g., cell proliferation represented here by the mitosis rate Λ_P) and those that promote infiltrative growth (e.g., cell migration, quantified by the chemotaxis migration coefficient Λ_M). One prediction of this model is that the tumor interface may be used to differentiate tumors that have highly proliferative cancer cells from those that are more indolent. We tested this idea by characterizing the tumor interface of PDAC on CT scans.

We developed our method to characterize the tumor interface by using pretreatment CT images from 12 patients who underwent upfront resection of PDAC (Supplementary Table S1; ref. 3). We quantified the “delta” measurement as the volumetric difference in the mean of the Hounsfield unit (HU) distributions for the observed PDAC tumor at the border and the surrounding pancreatic parenchyma (Fig. 1A and B).

Next, two radiologists (E.P. Tamm and P. Bhosale) ranked the conspicuity of the PDAC interface on a 1-to-5 scale (Supplementary Fig. S2) for these 12 patients (Supplementary Table S1). We found that a cutoff for the delta measurement of 40 HU distinguished conspicuous from inconspicuous lesions (Supplementary Fig. S2). We validated this result with presurgical CT scans of 101 patients who underwent upfront surgery for resectable PDAC (Supplementary Table S2). We performed receiver-operating characteristic (ROC) area under the curve (AUC) analysis over a range of 35 to 45 HU, giving an AUC that ranged from 0.81 to 0.95, and found that the cutoff of 40 performed best in both arterial and portal venous phases. Supplementary Figure S3A shows the ROC curves for the different delta cutoffs in the arterial phase that best distinguished conspicuity of the tumors read by both radiologists. In the portal venous phase, the delta cutoff of 40 also had the best AUC (0.88 for ET's score and 0.83 for PB's score). Thus, we used the cutoff of 40 to define high-delta and low-delta groups. There was high concordance for the delta classification in 40 patients from the cohort in Supplementary Table S2 for whom pancreatic protocol CT scans were performed both outside and inside our institution (kappa = 0.8; Supplementary Figs. S3E and S4).

The same 12 patients who provided the CT scans for developing the delta measurement (Supplementary Table S1; ref. 3), also provided tissue samples that were scored by pathologists (H. Wang and D. Chatterjee) for the amount of stroma in the resected tumors (Fig. 1C). We found a statistically significant association between the amount of stroma and the 40 HU cutoff for the delta measurement (Fig. 1D). We validated this association in another 33 patients who underwent upfront resection (Supplementary Fig. S5).

We then used a machine-learning algorithm (Supplementary Fig. S6) to identify different cellular populations on slides that were stained with H&E. For the 12 patients who underwent upfront resection (Supplementary Table S1), the quantitative pathology analysis demonstrated that high-delta tumors had a lower percentage of stroma cells (26% ± 10%) than did the low-delta tumors (40% ± 5%, Fig. 1E). The percentage of lymphocytes was also significantly different between low- and high-delta tumors (Fig. 1F). IHC staining showed high-delta PDAC had a higher proportion of cells that coexpressed markers for T regulatory cells (Tregs) than low-delta PDAC (Fig. 1G). We found the same associations for proportions of stroma cells and lymphocytes in a second group of patients who underwent upfront surgery for localized PDAC (Supplementary Figs. S7 and S8; Supplementary Table S2).

We reviewed 64 patients with PDAC who underwent genomic profiling (31). Fifty-two (81%) of these patients had a high-delta pancreatic tumor. In the entire group, the most common somatic mutations were in KRAS (81%), TP53 (55%), SMAD4 (8%), and PIK3CA (4%). Notably, KRAS and TP53 mutations frequently occurred together (48%), followed by KRAS and SMAD4 (11%) and KRAS and PIK3CA (4%). Patients with high-delta tumors were more likely to have cooccurring somatic mutations (60%) in these common genes than patients with low-delta tumors (17%, Fisher exact test, P = 0.01; Fig. 2A). We identified 15 patient-derived xenografts (PDX) from patients who underwent upfront resection of localized PDAC (23), from which we performed whole-exome sequencing and mutation profiling. Considering missense and splice site mutations, of the nine high-delta PDAC samples, seven of them had both KRAS and TP53 mutations. Of the six low-delta samples, only one of them had mutations in both KRAS and TP53. Harboring mutations in both KRAS and TP53 was associated with high-delta PDAC (P = 0.04), where such tumors had 13.6 times higher risk of being high delta (95% CI: 0.86–934; Fig. 2B).

We isolated and characterized five PDAC cell lines from patients who underwent upfront surgery (ref. 23; Fig. 3A and B). We applied the delta measurement to the preoperative CT scans from these patients. We took pictures of the cell lines at low passage in monolayer culture (less than P10), and used a membrane axis measurement to quantify cell morphology (defined as shortest axis measurement divided by longest axis measurement). We found that cancer cell lines derived from low-delta tumors had significantly higher membrane axis ratios (i.e., were rounder) than those derived from high-delta tumors (more elongated; Fig. 3C). We further observed higher expression of mesenchymal marker proteins in the cell lines derived from high-delta tumors than in those from low-delta tumors (Fig. 3D).

Because nuclear morphology is directly associated with cellular shape (32), we further quantified cancer cell nuclear morphology by using the machine-learning algorithm on the pathology samples represented in Fig. 1. The nuclei of the cancer cells in the high delta tumors were more elongated than the cancer cells in low-delta tumors, whose nuclei were rounder (Fig. 3E). This was also validated in another group of patients who underwent upfront surgery (Supplementary Fig. S7C).

We previously found that higher volumetric measurements of enhancement derived from CT scans (volumetric area under the enhancement curve, VAUC) were associated with longer survival outcomes (3). Here, analysis of 101 patients who underwent upfront resection for localized PDAC (Supplementary Table S2) and 106 patients who underwent neoadjuvant chemoradiation (Supplementary Table S3) demonstrated that patients with low-delta tumors had significantly higher VAUC (0.61 ± 0.21 and 0.65 ± 0.22, respectively) than did patients with high-delta tumors (0.40 ± 0.20 and 0.37 ± 0.15, respectively, P < 0.0001 for both groups; Fig. 3F).

We next measured tumors from the patients who underwent neoadjuvant chemoradiation (Supplementary Table S3). Analysis of the distributions of the tumor axis ratios demonstrated that the high-delta tumors were more likely to be round, whereas the low-delta tumors were less likely to be round (Fig. 3G).

We analyzed 84 consecutive patients who presented with DM at diagnosis and received upfront chemotherapy (Supplementary Table S4). Notably, 65 (77%) of those 84 patients had a high delta in the primary pancreatic tumor (Fig. 4A, left). Survival results indicated prolonged overall survival (OS) for patients with low-delta tumors (Supplementary Fig. S9A and Supplementary Table S5).

We extended this imaging-based classification to the CT scans of 101 patients with early-stage, resectable PDAC (Supplementary Table S2). In contrast to the patients with stage IV disease, a much higher percentage of patients in this early-stage PDAC group had low-delta tumors (54%, Fig. 4A, right). Moreover, patients with high-delta tumors had worse distant metastasis-free survival (DMFS) at 2 and 5 years (46% and 42%, respectively) than did patients with low-delta tumors (66 and 44%, respectively, log-rank P = 0.07, Fig. 4B, left). The earlier time to DM translated into poorer OS rates at 2 and 5 years for patients with high-delta tumors (44% and 8%, respectively) than for patients with low-delta tumors (61% and 31%, log rank P = 0.0013, Fig. 4C, left).

The imaging delta classification was independently significantly associated with OS in multiple regression analysis, and in comparison with CA19–9, the delta classification improved the model's predictive power (c-index: 0.679 in model B vs. 0.667 in model C; Supplementary Table S6). Exploratory analysis confined to patients with T3 disease, stratified by nodal disease, revealed a significant improvement in the c-index when the delta classification was included in the multivariate model (Supplementary Table S7). Similar analysis for CA19–9 revealed a lower c-index in the N0 group compared with the multivariate model using the delta classification (Supplementary Table S7). Sensitivity analyses of the delta measurement demonstrated significant differences in the outcomes of the patients with delta measurements in the range between 35 and 45 HU (Supplementary Fig. S10). The visual classification (i.e., scoring by radiologists; Supplementary Fig. S2) was also prognostic (Supplementary Fig. S11).

We performed a retrospective–prospective validation of the results in Fig. 4B and C (left). The same correlations for actuarial outcomes were seen in a second dataset of 106 patients with resectable pancreatic cancer who received neoadjuvant gemcitabine-based chemoradiation (19). The delta classification was highly associated with DMFS and OS (Fig. 4B and C, right). Furthermore, after neoadjuvant therapy, patients with high-delta tumors (20/58) were more likely than those with low-delta tumors (8/48) to have progressive disease by Response Evaluation Criteria in Solid Tumors (RECIST) v1.1 (P = 0.035), primarily due to distant progression. There was no difference in local progression between those with high-delta tumors (8/58), compared with low delta (4/48, P = 0.35). The delta classification was independently highly associated with OS in multivariate analysis (Supplementary Tables S8 and S9) and again effectuated improved predictive power when compared with CA19–9.

Because we previously found VAUC to correlate with patient outcomes and identified a cutoff of 0.6 as a way to stratify patients in terms of outcome, we explored a combined classification based on the delta score and VAUC. This combined classification seemed to subdivide patients further into different prognostic groups (Supplementary Fig. S12).

Motivated by recent studies suggesting that the stroma acts to restrict the metastatic spread of PDAC (7, 8), we hypothesized that our mathematical theory of tumor growth could be scaled up from the microscopic level to the macroscopic tissue level by considering the stroma as a global relaxation mechanism that strongly influenced gross tumor morphology (described in Methods).

We introduce a key, dimensionless model parameter (Λ = mitosis rate Λ_P/chemotaxis migration coefficient Λ_M) that describes the stability of the tumor interface and represents a ratio of biological factors that influence proliferation and migration. We performed parametric analyses of this key model parameter to ascertain its effect on gross morphology of the tumor (Fig. 5). When the proliferation rate is slower than the migration rate (i.e., a low Λ < 1), then the model predicts that the cancer cell clusters will intermingle with the stroma, resulting in an indistinct interface for the tumor (Fig. 5A, right; Λ = 0.2) and generating what have been described as “low-mode” interface instabilities (12). Conversely, when the proliferation rate is higher than the migration rate (i.e., a high Λ greater than the critical value of unity), the model predicts that the tumor cells will grow with a distinct interface, as proliferation overcomes any attempts of isolated migrating cells to separate from and leave the main tumor bulk (Fig. 5B, right; Λ = 1.5). Representative movies of these simulations are provided in Supplementary Materials. These simulations show direct analogy to the macroscopic features of PDAC tumors on patient CT scans (Fig. 5C and D).

Our measurements of the pathologic specimens produce surrogates for the determinants of the stability parameter (Λ) from our model of tumor growth. To illustrate, we used the average stromal content in high-delta (25% ± 10% and 30% ± 18%; Fig. 1D; Supplementary Fig. S5) and low-delta tumors (46% ± 17% and 59% ± 17%; Fig. 1D; Supplementary Fig. S5) and calculated a low-delta to high-delta stromal ratio of approximately 2. This ratio in stromal content would be inversely correlated with cell proliferation rates according to the model assumptions, which provides an estimate for the ratio of proliferation rates Λ_P between high and low-delta tumors. Similarly, the ratio of membrane axis ratios of cancer cells derived from low- and high-delta tumors is approximately 3, indicating that cancer cells from high-delta tumors are more elongated and have a mesenchymal-like phenotype (Fig. 3A–E). This number thus provides a surrogate for the ratio in migration rates Λ_M between high- and low-delta tumors.

Therefore, the estimated ratio of stability parameters (i.e., Λ_high-delta/Λ_low-delta) is approximately 2 of 3, which suggests that these parameters have similar values across the two phenotypes. Accordingly, the range of the stability parameter used in the computer simulations was narrowly centered around the critical value of unity (Fig. 5 used Λ = 1.5 and Λ = 0.2). Notably, the computer simulations correctly predict the observed, different phenotypes around the critical value of unity.

The current clinical practice for PDAC largely involves a one-size-fits-all approach. Although genetic subtypes have been identified (33, 34), clinicians continue to lack reproducible, biologically meaningful, and clinically relevant biomarkers to stratify patients into prognostic groups at initial diagnosis. We have focused on the pancreatic protocol CT scan, which is a ubiquitous, standard-of-care test, to identify such a biomarker. We developed a classification of PDAC based on the conspicuity of the tumors on CT. This resulted in pathologic, genetic, cellular, and clinical correlations that consistently indicate that high-delta tumors are more aggressive than low-delta tumors. Furthermore, our data support a mathematical theory that describes tumor morphologic patterns as a balance between tumor cell proliferation and migration, providing mechanistic insight into the observed phenomena.

The advantages of this imaging-based classification of PDAC are that it can be done at the point of care before treatment, and applies to both localized and metastatic PDAC. Molecular classifications of PDAC can be challenging because they require direct pathologic sampling of the tumor (33, 34). Our approach complements molecular characterizations by providing a noninvasive metric of tumor behavior. We anticipate that our imaging-defined high-delta group (median survival time 13 months in resected PDAC) largely comprises the “squamous” subtype of PDAC (median survival time 13.3 months; ref. 35), whereas the low-delta group (median survival time 24 months) is a mix of molecular subtypes with better prognosis (median survival times 23.7−30 months). Our results support a molecular basis of the imaging-defined subtypes, as high-delta tumors were more likely than low-delta tumors to harbor multiple poor prognosis mutations in combination with KRAS, such as TP53, SMAD4, and PIK3CA (Fig. 2).

Notably, mutant KRAS (36), TP53 (37, 38), SMAD4 (39), and PIK3CA (40, 41) have been linked to immunosuppressive states through immunomodulation and recruitment of T regulatory cells. Here, we found that high-delta tumors had a lower degree of stromal cell infiltrate than did low delta tumors (Fig. 1). Our findings also showed differences in the proportions of lymphocytes in high and low delta tumors. Notably, our data indicated that high-delta tumors have higher proportions of Tregs in the microenvironment compared with low delta tumors. Tregs suppress immune activity and have been associated with poor outcomes in multiple cancer types, including PDAC (42). This finding agrees with our observation that patients with high-delta tumors have poorer outcomes and may lead to rational methods to stratify and select patients for stromal-modifying therapies and immunotherapy. Furthermore, the findings in this work support the idea that the multiscale properties of cancer are intimately related to the underlying biology, which is central to the concept of oncophysics (4, 43).

In addition to differing immunologic and stromal properties, we noted variability in the enhancement properties of these tumors based on the imaging. This work builds directly on our previous imaging studies of PDAC (2, 3). While this study supports that low-delta PDAC has a better prognosis than high delta PDAC, our previous work on radiographic changes has shown heterogeneity in responses with regard to baseline conspicuity on CT (44). A classification approach of the VAUC and delta measurements at baseline, combined with the tumor's radiographic response to therapy may further stratify patients (Supplementary Fig. S12). The heterogeneity of responses indicates differing sensitivities within each delta subtype, and in vitro drug testing may help clarify the molecular determinants related to the observed heterogeneity (45). Because low- and high- delta PDAC exhibit heterogeneity in responses to neoadjuvant therapy in the primary tumor, the difference in prognosis in these groups is likely due to differences in development of distant metastasis. Indeed, our genetic, cellular, and pathologic data indicate that high-delta tumors have more aggressive biology leading to a higher propensity for distant metastasis, compared with low delta PDAC.

We acknowledge that this study represents a mixture of retrospective and prospective datasets, which were not specifically designed to investigate the imaging characteristics of PDAC. Furthermore, the method of CT acquisition evolved during the study period. However, the key aspect that enables reproducible measurements of the delta is our approach to push the contrast bolus during the pancreatic protocol CT scan (Methods). This has not significantly changed during the study period, and all the scans were reviewed with our radiologists to ensure a uniform quality of scans. A future direction for this work is prospective validation of our findings in a phase II study for borderline resectable PDAC (NCT01560949), SWOG S1505 (NCT02562716), and Alliance A021501 (NCT02839343). These trials may provide further evidence that the delta classification may help predict benefit from neoadjuvant therapy.

In summary, we identified a physical feature on pre-therapy, standard-of-care CT images—the presence or absence of a distinct border between PDAC tumors and the surrounding parenchyma—which associates with cellular, pathologic, biological, physical, and clinical characteristics. Each of these correlations indicates that patients with high-delta tumors generally have more aggressive disease and poorer outcomes than do patients with low-delta tumors (Table 1). These results can be used in basic and translational studies to identify how cancer/host interactions may be regulated and manipulated to alter the natural course of disease. Furthermore, our imaging-defined subtypes can be incorporated into ongoing and future clinical trials—as well as reviews of past “failed” therapies—to identify methods to personalize and improve the care of patients with this lethal disease.

K.A. Reiss reports receiving commercial research support from Bristol-Myers Squibb, Clovis Oncology, Eli Lilly Oncology, and TesaroBio. R.A. Wolff receives royalty payments from McGraw-Hill for acting as co-editor of the MD Anderson Manual of Medical Oncology. C.H. Crane reports receiving speakers bureau honoraria from Celgene. E.P. Tamm reports receiving commercial research grants from General Electric. No potential conflicts of interest were disclosed by the other authors.

Conception and design: E.J. Koay, Y. Lee, V. Cristini, P.C. Park, N. Garg, A. Maitra, P. Bhosale, J.B. Fleming

Development of methodology: E.J. Koay, Y. Lee, V. Cristini, J.S. Lowengrub, Y. Kang, H. Yan, N. Garg, M. Javle, H. Wang, P. Bhosale, E.P. Tamm, J.B. Fleming

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): E.J. Koay, Y. Lee, V. Cristini, J.S. Lowengrub, M. Almahariq, A.M. Amer, S. Yu, M.V. Rios Perez, D. Li, N. Garg, K. Reiss, A. Chauhan, N. Nikzad, R.A. Wolff, G.R. Varadhachary, R.T. Shroff, P. Das, C.M. Taniguchi, C.H. Crane, H. Wang, P. Bhosale, E.P. Tamm, J.B. Fleming

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): E.J. Koay, Y. Lee, V. Cristini, J.S. Lowengrub, F.A. San Lucas, B.P. Hobbs, R. Ye, M. Almahariq, D. Chatterjee, H. Yan, M.V. Rios Perez, D. Li, N. Nikzad, C.M. Taniguchi, H. Wang, P. Bhosale, E.P. Tamm, J.B. Fleming

Writing, review, and/or revision of the manuscript: E.J. Koay, Y. Lee, V. Cristini, J.S. Lowengrub, B.P. Hobbs, D. Elganainy, M. Almahariq, N. Garg, K. Reiss, M.M. Zaid, R.A. Wolff, M. Javle, G.R. Varadhachary, R.T. Shroff, P. Das, J.E. Lee, A. Maitra, C.M. Taniguchi, M.P. Kim, C.H. Crane, H. Wang, P. Bhosale, E.P. Tamm, J.B. Fleming

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): E.J. Koay, P.C. Park, N. Garg, N. Nikzad, M.H.G. Katz, E.P. Tamm, J.B. Fleming

Study supervision: E.J. Koay, V. Cristini, M. Javle, J.E. Lee, J.B. Fleming

Other (development of mathematical model and hypothesis of mechanistic correlation of lesion morphology and morphologic stability to cellular phenotypes and microenvironmental/stromal properties): V. Cristini

Other (research consultation): M. Ferrari

We thank Chris Wogan for her thoughtful editing of the manuscript. We also thank Mark Hurd, Amalia Gonzalez, Alexander Ondari, Christopher Bristow, and Stephanie Kerps for their research support. We gratefully acknowledge partial support from the Andrew Sabin Family Fellowship, Center for Radiation Oncology Research, the Sheikh Ahmed Center for Pancreatic Cancer Research, institutional funds from The University of Texas MD Anderson Cancer Center, equipment support by GE Healthcare and the Center of Advanced Biomedical Imaging, Philips Healthcare, Project Purple, and Cancer Center Support (Core) Grant CA016672 from the National Cancer Institute to MD Anderson. E. Koay was also supported by NIH grants U54CA210181-01, U54CA143837, and U01CA196403, the Pancreatic Cancer Action Network (14-20-25-KOAY), and the Radiological Society of North America (RSD1429). J. Fleming was supported by the Lustgarten Foundation (989161) and Viragh Family Foundation. This work was also supported by NIH T32CA009599. V. Cristini was supported by NSF DMS- 1716737, NIH 1U01CA196403, 1U01CA213759, 1R01CA226537, 1R01CA222007, the Rochelle and Max Levit Chair in the Neurosciences, and the University of Texas System STAR Award. J. Lowengrub was supported by NSF DMS-1714973 and NIH grants 1U54CA217378-01A1, P50GM76516, and P30CA062203. This research was partly performed in the Flow Cytometry & Cellular Imaging Facility.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.