Systems Analysis of a Mouse Xenograft Model Reveals Annexin A1 as a Regulator of Gene Expression in Tumor Stroma

Annexin A1 is a multi functional molecule which is involved in inflammation, innate and adaptive immune systems, tumor progression and metastasis. We have previously showed the impaired tumor growth, metastasis, angiogenesis and wound healing in annexin A1 knockout mice. While tumor is a piece of heterogeneous mass including not only malignant tumor cells but also the stroma, the importance of the tumor stroma for tumor progression and metastasis is becoming increasingly clear. The tumor stroma is comprised by various components including extracellular matrix and non-malignant cells in the tumor, such as endothelial cells, fibroblasts, immune cells, inflammatory cells. Based on our previous finding of pro-angiogenic functions for annexin A1 in vascular endothelial cell sprouting, wound healing, tumor growth and metastasis, and the previously known properties for annexin A1 in immune cells and inflammation, this study hypothesized that annexin A1 is a key functional player in tumor development, linking the various components in tumor stroma by its actions in endothelial cells and immune cells. Using systems analysis programs commercially available, this paper further compared the gene expression between tumors from annexin A1 wild type mice and annexin A1 knockout mice and found a list of genes that significantly changed in the tumor stroma that lacked annexin A1. This revealed annexin A1 to be an effective regulator in tumor stroma and suggested a mechanism that annexin A1 affects tumor development and metastasis through interaction with the various components in the microenvironment surrounding the tumor cells.


Introduction
We recently showed the pro-angiogenic functions in tumor development for annexin A1 [1] which was previously known as an inflammatory protein. Tumor growth and metastasis were significantly decreased when tumors grew in host animals unable to express annexin A1 [1]. While tumor is a piece of heterogeneous mass including not only malignant tumor cells but also the stroma, the importance of the tumor stroma for tumor progression and metastasis is becoming increasingly clear. The tumor stroma is comprised by extracellular matrix and non-malignant cells in the tumor, which include endothelial cells, fibroblasts, immune cells, inflammatory cells such as macrophages [2,3]. Although the tumor stroma has been sought after as a therapeutic target, the mechanisms for how the interaction of various components in the complex tumor stroma contributes to the tumor development are still poorly understood. Based on our recent finding of proangiogenic functions for annexin A1 in vascular endothelial cell sprouting, wound healing, and tumor growth and metastasis [1], and the previously known properties for annexin A1 in immune cells and inflammation [4], this study hypothesized that annexin A1 is a key functional player in tumor development, linking the various components in tumor stroma by its actions in endothelial cells and immune cells. Here, we used systems biology approach to analyze the tumors in whole animal models on the background of absence or presence of annexin A1 to show the global effects of annexin A1 on tumor stroma. It is impossible to study tumor stroma with isolated proteins, cultured cell homogenates or even whole live cells. Tumor cells in culture lack the microenvironment or the stroma provided by the body that hosts the tumor. Our whole animal models offered significant advantages over any other conditions. Of particular importance our systems based approaches provided thorough analysis and generated interesting and thought-provoking results.

Mice and Tumor Models
Annexin A1 knockout homozygous and congenic wildtype counterpart homozygous mice, tumor cell culture, B16 melanoma cell line and subcutaneous tumor model were as described previously [1]. Briefly, to grow the xenograft tumors on the mice, B16 melanoma cells were grown in continuous culture for no more than 3 consecutive passages, the actively growing cells were then trypsinized (0.25% trypsin/1 mM Na-EDTA; Gibco/BRL), resuspended in DMEM (Dulbecco modified Eagle medium), counted, and injected subcutaneously into the right flanks of the mice at 5610 6 cells in 200 mL DMEM for each mouse. The B16 cell line was an established cell line and was obtained and used in Dr. Schnitzer's lab in Sidney Kimmel Cancer Center and Proteogenomics Research Institute for Systems Medicine, San Diego, California, United States of America, and the B16 cell line related work was previously published [1].

Systems Analysis and Gene Microarray
The GeneChip mouse genome 430 2.0 array (Affymetrix) was used to analyze tumors from annexin A1 knockout or wildtype mice as described previously [1]. The dataset containing genes that are differentially expressed in the tumors from annexin A1 knockout mice and wildtype mice was analyzed using Genomatix software (Genomatix Software) and Gene Ontology database [5] and IPA software (Ingenuity Systems).
Study Approval. This study was approved by Sidney Kimmel Cancer Center and Proteogenomics Research Institute for Systems Medicine, San Diego, California, United States of America.

Results
The tumors were grown in both annexin A1 null and wild type mice using B16 melanoma cells as described previously [1]. In this system, only the stroma in the tumor differs because the nontumor cells are from the host body with or without annexin A1. From the raw data obtained previously [1] using mouse whole genome microarrays to run tumor samples from annexin A1 knockout and wildtype mice, this study performed thorough systems analysis. The gene expression levels in the whole tumors from annexin A1 knockout and wildtype mice were analyzed to obtain differential expression profiling (Fig. 1). For each gene probe in the mouse genome, the expression levels of the gene probe in tumors from annexin A1 null mice (ko) and from annexin A1 wild-type mice (wt) were compared and expressed as ratio of ko over wt, which reflects from three independent repeat experiments. In tumors from annexin A1 null mice, those genes with highest values were over-expressed (red colored region in Fig. 1) while those genes with lowest values were under-expressed (blue colored region in Fig. 1). The majority of the genome was expressed at similar extent in tumors from annexin A1 null and wild-type mice with the expression ratio values approaching one or the logarithm of the expression ratio values close to zero (Fig. 1).
Total about 1000 genes with either highest or lowest ratio values were identified as significantly changed genes in the tumors.
We then analyzed this list of significantly changed genes/ proteins in a systems biology approach with Genomatix and IPA software to reveal the functional connection with annexin A1. The data were examined in terms of biological processes category, which are the Gene Ontology (GO) structured networks of defined terms to describe gene product attributes. The data were statistically analyzed and expressed as the z-score of a term to check for over-or under-represented groups of genes. Integrated statistical rating allows for the identification and selection of clusters of functionally related genes. All the biological processes in the GO database were obtained from top of the hierarchial tree to further mining down of subcategories of each branch levels by levels. The top most general biological processes are shown with their z-scores (Fig. 2). Interestingly, most of the biological processes were positively enriched, suggesting that annexin A1 may act as a global regulator which affects many biological processes and some of these processes may balance each other to result in minimum changes, which could contribute to the annexin A1 knockout mice appearing normal overall. The z-score of a term shows whether a certain gene, or group of genes, is over-or under-represented in the set of the genes being analyzed. Among the total 22 categories, the immune system process was most significantly over-represented (z-score, 20.63). This indicates that a group of genes, which are functionally related to immune system, were over-represented in the list of the differentially expressed genes identified above. This is consistent with the previously known functions of annexin A1 in immune system [4], thus the other genes related to immune system changed their expression levels in the tumors grown in the host microenvironment where annexin A1 was knocked out. Breakdown of this immune system process into its subcategories is shown in detail with the biological processes down levels by levels ( Fig.S1) and shown in summary with representative biological processes (Fig. 3). The process of immune response to tumor cell was not much different between tumors in annexin A1 ko and wild-type conditions, consistent with our notion that it was the tumor stroma that made the tumor development significantly different in annexin A1 negative condition. Among the immune cells that were present in tumor stroma, T cell related processes stood out, which is consistent with the recently identified effects of annexin A1 on T cell functions [6].
Besides immune system process, other significantly overrepresented general biological processes include cell killing, cellular process, biological adhesion, multicellular organismal process, developmental process, locomotion, response to stimulus, biological regulation, positive regulation of biological process, negative regulation of biological process, and regulation of biological process (Fig. 2). All of these processes were further mined down to their subcategories and the summary of the analyses were shown (Figs. 4, 5, 6, and 7) with the detailed analyses for each subcategory shown in supplementary information section (Figs.S2, S3, S4, S5, S6, S7, and S8). Mining down in the cell killing process (Fig. 4A), again, the category T cell mediated cytotoxicity was over-represented. Further highlighting the involvement of the immune, inflammatory components in the tumor stroma, the cellular process brought out the enriched process of macrophage differentiation (Fig. 5). Biological adhesion process (Fig. 4B) showed an over-representation of leukocyte adhesion, in line with the interaction between extracellular matrix and infiltrated non-tumor cells in tumor stroma. The process of tumor necrosis factor production was over-represented (Fig. 4C), consistent with our previous findings of increased necrosis in tumors grown in annexin A1 null mice [1], the change in tumor Figure 1. Gene differential expression profile of tumors from annexin A1 knockout and wild type mice. B16 melanoma cells were implanted in annexin A1 knockout and wild type mice to grow up tumors. Whole tumors were analyzed using mouse whole genome microarray, in which for each gene, one or more than one probeset(s) were spotted for the gene. For each probeset on the microarray, the difference in gene expression between tumors from annexin A1 knockout mice and from wild type mice was expressed as log 2 ratio of the gene expression level in annexin A1 knockout (ko) over the gene expression level in wild type (wt). doi:10.1371/journal.pone.0043551.g001 necrosis factor production may contribute to the amount of necrosis in tumors. As shown previously [1], angiogenesis was over-represented ( Fig. 4D), the process involving endothelial cells in tumor stroma. The locomotion process (Fig. 4E) indicated an over-representation of leukocyte chemotaxis and neutrophil chemotaxis, both processes involving immune cells in tumor stroma. Again in the response to stimulus process (Fig. 6), the categories related to inflammatory response were significantly over-represented. It is noteworthy that the acute inflammatory response was also over-represented but not the chronic inflammatory response. Finally, biological processes related to regulation ( Fig. 7), with overlapping subcategories showed a general overrepresentation of regulation for those over-represented biological processes discussed above. The positive regulation of biological process was more over-represented than the negative regulation of biological process.
Furthermore, we identified the molecules that were responsible for those over-represented biological processes to reveal an interactome in tumor stroma with annexin A1 being a key functional player. Using Genomatix and IPA software, we obtained lists of biological processes or functions with a list of molecules in each process or function. We combined the lists of Figure 2. Biological process. The list of significantly changed genes was analyzed with Genomatix software and examined in terms of biological processes in Gene Ontology database and expressed as the Z score of a term for each biological process within the category. The higher value of Z score means more genes in that function of biological process were affected in this experiment, thus this biological process was over-represented. The Z scores for all the biological processes in top most general level of Gene Ontology structural networks of hierarchial tree of categories of biological process were shown here. doi:10.1371/journal.pone.0043551.g002 molecules from the processes or functions that involve various components in tumor stroma to generate the list of molecules that comprises the tumor stroma interactome which interacts with annexin A1 (Table 1). These molecules showed most significant changes, either most down-regulated (log2 (ko/wt),20.999), or most up-regulated (log2 (ko/wt).0.956), in the tumor stroma that lacked annexin A1. Based on the current literature data so far, most of the molecules on this list formed an interaction network built by IPA software (Fig. 8), of which, the interaction between annexin A1 and integrin beta 2 and the interaction between annexin A1 and vascular cell adhesion molecule 1 were identified based on current data. This interaction may occur in cell plasma membrane (Fig. 9). The expression ratio of integrin beta 2 in tumors from annexin A1 ko versus wt was 0.378 (Table 1), so the integrin beta 2 was down-regulated in the absence of annexin A1; and the vascular cell adhesion molecule 1 was even downregulated (ratio of ko versus wt was 0.214).
While the tumor cells were originally same before being implanted into annexin A1 ko and wt mice, therefore, the tumor stroma in annexin A1 ko and wt mice exhibited significant differences in the gene expression profiles between the host body with and without annexin A1, although the genes in tumor cells  Figure 2, the most over-represented immune system process was further mining down levels by levels into its subcategories with all biological processes in each level shown in Figure S1 and from each level, the over-represented biological processes are summarized here. Each color represents each level in the order from top to bottom as mining down the levels from top level of the hierarchial tree in the immune system process. doi:10.1371/journal.pone.0043551.g003  Figure 3, the most general processes were further mining down levels by levels into its subcategories with all biological processes in each level shown in Figures S2, S3, S4, S5, and S6 and from each level, the overrepresented biological processes are summarized here. Each color represents each level in the order from top to bottom as mining down the levels from top level of the hierarchial tree in these biological processes. doi:10.1371/journal.pone.0043551.g004 may change their expression levels due to the interaction between tumor cells and the tumor stroma. Annexin A1 is abundant in neutrophils, monocytes, macrophages in wild type animal models [4]. In the annexin A1 null animal model, these non-tumor cells infiltrated into tumor stroma do not express annexin A1, this lack of annexin A1 protein function reflected as the gene expression level changes of several other proteins. This network of changes indicated these important genes/proteins ( Table 1) that are interacting with annexin A1 directly or indirectly. The interaction of annexin A1 with these proteins was previously unknown.

Discussion
Through systems biology analysis, we found this list of genes ( Table 1) that significantly changed in the tumor stroma that lacked annexin A1. This showed a mechanism that annexin A1 affects tumor development and metastasis through interaction with the various components in the microenvironment surrounding the tumor cells. Literature has also showed that these genes are associated with tumor and/or tumor stroma: IGH-6 (ko/ wt = 0.065) with thymic lymphoma and lymphoblastic leukemia [7][8], TGFBI (ko/wt = 0.108) with neuroblastoma, lung carcinoma, ovarian and prostate cancers, renal, gastrointestinal and brain tumors, colon cancer [9][10][11][12][13][14], PTPRC (ko/wt = 0.143) with human Figure 5. Cellular process. Similarly, from Figure 2, the cellular process was further mining down levels by levels into its subcategories with all biological processes in each level shown in Figure S7 and from each level, the over-represented biological processes are summarized here. Each color represents each level in the order from top to bottom as mining down the levels from top level of the hierarchial tree in the cellular process. doi:10.1371/journal.pone.0043551.g005 breast tumor stroma [15], MMP13 (ko/wt = 0.147) with skin tumor stroma [16], CCK (ko/wt = 3.347) with human pancreatic cancer [17], TFRC (ko/wt = 2.967) with colon cancer, human esophageal squamous cell carcinoma, human brain tumors, mouse mammary adenocarcinoma and rodent liver tumor [18][19][20][21][22], TIMP2 (ko/wt = 2.224) with human colorectal cancer, human pancreatic carcinoma, human renal cell carcinoma and human squamous cervical carcinoma [23][24][25][26], BRCA1 (ko/wt = 1.917) with human breast tumor stroma [27].
It is worth noting that TGFBI, as an extracellular matrix protein, has been reported differently either as a tumor suppressor [12,28] or over-expressed in tumors [13][14]. TGFBI was shown to promote metastasis by acting on tumor stroma [14]. Here we showed that TGFBI was under-expressed in tumors from annexin A1 ko mice. In tumors grown on annexin A1 ko mice, the metastasis was reduced [1]. Thus, it suggests that annexin A1 interacts with or otherwise affects TGFBI and this interaction affects the metastatic ability of tumors. In addition, loss of TGFBI can induce resistance to chemotherapeutic agent paclitaxel in ovarian cancer cells [29], while human tumor cells with overexpressed levels of TGFBI showed an increased sensitivity to etoposide, paclitaxel, cisplatin and gemcitabine [10]. Similarly, stroma-derived annexin A1 has been shown to play a role in cirradiation-induced T-cell lymphoblastic lymphoma development, Figure 6. Response to stimulus. From Figure 2, the response to stimulus process was further mining down levels by levels into its subcategories with all biological processes in each level shown in Figure S8 and from each level, the over-represented biological processes are summarized here. Each color represents each level in the order from top to bottom as mining down the levels from top level of the hierarchial tree in the response to stimulus process. doi:10.1371/journal.pone.0043551.g006   with annexin A1 being a candidate resistance gene against cradiation exposure [30]. It has been reviewed that the cellular and non-cellular components of the tumor microenvironment contribute to the chemoresistance [31]. Therefore, annexin A1 and TGFBI may be further studied together to elicit a mechanism for drug resistance. Annexin A1 and TGFBI could be important therapy targets to manipulate to improve the efficacy of radio-and chemo-therapies for cancer patients, to reduce resistance to treatments and decrease recurrence of cancer. Protein tyrosine phosphatase gene expression was shown to be up-regulated in stromal fibroblasts from human breast tumors [15], while our study showed PTPRC down-regulated in tumor stroma without annexin A1 (Table 1), this suggests that inhibiting annexin A1 may disrupt tumor stroma to hamper tumor progression via decreasing expression level of protein tyrosine phosphatase. MMP13 was shown to promote angiogenesis [16], which is consistent with the finding in this study of down-regulated MMP13 expression (Table 1) and impaired tumor growth and angiogenesis in tumors on annexin A1 knockout mice [1].
The role of endogenous cholecystokinin (CCK) in human pancreatic cancer is not clear. The pancreatic cancer cells produced both CCK and gastrin, however the CCK level was lower than the gastrin. It seemed gastrin played a more dominant role than CCK in stimulating tumor growth. While downregulation of gastrin inhibited growth of pancreatic tumor, change in the CCK level did not affect the tumor growth [17]. Here we showed a huge increase of CCK mRNA expression in tumors from annexin A1 ko mice, the further study of interaction of annexin A1 and CCK will shed light on understanding the role of CCK in human cancers.
Also, this study showed a possible interaction between annexin A1 and BRCA1. It has been shown that in the breast cancers containing mutated BRCA1, a common breast cancer susceptibility gene, changes in the tumor stroma facilitate the malignant transformation of the tumor cells [32]. BRCA1 is also believed to be a regulator in mammary stem cell differentiation and associated with cancer stem cells. One of the identified markers for selection of human stem cells and cancer stem cells is aldehyde dehydro- genase 1 (ALDH1), however, the exact function of ALDH1 in stem cells is still mostly unknown. ALDH1 expression is significantly increased in both tumor cells and tumor stroma in breast cancers carrying BRCA1 mutations [33]. The cancer stem cells possess tumor initiating capacity and therapy resistance. Therefore, consistent with the above mentioned annexin A1 being candidate resistance gene, an association between annexin A1 and stem cells and cancer stem cells is strongly possible and worth further investigation to help understand the exact functions of stem cell markers.  Table 1    Similarly as Figure S1, the top level category, multicellular organismal process, labeled (A), was further mining down levels by levels into its subcategories labeled alphabetically with each letter for each down level and for each level, representative categories were further broken down into all its subcategories shown here. Similarly as Figure S1, the top level category, cellular process, labeled (A), was further mining down levels by levels into its subcategories labeled alphabetically with each letter for each down level and for each level, representative categories were further broken down into all its subcategories shown here.  Figure S1, the top level category, response to stimulus, labeled (A), was further mining down levels by levels into its subcategories labeled alphabetically with each letter for each down level and for each level, representative categories were further broken down into all its subcategories shown here.