Integrated Bioinformatics Approach Reveals Crosstalk Between Tumor Stroma and Peripheral Blood Mononuclear Cells in Breast Cancer

Breast cancer is a major cause of morbidity and the leading global cause of cancer death in women. Breast cancer constitutes 16% of all female-related cancers yet resulting in approximately half a million deaths every year (2015a). Clinical management relies on known prognostic factors, such as hormone receptor (HR) and HER2 status, for predicting responsiveness to therapies (Kwast et al., 2014; Song et al., 2015). However, Even when morphological characteristics and ER phenotypes are similar, patients have varying prognosis or chemotherapy response (2015b). Gradually increasing resistance against conventional therapies requires improving strategies to diagnose and treat breast cancer. Recent insights into the cancer biology has provided evidence that cancer progression is driven not only by a tumor’s genetic alterations and interactions within the tumor microenvironment (TME), but also by complex systemic processes (McAllister and Weinberg, 2014a). Tumor-derived factors may mobilize host cells from distant tissues, for instance, host circulating cells, the bone marrow and spleen (Shaked et al., 2014), recruiting various peripheral blood cells from the circulation into the TME, resulting impingement on cancer progression


Introduction
Breast cancer is a major cause of morbidity and the leading global cause of cancer death in women. Breast cancer constitutes 16% of all female-related cancers yet resulting in approximately half a million deaths every year (2015a). Clinical management relies on known prognostic factors, such as hormone receptor (HR) and HER2 status, for predicting responsiveness to therapies (Kwast et al., 2014;Song et al., 2015). However, Even when morphological characteristics and ER phenotypes are similar, patients have varying prognosis or chemotherapy response (2015b). Gradually increasing resistance against conventional therapies requires improving strategies to diagnose and treat breast cancer.
Recent insights into the cancer biology has provided evidence that cancer progression is driven not only by a tumor's genetic alterations and interactions within the tumor microenvironment (TME), but also by complex systemic processes (McAllister and Weinberg, 2014a). Tumor-derived factors may mobilize host cells from distant tissues, for instance, host circulating cells, the bone marrow and spleen (Shaked et al., 2014), recruiting various peripheral blood cells from the circulation into the TME, resulting impingement on cancer progression RESEARCH ARTICLE
In addition to tumor-driven systemic perturbations , the stroma cells of TME be educated and sculpted by tumor cells via paracrine and juxtacrine, also have capacity of promoting tumor progress through systemic effects (Quail and Joyce, 2013). For example, Cancer-associated fibroblasts (CAFs) derived circulating CXCL12 was shown to mobilize the progenitor cells into the circulation, subsequently leading to their recruitment into the TME to promote angiogenesis (Hattori et al., 2001;Orimo et al., 2005). These TME-driven systemic perturbations , which can also influence disease outcome, remain ambiguous.
However, the composition of TME is very complicated, including immune cells, CAF, blood vessels, lymphatic vessels. This cell-type heterogeneity makes it very difficult to study their crosstalk with peripheral blood cells via coculture. The genome-wide ''omics'' seems more important due to the absence of such an ''experimental testing ground''. Re-analysis of large data sets with advanced analysis tools will offer valuable clues.
Herein, we analyzed two mRNA expression profiling datasets, including peripheral blood mononuclear cells (PBMCs) and tumor stroma samples from breast cancer patients respectively. To provide biologically meaningful results, we propose comprehensive bioinformatics approaches to construct crosstalk network between the breast tumor stroma and the peripheral blood cells.

Data sources
Microarray data were downloaded from Gene Expression Omnibus (GEO) (Barrett et al., 2013), Accession number GSE9014, includes 53 cases of tumor stroma laser-capture microdissected from IDC breast cancer cases, and 31 cases of individual-matched normal adjacent stroma. Genes whose expression varied most between tumor tissue and normal stroma for the 31 tissue-matched pairs of this data had been identified to be independent of ER, HER2 and lymph node status, as well as age, grade and tumor size (Finak et al., 2008). Accession number GSE27562 (LaBreche et al., 2011), includes gene expression analysis of PBMCs from 57 women with a diagnosis of breast cancer and 31 healthy women with normal mammograms.

Basic stages of comprehensive bioinformatics approach
I. Data preparation: Obtain transcriptional profiling data from peripheral blood and stromal tissue compartments from samples of breast cancer patients.
II. Crosstalk components analysis: connect network of ligands-receptors-interacting proteins. Identify the key components as follows: i.Upregulated ligands in tumor stromal cells compared to normal counterparts; ii.Their corresponding receptors expressed in peripheral blood cells of breast cancer patients; iii. Intracellular interacting proteins of these receptors expressed in peripheral blood cells of breast cancer patients; III. Network analysis: Perform network analyses to identify network hubs, differential network features, and differential functional enrichment.

Data Preprocessing and Identification of differentiallyexpressed genes (DEGs)
The raw data were analyzed using the integrated software package BRB-Array Tools version 4.0 (2015c). BRB-Array Tools performed a series of preprocessing steps including filtering and computing hybridization data of each probe in the probe set. Probes with the following conditions were removed: 1) greater than 20% of expression data values had more than 1.5-fold change in either direction from the median value; 2) more than 50% of the gene expression data were "absent". The Robust Multiarray Average (RMA) method (Irizarry et al., 2003) was used for normalization. The SAM (Statistical Analysis of Microarray) method (Tusher et al., 2001) was used to analyze the transcription profiles and screen for DEGs (False Discovery Rate [FDR]=1%, fold change [FC]≥1.5).

Screening of DEGs for ligands and their cognate receptors
The significantly upregulated ligands of breast cancer stroma and their cognate receptors were sorted based on two databases, the Database of Ligand-Receptor Partners (DLRP) (Graeber and Eisenberg, 2001), a database of protein ligand and protein receptor pairs that are known to interact with each othe, and The Universal Protein Resource (UniProt) (2015d), UniProtKB entry: Cellular component/Subcellular locations/Secreted. The cognate receptors of these ligands were identified by analyzing protein interactions from the two databases simultaneously. To verify interaction relationships the ligand-receptor pairs were evaluated manually using published literature.

Selecting intracellular proteins interacting with receptors
The interaction proteins of the receptors was selected based on the UniProt database and the Search Tool for the Retrieval of Interacting Genes (STRING) database (Franceschini et al., 2013) (combined score >0.4, proteins with subcellular localizations in ''extracellular'' and ''secreted'' were excluded). The crosstalk network of ligands-receptors-interacting proteins was integrated and visualized by Cytoscape software (Shannon et al., 2003). The property of the network was evaluated using the network analyzer plugin of Cytoscape.

GO and KEGG pathway enrichment analysis
To uncover further insights into the precise biological function and signaling pathways involved with the crosstalk genes identified in the present study, Gene Ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis were performed for the genes in the crosstalk network. The online based tool of the database for annotation, visualization and integrated discovery (DAVID) (Huang et al., 2009) was used to perform these analyses with P<0.05.

Differential up-regulation ligands in breast cancer stroma samples
A total of 3420 transcripts, corresponding to 2737  DOI:http://dx.doi.org/10.7314/APJCP.2016.17.3.1003 Integrated Bioinformatics -Crosstalk Between Tumor Stroma andPeripheral Blood Mononuclear Cells in Breast Cancer DEGs, including 1122 upregulated DEGs and 1615 downregulated DEGs were identified between the breast cancer stroma samples and the controls. Among these upregulated DEGs, 19 were define as ligand according to the two databases, DLRP and UniProt. There were 26 corresponding receptors expressed in PBMCs were retrieved as potential crosstalk components. The complete list of ligand-receptor pairs were show in Table 1.

GO and pathway enrichment analyses
We used DAVID to analyze all proteins in the interaction network. A total of 359 GO functional terms were enriched, the top 20 functional terms ranked by statistical significance were listed in Table 2. We found that they were mainly involved in wound healing, inflammation, cell adhesion and migration. Top 10 KEGG pathway are shown in Table 3. ECM-receptor interaction, focal adhesion, chemokine signaling pathway were included in the list.

Discussion
Over the past two decades, the TME have collectively risen in prominence, embodied in the view that, as a systemic disease, cancer progression may be directed by the body's systemic responses to malignancy and by the involvement of organ systems located at sites distant from the site of primary tumor growth (McAllister and Weinberg, 2014b). Up to now, most analyses of such systemic instigation by tumor have concentrated on either the tumor tissues or the tumor cells, few have examined the dialogue between the tumor stroma (the local tumor environment) and the systemic tumor environment. The primary aim of this study was to decipher the secretome of breast tumor stroma cells, and the regulation of circulating cell function by stroma cell-secreted proteins.
In the present study, using mRNA expression profile, comparison of the transcriptomes of breast tumor stroma and normal control tissue was made, candidate genes potentially involved in the stroma-circulating cell dialogue were identified. There were 19 ligands that had up-regulated in breast tumor stroma, for which 26 respective receptors expressed on the peripheral blood cells of breast cancer patients. In order to explore the potential downstream target proteins of the receptors in PBMCs, we predicted PPIs of the receptors using the STRING and UniProt databases. A total of 73 interactive proteins were identified. Finally, a crosstalk network was constructed by integrating ligands-receptors-interacting proteins relationships, and was visualized using Cytoscape.
Based on the GO term and KEGG pathway analyses, the present study obtained a comprehensive understanding of the crosstalk proteins (genes). The most significant biological processes and pathway indicated an abundance of molecules contributing to a proinflammatory environment, invovled response to wounding, leukocyte activation, adhesion and migration. In cancer, TME consists of different entities, such as infiltrating immune cells (IICs), cancer-associated fibroblastic cells (CAFs), and mesenchymal stem cells (MSC) (Hanahan and Coussens, 2012;Maenhout et al., 2014). The mobilization and subsequent trafficking of these cells to tumors is thought to be due to inflammatory signaling in a tumor resembling that of an unresolved wound (Spaeth et al., 2008). For instance, among the leukocytes that infiltrate the tumor microenvironment, myeloid-derived suppressor cells (MDSCs) can be recruited to the tumor site from the peripheral blood, by a number of chemokines produced by tumors (Ley et al., 2007;Ding et al., 2015;. Our result indicated that, in addition to tumor cell itself , stromal cells may also contribute significantly to recruitment, activation, programming of those cells. The hub proteins of the crosstalk network with degree more than 10 were ITGB1, ITGB3, ITGA2B, CXCR4. The first three were integrins. The integrins are expressed in many leukocytes such as T-lymphocytes, monocytes, and granulocytes and play an integral role in leukocyte migration by binding to endothelial cells and stimulating extravasation (Alon and Feigelson, 2012;Kuehn et al., 2014). Activation of integrin α4β1, leading to the extravasation of myeloid lineage cells from the circulation, and recruit to the tumor microenvironment (Schmid et al., 2013). CXCR4 is a receptor specific for CXCL12 (also known as stromal-derived factor 1α). Increasing evidence supports a critical role of the CXCL12/CXCR4 axis for cells trafficking and recruiting to tumor site (Domanska et al., 2013). These different cell types including Endothelial Progenitor Cells (de la Puente et al., 2013), myeloid bone marrow-derived cells (BMDCs) (Hiratsuka et al., 2011), neutrophil (Seubert et al., 2015), T-regulatory cells and M2-type macrophages (Gil et al., 2014), mast cell (Ellem et al., 2014), mesenchymal stem cell (Lourenco et al., 2015), B cell (Shetty et al., 2012).
Overall, we carried out an integrated bioinformatics analysis of crosstalk between tumor stroma and the peripheral circulation. Our results provide new insight into the communication between local tumor environment and the host systemic environment. We hypothesize that secretion of these proteins by the tumor stroma and binding to their receptors in the PBMCs may play a role in exerting systemic changes which lead to the mobilization