Genetic profiling of IGF1R signaling-mediated mammary gland tumor formation reveals cellular processes associated with decreased tumor latency
Previously, we chronically treated p53R270H/+WAPCre mice with insulin and insulin analogues to study the effect of chronic IGF1R stimulation on tumorigenesis (Fig. 1a). We found that MG tumor latency time in X10 and IGF1 (compounds with a high affinity towards the IGF1R) treatment groups was significantly decreased (Fig. 1b). Tumors in mice treated with glargine (an insulin analogue that induces only a mild IGF1R activation in the MG) also developed earlier compared to vehicle- and insulin-treated animals, but this trend was not significant [7]. To further investigate this, we evaluated receptor gene expression levels. We found that gene expression levels of the insulin receptor (both the A and B isoform; ira, irb, respectively) and insulin-like growth factor 1 receptor (igf1r) are significantly upregulated in pre-neoplastic MG tissue of old (50 weeks on average) p53R270H/+WAPCre mice compared to healthy 8-week-old MGs (Fig. 1c). Also, ira and igf1r levels were slightly upregulated in MG tumors. Interestingly, we found a significant and sevenfold downregulation of the B-isoform of the insulin receptor in MG tumor tissue compared to normal MG tissue. This effect was not treatment specific (Fig. 1c right). Levels of igf1r were only decreased after chronic treatment with IGF1. These data are indicative for an involvement of INSR/IGF1R signaling in tumor development and/or progression. However, the different receptor distributions cannot explain the differences in tumor latency time, since this was not a treatment-related effect.
To gain more insight in to overall differences between various insulin analogues in tumor development and progression, we performed RNAseq NGS of 50 tumors (10 tumors per treatment group). We reasoned that transcriptome analysis would shed light on the cause for the differences in tumor latency time by chronic insulin analogue treatment. From the 40 tumors per treatment group, 10 epithelial to mesenchymal transition (EMT)-like tumors were selected per treatment condition for the transcriptome analysis; these 10 EMT tumors had a tumor latency time that was closest to the median latency time of that specific entire treatment group. The selected tumors have been indicated with colored dots in the Kaplan-Mayer curve (Fig. 1b). In Fig. 1d, a correlation analysis shows the sample-to-sample distance. All transcriptomics data were used to define the correlation within and between the different treatment groups. In spite of the relatively high variation in this in vivo experiment, there was a higher correlation within each treatment group than between the different treatment groups. This suggests the existence of a treatment-specific response under each treatment condition. To define these treatment-specific responses, we focused on four particular hallmarks of cancer (Fig. 1e) that might have been involved in the decreased tumor latency, either by tumor initiation or progression. These pathways were selected in close relation to the signaling pathways downstream of IR and IGF1R, and included: i) sustained growth signaling; ii) tissue invasion and metastasis; iii) deregulated cellular metabolism; and iv) genetic instability and mutations.
Chronic IGF1R activation of the MG results in tumors with a ‘sustained growth signaling’ signature
Proliferative signaling is normally a highly regulated process. In a tumor, the replicative cell homeostasis is deregulated causing sustained growth signaling [9]. Since mice were chronically exposed to insulin analogues that activate (to a lesser or greater extent) the IGF1R, we postulated that the decreased tumor latency time (after X10/IGF1 stimulation) would be a direct result of an upregulated IGF1R signaling pathway. We incorporated mouse orthologs of all human genes (~60 genes) that are directly involved in “INSR and IGF1R signaling” according to Ingenuity Pathway Analysis (IPA) (Fig. 2a). A hierarchical clustering analysis was performed for all 50 MG tumors on the relative expression levels of these genes (Fig. 2b). Two clusters could be defined. In cluster 1, tumors from low mitogenic treatments are enriched, whereas cluster 2 predominantly consists of the X10/IGF1 (high mitogenic compound) treated tumors. These two different clusters could not be linked to specific signaling pathways (e.g., PI3K, MAPK, JAK signaling cascade).
To predict the proliferative potential of all the individual tumors, we next built an SVM classifier using transcriptomic signatures associated with proliferation (Fig. 2ci). No significant effect was observed; however, a trend for the tumors in the IGF1 treatment groups to have an increased proliferative potential could be seen. As a positive control for the SVM simulation, similarly we ran the SVM simulation on our earlier micro-array data of the MCF7-IGF1R cells stimulated with the different insulin analogues. We observed that 1 h after glargine, X10, and IGF1 stimulation a strong and significant proliferative genetic signature could be detected. These results are in line with the proliferative potential as determined with functional in vitro assays using this cell line and insulin analogue exposures [5], which makes us believe that the transcriptomic profiles and used SVM simulation is predictive for the biological effect.
An effort was made to determine the growth rate of the individual tumors by calculating the tumor volume from time of detection until the time of sacrifice (data not shown). The EMT tumors in general were highly proliferative resulting in a very small timeframe (a couple of days to weeks) in which tumor volume could be assessed. Unfortunately, this resulted in just a few data points per tumor. Knowing that tumor growth is very likely not linear during tumor development we did not feel comfortable sharing these results. The proliferative potential was also assessed with ki67 staining on several tumors, but as expected there appeared to be a very small window and all tested tumors had a ki67 score >45% indicative of being highly proliferative (data not shown). Next, we attempted to isolate cells from each tumor; we anticipated that monitoring the growth rate of each mouse mammary gland tumor-derived cell line would resemble the growth rate of its originating tumor. Although it appeared to be rather easy to generate a cell line from each individual tumor, it was quite evident that just a small subset of cells would grow out from these tumors. We feared that this clonal expansion would give a skewed image of true tumor growth, and therefore we decided not to include these data.
Chronic IGF1R activation induces tumors with a more mesenchymal phenotype and a higher migration potential
Next, we assessed the hallmark “tissue invasion and metastasis”. We calculated the average number of tumors per mouse per treatment group. Interestingly, we found a significantly higher number of tumors in the IGF1 group (and a non-significant increase for the chronic X10 and glargine treatment groups) compared to the vehicle-treated mice (Additional file 1). This effect can be both due to an increased incidence of multiple primary tumors as well as an increased formation of metastasis in the IGF1 group. To assess the underlying genetic pathways that might have induced metastasis formation we incorporated all genes that are directly associated with the EMT or the mesenchymal to epithelial transition (MET) according to IPA. The gene expression of the mouse orthologs of these genes is presented in an unsupervised hierarchical clustering (Fig. 3a). Two clear gene groups were observed: group A containing genes associated with epithelial phenotype (Elf5, Serpinb5, Cdh1, Grhl2, Elf3, and Wnt4), and group B containing twelve genes associated with mesenchymal cells (Pdgfrb, Six1, Snai2, Tcf4, Zeb1, Klf8, Wnt5a, Snail, Vim, Twist1, and Cdh2). From this hierarchical clustering, three separate treatment clusters could be discriminated. Cluster 1 consists of tumors that show high expression of epithelial- and low expression of mesenchymal-associated genes; three out of six tumors originate from the vehicle-treated animals. Cluster 2 shows moderate expression levels of epithelial- and mesenchymal-associated genes; interestingly the majority of these tumors were from the glargine treatment condition. Cluster 3 consists of tumors that have a low expression of epithelial and high expression of mesenchymal markers. IGF1 treatments were not present at all in cluster 1 and 2, while IGF1 and X10 treatment groups were over-represented in cluster 3. This suggests that IGF1R-mediated signaling is a driver of tumors with an EMT phenotype. We substantiated the phenotypes for the three different clusters (Fig. 3b). The tumors in cluster 1 predominantly expressed epithelial cells that contain many E-cadherin-positive cells. Tumors of the cluster 2 and 3 phenotype were characterized as EMT tumors predominantly consisting of mesenchymal cells, lack of E-cadherin staining, and clear smooth muscle actin (SMA) staining (data not shown) [7].
To link the transcriptomics data further to the phenotype we performed a SVM simulation to derive predictions for the cell migratory potential of each tumor (Fig. 3ci). Overall, no significant difference in the migration potential was observed for the different treatment groups; potentially the intrinsic high variation of in vivo experiments and smaller group sizes can explain lack of power in this analysis. To evaluate whether IGF1 signaling itself is a strong activator in vitro, we also applied our SVM classifier to our previous MCF7-IGF1R transcriptomics data exposed to the different insulin analogues (Fig. 3cii). IGF1, X10, and glargine caused a significant increase in the migration potential. Altogether, these data suggest that chronic exposure to X10, and especially IGF1, induces tumors with a more mesenchymal phenotype that are possibly more aggressive in terms of their migratory potential.
Mice receiving a chronic X10 and IGF1 treatment develop tumors with a higher Warburg potential
Through the activation of the INSR, insulin and insulin analogues directly affect cell metabolism. Therefore, we wondered whether the chronic treatment with insulin analogues induce tumors with a deregulated cellular metabolism or a higher Warburg potential. To test this we used all genes directly involved in “glycolysis” and in “oxidative phosphorylation” according to IPA. An unsupervised hierarchical clustering of the expression profiles of each tumor revealed two very clear gene groups (Fig. 4a). Strikingly, the first gene group consists of eleven genes that are all directly involved in glycolysis; all the genes in the second cluster are involved in oxidative phosphorylation. Four clear tumor clusters could be detected: cluster 1 with high expression of glycolysis genes and low expression of oxidative phosphorylation genes; cluster 2 with moderate glycolysis; cluster 3 with little glycolysis and moderate oxidative phosphorylation; and cluster 4 with high oxidative phosphorylation and little glycolysis. For the insulin group almost all tumors were related to cluster 3 (little glycolysis and moderate oxidative phosphorylation: 90%). For glargine-related tumors, 70% showed enrichment for genes regulating oxidative phosphorylation. IGF1 and X10 tumors showed increased expression of genes involved in glycolysis (60% and 50%, respectively) compared to the vehicle control. This suggests that chronic insulin signaling drives tumors in a more oxidative phosphorylation mode.
Since IGF1 and X10 were mostly related to glycolysis we performed an independent bioinformatics analysis on the activity of various metabolic processes. This revealed several significantly down- and upregulated metabolic programs in IGF1/X10-derived tumors compared to the vehicle/insulin/glargine tumors (Fig. 4b). Importantly, as expected, oxidative phosphorylation was significantly downregulated in the IGF1/X10 group (P = 0.000107) while glycolysis was significantly upregulated. Again, this analysis was also performed on the microarray data of MCF7-IGF1R cells exposed to the different insulin analogues. Similar to the in vivo data, IGF1 and X10 treatment significantly downregulated the citric acid cycle as well as oxidative phosphorylation (both P < 0.0005) and upregulated glycolysis (P < 0.00001) (Fig. 4c). These data suggest that both IGF1 and X10 promote a Warburg effect in mammary gland tumors (chronic exposure) and human tumor cell lines (a single 1 h exposure).
Biomass production rate is highly upregulated in tumors of chronic X10/IGF1-treated mice
Tumor mice receiving a chronic X10 and IGF1 treatment showed increased expression of glycolysis related genes (Fig. 4) and are associated with a decreased tumor latency time (Fig. 5a). We evaluated whether this difference in metabolic capacity is related to enhanced growth. To test this hypothesis, we predicted the rate of the accumulation of biomass in all these different tumors based on the gene expression profiles of biomass-producing metabolic processes. In Fig. 5b, a bar graph of these results is presented. A highly increased mammary gland tumor biomass production rate was observed in chronic X10- and IGF1-treated tumors. This finding is in agreement with the observed increased proliferative potential in the X10/IGF1 tumors. We hypothesize that the decreased tumor latency time by X10 and IGF1 treatment was caused by enhancing the tumor development rather than interfering with the tumor initiation.
Ezh2 and Hras mutations are enriched in chronic IGF1R-activated tumors
Since IGF1 and X10 promote cell proliferation, this could lead to a manifestation of mutations and consequently modulation of tumor development and progression. To evaluate if chronic IGF1 or X10 treatment affects the number of mutations, a mutational analysis was performed on all tumors based on the 40 million 100 base pair reads for each tumor. The number of different mutations was ~3000 mutations per tumor and no significant difference between treatment groups could be detected. Furthermore, there was no correlation between the number of mutations per tumor and the tumor latency time (best fit correlation values, slope = 0.00001207, r
2 = 0.0000099, P value = 0.9827) (Fig. 6a).
Next, we focused on specific clinically relevant mutations that are part of the ~140 known human tumor drivers according to Vogelstein and colleagues [10]. In the 50 mouse tumors that we sequenced, 102 of these tumor drivers were mutated. On average, each tumor had about 35 of these tumor driver mutations and again it seems that the chronic treatment did not affect the overall number of tumor driver mutations per tumor (Fig. 6b). Also, for the tumor driver mutations, there was no correlation between the number of mutations and latency time of each tumor. Furthermore, no treatment-specific effects could be detected for the average number of point mutations, the average number of frame shifts, the average number of start CODON insertions, or the average number of stop CODON insertions (Additional file 2).
We next determined the mutations that are over-represented (>50%) in the X10/IGF1 treatments (Fig. 6c). All mutations are involved in cancer development and/or progression [10, 15]. Although no specific core pathway links these individual mutations, strikingly Ezh2 and Hras were strongly over-represented in IGF1/X10-treated tumors. In a similar way we looked at specific mutations that have been over-represented in insulin analogue treatment altogether and highly under-represented in the vehicle condition. No general mechanism was identified that could explain the enrichment of these mutations in the treatment groups. However, a set of five mutations was not detected at all in the vehicle control, yet these were specifically affected in the insulin treatment group. Overall, these data suggest that the IGF1/X10 enhanced tumor formation is associated with several key candidate cancer drivers that could contribute, including Ezh2 and Hras which are known drivers in human breast cancer [16, 17].