Gene Expression Profile of Metastatic Human Pancreatic Cancer Cells Depends on the Organ Microenvironment

At the time of diagnosis, human neoplasms are heterogeneous and consist of multiple subpopulations of cells with differing biological properties that include metastatic potential (1, 2). Cells with different metastatic properties have been isolated from primary neoplasms and from metastatic lesions, showing that not all tumor cells have the capacity to complete the arduous metastatic process (1). Card analysis using marker chromosomes showed that metastases originate from a single progenitor cell (3). Isolation of cells with different metastatic potential from heterogeneous neoplasms has been accomplished by in vivo selection or by random cloning of cultured parental tumor cells (1). Tumor cells with different metastatic capabilities have been shown to differ in expression of proteins, such as collagenases, E-cadherin, interleukin-8, basic fibroblast growth factor, and many others (4). In most cases, however, these differential expressions were most evident in tumor cells growing at specific orthotopic sites.

Recent advances in the development of microarrays provide a method to assess genome-wide modifications in gene expression that can differentiate between metastatic and nonmetastatic cells isolated from a single parental neoplasm. Many biological investigations of solid tumor cells often use cell lines growing in vitro as monolayer cultures. Although easy to accomplish, monolayer cultures are not subjected to any cross talk (e.g., paracrine signaling pathways associated with growth in vivo; ref. 5). Clinical observations of cancer patients and studies in rodent models of cancers have concluded that certain tumors tend to metastasize to certain organs (6). The concept that metastasis results only when certain tumor cells interact with a specific organ microenvironment was originally proposed in Paget's venerable “seed and soil” hypothesis (7). Indeed, spontaneous metastasis is produced by tumors growing at orthotopic sites, whereas the same tumor cells implanted into ectopic sites fail to produce metastasis (8).

The purpose of this study was to determine whether expression of genes associated with metastatic potential of tumor cells depends on the organ microenvironment. We have previously described an orthotopic model for metastasis of human pancreatic cancer used to select variant cell lines with increasing metastatic capacity (9). We used gene expression profiles generated from microarray analysis to compare the phenotypes of three pancreatic cell lines with different metastatic potentials growing in vitro as monolayer cultures or in vivo at ectopic (s.c.) and orthotopic (pancreatic) sites. We found that the differential gene expression profile in metastatic pancreatic cancer cells depends on growth in a biologically relevant orthotopic organ microenvironment.

Pancreatic cancer cell lines and culture conditions. The selection of the human pancreatic cell lines with different metastatic potentials was previously described (9). The low metastatic COLO 375 FG human pancreatic cancer cell line was injected into the spleen of nude mice. Cells isolated from liver lesions designated as L3.3 produced liver lesions at a higher incidence than the original COLO 375 FG cells did. Next, the L3.3 cells were injected into the pancreas of nude mice. Spontaneous liver metastases were harvested and cells were established in culture. The selection cycle was repeated thrice to yield a cell line designated L3.6pl (pancreas-to-liver metastasis; ref. 9).

The variant cell lines (FG, L3.3, and L3.6pl) were maintained in Eagle's MEM supplemented with 10% fetal bovine serum, sodium pyruvate, nonessential amino acids, l-glutamine, a 2-fold vitamin solution (Life Technologies, Inc., Grand Island, NY), and a penicillin-streptomycin mixture (Flow Laboratories, Rockville, MD). Adherent monolayer cultures were maintained on plastic and incubated at 37°C in a mixture of 5% CO₂ and 95% air. The cultures were free of Mycoplasma and the following pathogenic murine viruses: reovirus type 3, pneumonia virus, K virus, Theiler's encephalitis virus, Sendai virus, minute virus, mouse adenovirus, mouse hepatitis virus, lymphocytic choriomeningitis virus, ectromelia virus, and lactate dehydrogenase virus (assayed by M.A. Bioproducts, Walkersville, MD). The cultures were maintained for no longer than 12 weeks after recovery from frozen stocks.

Animals and production of tumors. Male athymic BALB/c nude mice (NCI-nu) were obtained from the Animal Production Area of the National Cancer Institute Frederick Cancer Research and Development Center (Frederick, MD). The mice were maintained under specific pathogen-free conditions in facilities approved by the American Association for Accreditation of Laboratory Animal Care and in accordance with current regulations and standards of the U.S. Department of Agriculture, Department of Health and Human Services, and NIH. The mice were used in accordance with institutional guidelines when they were 8 to 12 weeks old.

Tumor cell injection techniques. For in vivo injection, cells were harvested from subconfluent cultures by a brief exposure to 0.25% trypsin and 0.02% EDTA. Trypsinization was stopped with medium containing 10% fetal bovine serum and then the cells were washed once in serum-free medium and resuspended in HBSS. Only single-cell suspensions of >90% viability (assessed by trypan blue exclusion assay) were used for injection.

Orthotopic injection. The mice were anesthetized with pentobarbital sodium administered by i.m. injection. A small left abdominal flank incision was made and the spleen exteriorized. Tumor cells (1 × 10⁶ per 50 μL HBSS) were injected subcapsularly in a region of the pancreas just beneath the spleen. We used a 30-gauge needle with a 1-mL disposable syringe. A successful subcapsular intrapancreatic injection of tumor cells was identified by the appearance of a fluid bleb without i.p. leakage. To prevent leakage, a cotton swab was held for over the site of injection for 1 min. One layer of the abdominal wound was closed with wound clips (Auto-clip, Clay Adams, Parsippany, NJ). The animals tolerated the surgical procedure well and no anesthesia-related deaths occurred.

S.c. injections. Tumor cells (1 × 10⁶ per 50 μL HBSS) were injected into the subcutis of the lateral flank. Tumors were harvested when they reached 7 to 9 mm in diameter. Half of each tumor was snap frozen in liquid nitrogen for mRNA extraction; the other half was fixed in formalin and embedded in paraffin.

RNA sample preparation. Total RNA was extracted from cultured cells or tumor tissues by using Trizol reagent (Invitrogen, Carlsbad, CA). The extracted RNA was passed through an RNeasy spin column (Qiagen, Valencia, CA) with on-column DNase I treatment to remove any contaminating genomic DNA according to the manufacturer's protocol.

Affymetrix GeneChip hybridization. Human Genome U133 Plus 2.0 GeneChip arrays (Affymetrix, Santa Clara, CA) were used for microarray hybridizations. This GeneChip carries 54,675 probe sets. For microarray hybridization, we followed the protocol described in the manufacturer's eukaryotic one-cycle target preparation protocol. In short, 10 μg of total RNA were used to prepare antisense biotinylated RNA and single-stranded cDNA was synthesized using a T7-Oligo(dT) promoter primer followed by RNase H–facilitated second-strand cDNA synthesis, which was purified and served as a template in the subsequent in vitro transcription. The in vitro transcription reaction was carried out in the presence of T7 RNA polymerase and a biotinylated nucleotide analogue/ribonucleotide mix for cRNA. The biotinylated cRNA targets were then cleaned up and fragmented. The fragmented cRNA was used for hybridization to the U133 Plus 2.0 chip at 42°C for 16 h. The chips were washed and stained with Affymetrix GeneChip Fluid, then scanned and visualized with a GeneArray Scanner (Hewlett-Packard, Palo Alto, CA).

Statistical analysis. Perfect match and mismatch features on the scanned microarray images were quantified using the Affymetrix Microarray Suite 5.0 software (MAS 5.0). After quantification, quality control was done using the Affymetrix quality control metrics. Two arrays that failed quality control (i.e., percent present calls ≤40% when the median percent present call was 51.8%) were repeated with fresh extractions of RNA. The qualified data sets were then analyzed with DNA-Chip Analyzer software (dChip³

We filtered the genes to do hierarchical clustering using the most reliably measured genes. More precisely, we retained genes that were called present in at least 20% of genes and had an expression value >500 in at least half of the samples.

We did a separate two-way ANOVA for each probe set on the Affymetrix array using the linear model Y_ijk = μ + α_i + β_j + γ_ij + E_ijk (Eqn. 1), where Y_ijk is the observed expression of the gene in replicate k of cell line i under growth condition j; μ is the overall mean expression of the gene; α_i represents the effect due to the cell line (i = FG, L3.3, or L3.6pl); β_j represents the effect due to the growth condition (j = culture, s.c., or orthotopic); γ_ij represents an interaction between cell line and growth condition; and E_ijk ∼ N (0, σ) is the residual errors, which are assumed to be independent and identically distributed. We did an F test for each probe set to determine whether the model fit the observed data better than the null model Y_ijk = μ + E_ijk and computed the corresponding (unadjusted) P values. To account for multiple testing, we modeled the collection of all 54,675 P values using a β-uniform mixture model (12). We then used the β-uniform mixture model to estimate the false discovery rate (13).

Pattern identification. All genes that were found to be significant by ANOVA were clustered, using unsupervised hierarchical clustering in R, to identify patterns of gene expression across the sample set.

Select genes that are up-regulated in vitro, in vivo, and in metastasis. To select the 50 most highly up-regulated genes in each of these three categories, we focused on genes with significant ANOVA results. These genes were filtered to identify those expressed at a higher level in orthotopic tumors than in culture (in vivo set), higher in culture than orthotopically (in vitro set), or higher in L3.6pl than in FG (metastasis set). The filtered genes were ranked by the P values for the significance of the term in the model related to growth condition (in vivo and in vitro sets) or to cell lines (metastasis set).

Gene Ontology analysis. Patterns that were represented in the top 50 lists were analyzed to identify functional categories using the Database for Annotation, Visualization and Integrated Discovery (DAVID).⁵

We did 18 microarray experiments using human U133 Plus 2.0 GeneChip microarrays, which contained 54,675 probe sets. Samples came from several mouse xenograft models of human pancreatic cancer, in which three different cell lines were grown in culture or implanted s.c. or orthotopically. The experiment followed a complete 3 × 3 factorial design with two replicates of each combination. To determine differentially expressed genes between tumor models (cell culture, s.c., or orthotopic) and between cell lines based on their metastatic potential (L3.6pl, high; L3.3, low; FG, none), we used both hierarchical clustering and gene-by-gene ANOVA.

Hierarchical clustering. A total of 8,686 genes passed the filter (see Materials and Methods) and we did hierarchical clustering of the samples using the Pearson product-moment correlation coefficient to define distances and Ward's linkage rule (Fig. 1). All replicate pairs clustered as nearest neighbors, increasing our confidence in the reproducibility of the results. The results showed that the six experiments done on cell cultures clustered as a block, confirming that gene expression in vitro differs from the expression pattern of the same cells growing in vivo. At the next level, it seemed that clustering was driven more by the cell line (equivalently, by metastatic potential) than by the growth conditions of the tumor model. Reproducibility of the clusters was tested using a bootstrap cluster test (14); at the very least, the four main branches were robust to perturbations of the data (Supplementary Fig. S1).

ANOVA. To examine the two main effects (tumor model type and cell line or metastatic potential) along with their interactions, we did a separate two-way ANOVA for each of the 54,675 probe sets on the U133 Plus 2.0 GeneChip. β-Uniform mixture analysis indicated that substantial numbers of genes change expression between cell lines and in response to growth conditions (Supplementary Fig. S2). Overall, we found 5,495 statistically significant probe sets when false discovery rate <0.001, which corresponded to a cutoff of P < 0.00046 on the unadjusted P values (see Materials and Methods). A complete list of the significant probe sets is available on our website.⁶

Using a global F test to compare the complete linear model (Eqn. 1) to the null model avoids additional multiple testing, but it requires an additional step to determine which factors drive the differential expression of individual genes. TSF2o address this issue, we clustered all 5,495 significant probe sets based on their standardized expression profiles in the 18 samples using the Pearson product-moment correlation coefficient to define distances and average linkage (Fig. 2). We cut the dendrogram to produce 30 different clusters; using 30 for the number of clusters was an arbitrary choice but it seems to be large enough to identify all the distinct patterns in Fig. 2.

To relate the cluster patterns of Fig. 2 to the factors in the two-way balanced design of the experiment, we computed the average expression pattern of the genes in each cluster and displayed each resulting pattern in a 3 × 3 layout (Fig. 3). Differences that are attributable to the local environment show up in these graphs as horizontal bands. For example, patterns 1, 2, 7, 13, and 14 seem to represent in vivo patterns. These genes were expressed at significantly lower levels in all three cell lines when they were growing in vitro and were expressed at higher levels when growing in vivo. In the same way, patterns 20, 21, 28, and 29 seem to represent an in vitro pattern.

We examined the top 50 genes (among the 5,495 significant ones) ranked by the significance of the “growth condition” term in the ANOVA model (see Materials and Methods) for the contrast between growth in vitro and in vivo. Most of the 50 genes up-regulated in cells growing in culture came from patterns 21 and 28. On the other hand, genes expressed in cells growing in vivo came primarily from patterns 1, 2, and 7 (Tables 1 and 2). To better understand the functional differences among those patterns, we placed the genes from each pattern into DAVID, which identified statistically significant functional categories using Gene Ontology (Table 3). Multiple categories were statistically significant. The patterns up-regulated in vitro were enriched for cell cycle (pattern 21), metabolism (patterns 21 and 28), and translation (pattern 28). The patterns up-regulated in vivo were enriched for regulation of transcription (pattern 1), antigen processing (pattern 2), and protein modification (pattern 7). These data suggest that gene expression of cell lines growing under in vitro conditions was optimized for cell proliferation, whereas gene expression of the same cell lines growing in vivo was affected by many other factors related to the microenvironment.

Differences between the cell lines, which may reflect their metastatic potential, show up as vertical bands in these graphs. For example, patterns 4, 5, and 6 seem to include genes that are highly expressed in L3.3, and patterns 22 and 26 are expressed at high levels almost exclusively in FG. Patterns 10 and 16 seem to include genes that are highly expressed in L3.6pl growing both in vitro and in vivo, and thus may represent a metastatic pattern. In the same way, patterns 7 and 30 represent a lower metastatic pattern because they are expressed at higher levels in FG and L3.3 than in L3.6pl.

We also examined the top 50 genes up-regulated in metastasis ranked by the significance of the “cell line” term in the ANOVA model (see Materials and Methods). Most of these genes came from patterns 5, 10, and 18 (Table 4). We identified multiple functions for these patterns using DAVID (Table 5), including morphogenesis (pattern 5), signal transduction (pattern 5), protein phosphorylation (pattern 10), glycolysis (pattern 10), Wnt receptor signaling pathway (pattern 18), and secretory pathway (pattern 18). These results suggest that ligand production and receptor activation are more prominent in the metastatic cell lines.

Many of the patterns illustrate complex effects of both the local environment and the metastatic potential. Some of these effects seem to be additive, as in pattern 7. In that case, expression is high in FG and in L3.3 cells, but only when they are growing in vivo. The expression is low for all three of the cell lines in vitro and stays low for the highly metastatic cell line L3.6pl even when it is grown in vivo. In addition, highly metastatic L3.6pl cells growing in the pancreas expressed significantly higher levels of 226 genes than did the L3.3 or FG variant cells (in patterns 10, 11, 13, 17, and 19). In addition, the L3.6pl cells growing in the s.c. space expressed significantly higher levels of 98 genes than the L3.3 or FG variant cells (see Materials and Methods). We used DAVID to examine gene patterns that were differentially expressed when the highly metastatic L3.6pl cells growing in the pancreas (orthotopic site) were compared with the same cells growing in the subcutis (ectopic site; Table 6). The patterns of genes expressed in orthotopic tumors were enriched for protein amino acid dephosphorylation, ion transport, ATP binding, and actin cytoskeleton organization and biogenesis. By contrast, the patterns of genes expressed ectopically were enriched for tissue development, negative regulation of cell proliferation, regulation of cell differentiation, and negative regulation of apoptosis. (Lists of the top 50 up-regulated genes in orthotopic and s.c. L3.6pl tumors are provided in Supplementary Tables S2 and S3).

We determined that growth in a specific organ microenvironment is a prerequisite for identifying a differential gene expression pattern unique to metastatic human pancreatic cancer cells. Three variant lines of a single human pancreatic cancer with different metastatic potentials were grown in vitro, in the subcutis, and in the pancreas of nude mice. The pattern of gene expression was similar in the three lines growing in culture but not in vivo (Fig. 1). Cell cultures formed a distinct branch of the dendrogram, but there was a clear expression signature that differentiated between the three cell lines once they were growing in vivo. Specifically, the pattern of gene expression associated with the metastatic phenotype of the highly metastatic cell line L3.6pl was most different from that of the L3.3 or FG cell lines when the tumor cells were growing in the pancreas. (This assertion follows from the height of the main splits in the dendrogram of Fig. 1; it is also apparent from a principal component analysis contained in Supplementary Fig. S3.) We previously showed that implantation of tumor cells at orthotopic or ectopic sites in nude mice produces tumors with different metastatic potentials (1, 4, 8). The orthotopic tumors are invasive and metastatic, whereas the ectopic s.c. tumors are not, suggesting that different microenvironments may differentially influence the expression of metastasis-related genes (15–17). A significant difference in gene expression patterns was identified between orthotopic tumors and ectopic tumors of L3.6pl (Table 6). The orthotopic tumors expressed 226 unique genes and the ectopic tumors expressed 98 genes. Our data showing that the growth of the variant cell lines in the subcutis did not yield similar results indicate that the orthotopic microenvironment significantly influences specific gene expression in pancreatic cancer cells.

By jointly analyzing the data from the three cell lines in multiple environments, we have been able to identify expression patterns that are specific to the highly metastatic cell line L3.6pl when it is grown orthotopically. Specifically, the genes in patterns 10, 11, 13, 17, and 19 have their highest expression level in orthotopically grown L3.6pl cells, and the genes in patterns 6, 22, 23, 24, 29, and 30 have their lowest expression under the same conditions. The 50 most highly expressed genes of the metastatic cell line, regardless of the microenvironment, are listed in Table 4. The dominant pattern was pattern 10, which has its highest expression in orthotopically grown L3.6pl cells. The Gene Ontology analysis identified functional categories that were up-regulated in the metastatic cell line (Table 5). Many of the genes with catalytic (enzymatic) activities are involved in several phosphorylated signaling pathways.

Some of the genes have previously been identified as metastasis-related genes (18, 19). Tropomysin-related kinase B (Table 4) was overexpressed in the patients with pancreatic cancer and its expression correlated with perineural invasion, a positive retroperitoneal margin, and shorter latency to the development of liver metastasis (18). Axl receptor tyrosine kinase (Table 4) is highly expressed in a variety of tumors (19–22) and was reported to be expressed at 10-fold higher levels in a peritoneal metastatic nodule than in other normal and malignant tissues (19). Axl signaling has been shown to affect neovascularization in vitro and angiogenesis in vivo (23).

Many investigators have used mRNA microarrays to study the effects of in vitro and in vivo microenvironments on gene expression profiles in various cancers (24–28). Sandberg and Ernberg (29) reported on a meta-analysis of gene expression profiles from three laboratories comprising 60 cell lines and 311 tissues (135 normal tissues and 176 tumor tissues). They concluded that the genes involved in cell cycle progression, macromolecule processing and turnover, and energy metabolism were up-regulated in cell lines, whereas cell adhesion molecules and membrane signaling proteins were down-regulated. Our data reveal the same tendency for differences between in vitro and in vivo microenvironments. Differences in gene expression by tumor cells growing in vitro and in vivo were also reported by Camphausen et al. (30), who compared two different glioblastoma cell lines, U87 and U251, grown in cell culture, s.c., and intracerebrally. The gene expressions in the cells grown in the three environments differed for both cell lines. Unlike our data, the gene expression patterns of the U87 and U251 cells in that study were more similar in the orthotopic tumors than in ectopic s.c. tumors, and the greatest difference between the two cell lines was found in the in vitro cultures (30). A follow-up study that examined the effects of radiation on U87 and U251 cells reported few commonly affected genes in cultured cells and greater changes that were detected after orthotopic radiation (31).

The pathogenesis of metastasis selects for tumor cells that succeed in invasion, embolization, survival in the circulation, arrest in the capillary beds, and then extravasation and growth in the organ parenchyma (32). Thus, the outcome of cancer metastasis depends on multiple interactions of metastatic cells with a specific organ microenvironment (32, 33). Recent in situ hybridization and immunohistochemical staining studies have shown that the expression of various genes and proteins varies between different regions of a single tumor and between tumors implanted in different sites of nude mice (16, 17). Our data show that the organ microenvironment clearly influences the pattern of gene expression profiles of tumor cells that have different metastatic potentials. These data suggest that both the nature of tumor cells and the host microenvironment contribute to a variety of gene expressions necessary for the development of cancer metastasis.

Grant support: Cancer Center Support Core Grant CA16672 and Specialized Program of Research Excellence in Prostate Cancer Grant CA902701 from the National Cancer Institute, NIH.

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.

We thank Karen F. Phillips, ELS, for critical editorial review, and Lola López for expert preparation of this manuscript.