Yu Zhang1†
Sen Zhang- 1Department of Orthopaedics, The First People's Hospital of Chengdu, Chengdu, China
- 2Department of Stomatology, Shenzhen Hospital, Southern Medical University, Shenzhen, China
- 3State Key Laboratory of Bioactive Substances and Functions of Natural Medicines, Institute of Materia Medica, Chinese Academy of Medical Sciences & Peking Union Medical College, Beijing, China
Background: Bone metastasis frequently occurs in advanced breast cancer patients, and it is one of major causes of breast cancer associated mortality. The aim of the current study is to identify potential genes and related signaling pathways in the pathophysiology of breast cancer bone metastasis.
Methods: Three mRNA expression datasets for breast cancer bone metastasis were obtained from Gene Expression Omnibus (GEO) dataset. The differentially expressed genes (DEGs) were obtained. Functional analyses, protein-protein interaction (PPI) network, and transcription factors (TFs)-target genes network was constructed. Real-time PCR using clinical specimens was conducted to justify the results from integrated analysis.
Results: A 749 DEGs were obtained. Osteoclast differentiation and rheumatoid arthritis were two significantly enriched signaling pathways for DEGs in the bone metastasis of breast cancer. SMAD7 (degree = 10), TGFBR2 (degree = 9), VIM (degree = 8), FOS (degree = 8), PDGFRB (degree = 7), COL5A1 (degree = 6), ARRB2 (degree = 6), and ITGAV (degree = 6) were high degree genes in the PPI network. ETS1 (degree = 12), SPI1 (degree = 12), FOS (degree = 10), FLI1 (degree = 5), KLF4 (degree = 4), JUNB (degree = 4), NR3C1 (degree = 4) were high degree genes in the TFs-target genes network. Validated by QRT-PCR, the expression levels of IBSP, MMP9, MMP13, TNFAIP6, CD200, DHRS3, ASS1, RIPK4, VIM, and PROM1 were roughly consistent with our integrated analysis. Except PROM1, the other genes had a diagnose value for breast cancer bone metastasis.
Conclusions: The identified DEGs and signaling pathways may make contribution for understanding the pathological mechanism of bone metastasis from breast cancer.
Introduction
Breast cancer, a complex and heterogeneous disease, is the leading cause of cancer-associated death in the women worldwide (1). In China, breast cancer exhibits one of the fastest growing incidence of cancer, which has become the leading cause of mortality in females (2). Clinically, metastasis is a major challenge that leads to 80% of cancer-associated death (3). It is worth mentioning that most patients with breast cancer die from tumor metastases. It is found that bone is one of the most frequent metastatic site for breast cancer, with about 70% of patients with advanced breast cancer exhibiting bone metastasis (4). Furthermore, bone micrometastases are present in approximately 30% of breast cancer patients diagnosed with stage I, II, and III (5). These breast cancer patients with bone metastasis also have substantial complications, including pain, hypercalcemia, and increased risk of fracture (6).
It is reported that the pathology of bone metastasis is characterized by a vicious cycle where the balance between bone forming cells “osteoblasts” and the bone resorbing cells “osteoclasts” is imbalanced (7). For women with early breast cancer, younger age, primary tumor size >2 cm, estrogen receptor positive/progesterone receptor negative, the presence of significant nodal disease are risk factors for developing bone metastases (8, 9). Generally, it is incurable except palliative therapies when the tumor has metastasized to bone (10). Although significant progress in long—term survival has been made, early diagnosis, prognosis of breast cancer remain poor (11). Therefore, identifying potential pathological genes and signaling pathways influencing the bone metastasis of breast cancer process is very important.
In this study, we aimed to find differentially expressed genes (DEGs) in breast cancer bone metastasis by integrated analysis. Then, functional enrichment analysis including Gene Ontology (GO) and the Kyoto Encyclopedia of Genes Genomes (KEGG) was used to investigate the biological function of DEGs. Protein-protein interaction (PPI) and transcriptional factors (TFs)-target genes network construction was performed to further investigate the function of identified DEGs. Quantitative reverse transcription-polymerase chain reaction (QRT-PCR) was applied to validate the expression of candidate DEGs. Finally, receiver operating characteristic (ROC) analyses was applied to analyze diagnosis ability of identified DEGs. Our study may be helpful in understanding the pathogenic mechanism of bone metastasis from breast cancer, and most difference from previous studies is that current study focuses on bone metastasis of breast cancer and aim to find the differential expressed genes and pathways from other site metastasis of breast cancer.
Materials and Methods
Datasets
In this study, we searched datasets of breast cancer bone metastasis from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/) with the keywords (“breast neoplasms” [MeSH Terms] OR breast cancer [All Fields]) AND (“bone and bones” [MeSH Terms] OR bone [All Fields]) AND “Homo sapiens”[porgn] AND “gse”[Filter]. The dataset will be screened for the integrated analysis; (1) The study type is “expression profiling by array.” (2) All selected datasets were from breast cancer metastasis of bone and/other sites metastasis of tumor samples. (3) The standardized or primary datasets were included in this study. Finally, a total of 3 datasets including GSE14020, GSE54323, and GSE46141 was screened, which was shown in Table 1.
Table 1. Three datasets of selected breast cancer with bone metastasis.
Identification of DEGs
Raw expression data in this study was downloaded. Limma and Meta-analysis for MicroArrays (metaMA) packages were used to identify the DEGs, and the inverse normal method was used to combine the p-value in metaMA. The false discovery rate (FDR) was performed for multiple testing corrections of raw p-value through the Benjamin and Hochberg method (12, 13). The threshold of DEGs was set at p-value < 0.01.
Functional Enrichment Analysis of DEGs
In order to obtain the biological function and signaling pathways of DEGs, the GeneCoDis3 (http://genecodis.cnb.csic.es/analysis) software was used for Gene Ontology (GO, http://www.geneontology.org/) annotation and Kyoto Encyclopedia of Genes Genomes (KEGG, http://www.genome.jp/kegg/pathway.html) pathway enrichment of DEGs. The threshold of GO function and KEGG pathway of DEGs was all set as FDR < 0.05.
PPI Network Construction
It is useful to understand the molecular mechanism of bone metastasis from breast cancer by studying the interactions between proteins. In order to gain insights into the interaction between proteins encoded by DEGs and other proteins, the database of BioGRID (http://thebiogrid.org) was used to retrieve the predicted interactions between all proteins encoded by DEGs and other proteins. Then, PPI network was visualized by the Cytoscape Software (http://cytoscape.org/). A node in the PPI network denotes protein, and the edge denotes the interactions.
Analysis of Potential TFs to DEGs
TFs play crucial roles in regulating gene expression. We downloaded the TFs in the human genome and the motifs of genomic binding sites from the TRANSFAC. Moreover, the 2 KB sequence in the upstream promoter region of DEGs was downloaded from UCSC (http://www.genome.ucsc.edu/cgi-bin/hgTables). Target sites of potential TFs were subsequently distinguished. Finally, the transcriptional regulatory network was visualized by Cytoscape software.
QRT-PCR in vitro
Twenty-five patients with bone metastasis of breast cancer, 25 patients with lung metastasis of breast cancer, and 25 patients with liver metastasis of breast cancer were incorporated in our study. Meanwhile, 25 primary breast cancer patients without metastasis was enrolled, and their adjacent non-cancer breast tissues was considered as normal tissue to be used as control in the present study. Their tumor tissues from metastatic sites were used for RNA isolation. Ethical approval was obtained from the ethics committee of the cancer hospital and informed written consent was obtained from all of subjects.
Total RNA of the tumor tissues from metastatic sites was extracted using the RNA liquid Reagent (TIANGEN) according to the manufacturer's protocols. One microgram RNA was applied to synthesize cDNA by Fast Quant Reverse Transcriptase (TIANGEN). Then real-time PCR was performed in an ABI 7300 Real-time PCR system with SYBR® Green PCR Master Mix (TIANGEN). All reactions were carried out in triplicate. GAPDH was used for the internal reference. Relative gene expression was analyzed by fold change method.
Receiver Operating Characteristic Analyses
In order to assess the diagnostic value of IBSP, MMP9, TNFAIP6, DHRS3, RIPK4, CSPG4, and CD200 in bone metastasis of breast cancer, receiver operating characteristic (ROC) analyses were performed using pROC package in R language. The area under the curve (AUC) under binomial exact confidence interval was calculated to generate the ROC curve.
Cell Culture and Transfection With Plasmid and siRNA
Human breast cancer cell line MCF-7 was purchased from ATCC (Manassas, VA, USA) and was cultured in RPMI 1640 medium (Invitrogen, CA, USA) supplemented with 10% fetal bovine serum (Invitrogen, Carlsbad, CA, USA) and incubated at 37°C and 5% CO2. CD200 and TNFAIP6 were selected from up-regulated gene groups and DHRS3, ASS1, and RIPK4 were selected from down-regulated gene groups. Commercial RNA interfering sets targeting CD200 and TNFAIP6 were purchased from Ribobio technology (Shenzhen, Guangdong, China), and negative control of siRNA was also transfected into MCF-7 cells to avoid the influence of transfection itself. Exogenous human full-length DHRS3, ASS1, and RIPK4 plasmid vectors were purchased from Sino biological Inc. (Beijing, China); and empty vector was also transfected into MCF7 cells as negative control. Lipofectamine 2000 (Invitrogen) was used for all transfection studies by following the manufacturer's protocol. Transfection efficiency was examined using qRT-PCR. Briefly, total RNA was extracted from cells using Trizol (Invitrogen, CA, USA) and further purified by RNAeasy mini kit plus DNase I treatment (Qiagen, Germany). The relative mRNA level for each gene was quantified by real-time RT-PCR with SYBR Green (Applied Biosystems, CA, USA), using GAPDH as a control.
Wound Healing Assay
The migration ability of cells was measured by wound-healing assay. Briefly, 2 × 104 cells were inoculated in 6-cm tissue culture dishes, cultured overnight, scratch was performed when cell growth fusion reached to 90%. Migration images were captured 24 and 48 h after scratching, the migration rate (width) of cells were calculated.
Transwell Migration
Corning™ BioCoat™ Matrigel™ Invasion Chamber with Corning™ Matrigel Matrix (Thermo Fisher Scientific, Waltham, MA, USA) was used for cell invasion assay. One hundred microliters of cells with concentration of 5 × 105/mL were added in the upper chamber with (8 μm pore), 600 μL medium containing 10% serum was added in the lower chamber. After cultured for 6 h, medium in the chamber were removed and unmigrated cells were swabbed. Cells was fixed by 4% polyformaldehyde for 10 min, stained with crystal violet for another 10 min. The filter membrane was photographed under the inverted microscope (200×) after sealed with the neutral gum. Cells were counted by Image Pro Plus Version 6 software, 3 wells in each group and 5 vision fields of each well were randomly selected to calculate the average number of cells.
Results
DEGs Analysis
A total of 749 DEGs were identified with the threshold of p-value < 0.05, consisting of 469 up-regulated genes and 280 down-regulated genes. Top 10 up-regulated DEGs were integrin binding sialoprotein (IBSP), matrix metallopeptidase 9 (MMP9), matrix metallopeptidase 13 (MMP13), osteomodulin (OMD), cathepsin K (CTSK), TNF alpha induced protein 6 (TNFAIP6), SATB homeobox 2 (SATB2), olfactomedin like 2B (OLFML2B), EGF like repeats and discoidin domains 3 (EDIL3), and procollagen C-endopeptidase enhancer (PCOLCE). The top 10 down-regulated DEGs were dehydrogenase/reductase 3 (DHRS3), argininosuccinate synthase 1 (ASS1), receptor interacting serine/threonine kinase 4 (RIPK4), vimentin (VIM), chondroitin sulfate proteoglycan 4 (CSPG4), ring finger protein 114 (RNF114), NDRG family member 2 (NDRG2), integrin subunit beta 8 (ITGB8), solute carrier family 6 member 6 (SLC6A6), and formyl peptide receptor 1 (FPR1). Detailed information about top 10 up- and down-regulated DEGs was listed in Table 2. The heat map of top 200 DEGs was shown in Figure 1. From the heat map, we can see that these DEGs were distinguished between the bone metastasis and other site metastasis.
Table 2. Top 10 up- and down-regulated DEGs.
Figure 1. The heat map of top 200 DEGs. Diagram presents the result of a two-way hierarchical clustering of top 200 DEGs and samples. The clustering is constructed using the complete-linkage method together with the Euclidean distance. Each row represents a DEG and each column, a sample. The DEG clustering tree is shown on the right. The color scale illustrates the relative level of DEG expression: red, below the reference channel; green, higher than the reference.
Functional and Pathway Enrichment Analyses of DEGs
To investigate the biological function of the identified DEGs in breast cancer bone metastasis, GO term and KEGG pathway enrichment analyses was performed. In GO term and KEGG pathway enrichment analyses, cell adhesion, angiogenesis, and signal transduction were the most significant enrichment in biological process; protein binding, calcium ion binding and collagen binding were the most remark enrichment in molecular function; extracellular matrix, extracellular region, and plasma membrane were the most significant enrichment in cellular component. ECM-receptor interaction, focal adhesion, and osteoclast differentiation were the most significant enrichment in KEGG signaling pathways. Top 25 GO terms of DEGs were presented in Figures 2A–C, respectively. The top 25 KEGG enrichment signal pathways of DEGs were shown in Figure 2D.
Figure 2. The functional enrichment of GO terms and KEGG pathways for differential mRNA. (A) Top 25 significant enrichment biological processes of DEGs. (B) Top 25 significant enrichment molecular functions of DEGs. (C) Top 25 significant enrichment cellular components of DEGs. (D) Top 25 significant enrichment KEGG signal pathways of DEGs.
PPI Network
To obtain the interaction between the proteins encoded by DEGs and other proteins, PPI network was explored and visualized by Cytoscape. PPI networks of proteins encoded by DEGs (FDR < 0.05) and interacting proteins encoded by DEGs (p < 0.01) were shown in Figure 3. As shown in Figure 3, the network consisted of 173 nodes and 158 edges. The rose red and green ellipse presented the up- and down-regulated proteins encoded by DEGs, respectively. The rectangle presented DEGs-encoded proteins with high degree (degree > 5). These high degree proteins were SMAD7 (degree = 10), TGFBR2 (degree = 9), VIM (degree = 8), FOS (degree = 8), PDGFRB (degree = 7), COL5A1 (degree = 6), ARRB2 (degree = 6), and ITGAV (degree = 6).
Figure 3. The PPI networks. The node denoted protein, and the edge denoted the interactions. The rose red and green ellipse presented the up- and down-regulated proteins encoded by genes, respectively. The rectangle presented high degree (degree > 5) proteins encoded by DEGs. The solid line and dotted line presented direct interaction and colocalization, respectively.
Establishment of TFs-Target Genes Regulatory Network
In order to investigate the TFs-target genes regulatory network for bone metastasis of breast cancer, we utilized TRANSFAC to obtain 53 TFs (from proteins encoded by all DEGs) regulating other DEGs (Figure 4). In the network, the hexagon and rectangle presented the TFs and DEGs, respectively. The rose red and green ellipse presented the up- and down-regulation, respectively. In the end, we obtained transcriptional regulatory networks consisting of 77 nodes and 77 edges. In this network, the top 7 TFs that covered the most DEGs were ETS1 (degree = 12), SPI1 (degree = 12), FOS (degree = 10), FLI1 (degree = 5), KLF4 (degree = 4), JUNB (degree = 4), NR3C1 (degree = 4).
Figure 4. The TFs-target genes networks. In the network, the hexagon and rectangle presented the TFs and DEGs, respectively. The rose red and green ellipse presented the up- and down-regulation, respectively.
QRT-PCR
In this study, eight candidate genes were selected from the top 10 up- or down-regulated DEGs, which were IBSP, MMP9, MMP13, TNFAIP6, CD200, DHRS3, ASS1, RIPK4, VIM, and PROM1 for validation of integrated analysis results (Table 2, Figures 5A,B). CD200 and PROM1 was also selected for qRT-PCR validation as cluster of differentiation (CD) markers from DEGs. The relative mRNA expression from adjacent non-cancer breast tissues was normalized as 1. As shown in Figure 5, results showed that the relative expression of MMP9, MMP13, TNFAIP6, and CD200, were significantly up-regulated (P < 0.05), while DHRS3, ASS1, and VIM were significantly down-regulated in the bone metastasis compared with lung and liver metastasis (P < 0.05). It is noted that, although no statistical significance was found for IBSP, RIPK4, and PROM1, their expression trends were similar with bioinformatics data. All the mRNA level of tested genes from adjacent non-tumor breast cancer was significantly different from metastasis tumor tissues, except RIPK4 (Figure 5).
Figure 5. Identification of differential mRNA by QRT-PCR using metastatic bone, liver, and lung tissues. (A) Expression of IBSP, MMP9, MMP13, TNFAIP6, and CD200 in the tissue of patients with breast cancer bone metastasis, lung metastasis and liver metastasis by qRT-PCR. The x axis and y axis presented gene name and relative expression, respectively. (B) Expression of DHRS3, ASS1, RIPK4, VIM, and PROM1 in the tissue of patients with breast cancer bone metastasis, lung metastasis, and liver metastasis by qRT-PCR. The x axis and y axis presented gene name and relative expression, respectively. #P < 0.05, adjacent non-tumor breast tissue compared with three site metastases; *p < 0.05, bone metastasis compare with lung and liver metastasis.
ROC Curve Analysis
In order to access the discriminatory ability of the above 10 genes used for qRT-PCR validation in breast cancer bone metastasis, ROC curve analyses were conducted and AUC were calculated. As shown in Figure 6, the AUC of these genes was more than 0.7 except PROM1. For breast cancer bone metastasis diagnosis, the sensitivity and specificity of these genes were very applicable.
Figure 6. The ROC curve analyses of IBSP, MMP9, MMP13, TNFAIP6, CD200, DHRS3, ASS1, RIPK4, VIM, and PROM1 in breast cancer bone metastasis. The x axis and y axis presented specificity and sensitivity, respectively.
Migration and Invasion Ability of Breast Cancer Cell Regulated by Gene Knocking Down or Over-expression
To evaluate the impact of these predicted genes on cell migration and invasion, the wound healing assay and matrigel invasion assay were employed. As shown in Figure 7A, by RT-PCR, we confirmed that RNA interference significantly reduced the mRNA level of CD200 and TNFAIP6, and meanwhile exogenous expressions of DHRS3, ASS1, and RIPK4 were remarkably higher than empty-vector control cells. We found that overexpression of DHRS3, ASS1, and RIPK4 inhibited MCF-7 cell migration, as well as knock-down CD200 and TNFAIP6 (Figure 7B). Consistent with this finding, matrigel invasion assay showed that DHRS3, ASS1, and RIPK4 overexpression significantly reduced invasion capacity of MCF7, which is similar as knocking down CD200 and TNFAIP6 (Figure 7C). These observations suggested that these predicted gene played important roles in regulating migration and invasive potential of breast cancer cells.
Figure 7. Wound healing assay and transwell invasion assay. (A) RT-PCR analysis of selected gene expression in MCF-7 cells after knocking down and knocking-up treatment. (B) wound closure was shown at 48 h post scratch (200×). Cell migration was assessed by recover of the scratch. The area of the wound was measured at the two time points in every group, and % reduction of initial scratch area was compared with control; and the results were expressed as percentage relative to the corresponding blank control. (C) Transwell invasion assay. MCF-7 cells penetrating the membrane were fixed and 0.1% crystal violet stained after 24 h as described in experimental procedures, and cell number penetrating the membrane were counted for each group. *P < 0.05, **P < 0.01. Data were the means of three measurements and the bars represented SD of the mean.
Discussions
Among numerous malignancies, metastatic breast cancer is the second leading cause of woman death (14). Clinically, breast cancer can metastasize to various sites such as brain, liver, lung, and bone, and bone is the most common metastasis site (15). Therefore, maintaining bone health is an important clinical challenge in breast cancer patients. In this study, we found several DEGs and signaling pathways related to breast cancers bone metastasis, which may be helpful in understanding the bone metastasis mechanism of breast cancer.
In this study, a total of 749 DEGs were identified, consisting of 469 up-regulated genes and 280 down-regulated genes and eight DEGs were selected from the top 10 up- or down-regulated DEGs, such as IBSP, MMP9, MMP13, TNFAIP6, DHRS3, ASS1, RIPK4, and VIM for the following qRT-PCR validation. CD200 and PROM1 was also selected for qRT-PCR validation as cluster of differentiation markers. Although there is no statistical difference except CD200 and DHRS3, the tendency of expression was generally consist with the bioinformatics analysis results. In addition, their diagnostic ability was also evaluated, and it appeared that IBSP, MMP9, MMP13, TNFAIP6, CD200, DHRS3, ASS1, RIPK4, and VIM may be considered as specific prognostic markers for bone metastasis to identify primary breast cancers that have potential to metastasize to bone tissue.
There are some preclinical evidences of factors associated with breast cancer cells homing to bone tissue. Awolaran et al. pointed out that 15 proteins expressed by breast cancer cells mediated breast cancer bone metastasis: ICAM-1, CDH11, osteoactivin, bone sialoprotein, CCN3, IL-11, etc. (16). In a study by Mercatali et al., SPARC, IBSP, MMP9 was found to be over expressed in patient-derived xenograft of breast cancer bone metastasis compared to cell lines of MCF7 and MDA-MB-231 (17).
In our study, we detected some published potential biomarker of breast cancer bone metastasis, or their function has been reported in breast cancer bone metastasis such as SPARC, IBSP, MMP9, MMP13. SPARC was found as a DEG in breast cancer bone metastasis, while wasn't in the list of the up- or down-regulated DEGs. SPARC, also known as osteonectin, is involved in the chemoattraction of both breast and prostate cancer, and has been detected as a biomarker for bone metastasis (18).
IBSP, an osteoblast differentiation marker, is associated with osteoblastogenesis, bone formation and tumor-induced osteolysis (19, 20). It is reported that IBSP is up-regulated in luminal B2 breast cancer and breast cancer with bone metastasis (17), while the expression level of IBSP was relatively low in older patients with breast cancer compared with those young patients (21, 22). Furthermore, the expression of IBSP is associated with development of metastasis and poor prognosis in breast cancer (23, 24).
Matrix metallopeptidase 9 (MMP9) and matrix metallopeptidase 13 (MMP13), members of the peptidase M10 family of matrix metalloproteinases (MMPs), is normally expressed by resident bone cells primarily by osteoclasts and osteoblasts (25). It was no surprise to find high levels of it in the bone metastasis. MMP9 and MMP13 is up-regulated in breast cancer tissue, and were involved in bone metastasis of tumor progression (26–28). It is reported that up-regulated MMP9 (29) and MMP13 is both associated with poorer overall survival and can serve as a prognosis indicator of breast cancer (30, 31).
Besides some published potential biomarkers of breast cancer bone metastasis were identified in our study, we also detected some novel biomarkers or their relationship with breast cancer bone metastasis has never been reported, such as TNFAIP6, DHRS3, ASS1, RIPK4, VIM, CD200. However, their function in breast cancer has been investigated except RIPK4.
Cancer cells can create an inflammatory microenvironment to enhance the tumor metastatic process (32). It is reported that excessive and chronic activation of the immune system is involved in breast cancer metastasis (33). TNF alpha induced protein 6 (TNFAIP6) encodes inflammatory response factors and regulates anti-inflammation response and immunosuppression (34, 35). It is indicated that TNFAIP6 is a prognostic-related cytokine in breast cancers (36).
Dehydrogenase/reductase 3 (DHRS3) has been found to be involved in the tumor suppressive pathway and constitutively expressed in breast cancer cell lines (37, 38). It is found that DHRS3 is down-regulated in breast cancer and breast cancer brain metastases (39, 40). CD200 molecule (CD200) serum levels are significantly higher in the early and the advanced stage breast cancer patients (41). Gorczynski et al. found that over expression of CD200 increased breast cancer lymph node metastasis (42). In addition, it has been demonstrated that the enhancement of the CD200-CD200 receptor interaction inhibits the metastasis of breast cancer (43). It is reported that CD200 analogs may have therapeutic potential in treating aggressive breast carcinoma (33). ASS1 (44) and VIM (45) are both prognostic indicators in breast cancer patients. The research of TNFAIP6, DHRS3, ASS1, RIPK4, VIM, and CD200 in breast cancer may indicate that those genes may play a crucial role in the development of breast cancer bone metastasis.
For RIPK4, their relationship to breast cancer still remains unreported. Receptor interacting serine/threonine kinase 4 (RIPK4) is involved in a number of signaling pathways and plays an important role in regulating cells proliferation, cell differentiation and cell apoptosis (46). More and more evidences have shown that RIPK4 is up-regulated in tumor tissues and promotes the occurrence and progression of cervical cancer and ovarian cancer (47, 48).
According to the KEGG pathway enrichment analysis, we found two important signaling pathways including osteoclast differentiation and rheumatoid arthritis in the bone metastasis of breast cancer. Bone metastasis of breast cancer induces severe osteolysis with increased bone resorption (25). Interestingly, osteoclast differentiation is involved in the process of bone resorption. It is reported that bisphosphonate, an inhibitor of osteoclastic bone resorption, is a valuable agent in treating bone metastasis of breast cancer patients (49). Previous epidemiologic study indicated that breast cancer patients with rheumatoid arthritis had poor prognoses and higher mortality compared with breast cancer patients without rheumatoid arthritis (50). In identified DEGs, we found two high degree gene in both PPI and TF-target gene networks, which were fos proto-oncogene, AP-1 transcription factor subunit (FOS) and one of top 10 up-regulated gene cathepsin K (CTSK) was involved in osteoclast differentiation and rheumatoid arthritis. FOS is an important osteoclasts' differentiation maker involved in osteoclastogenesis (51). Sotiriou et al. found that FOS was a confident driver gene for breast cancer metastasis (52). CTSK is an osteoclast phenotypic marker. It is found that CTSK is secreted into the isolated environment, which leads to the bone resorption (53). Le Gall et al. found that CTSK inhibitors have been used in the preclinical studies, which have shown the effectiveness of CTSK in the reduction of breast carcinomas bone metastases in experimental animal model (54).
To further confirm the results of bioinformatics prediction, the wound healing assay and matrigel invasion assay were performed. Total five genes were selected, two are up-regulated in breast cancer and three are down-regulated in breast cancer. As shown in Figure 7, altering expression of these genes in breast cancer cells significantly influencing the migration and invasion ability.
In summary, we found several DEGs including IBSP, MMP9, MMP13, TNFAIP6, DHRS3, ASS1, RIPK4, VIM, CD200, and PROM1 may play important roles in the process of bone metastasis from breast cancer. Among which, IBSP, MMP9, TNFAIP6, DHRS3, RIPK4, and CD200 had a diagnose value for patients with breast cancer bone metastasis. Additionally, osteoclast differentiation and rheumatoid arthritis signaling pathways and related DEGs including FOS and CTSK may be involved in the development of bone metastasis from breast cancer. However, there are limitations to our study. Firstly, sample size in the QRT-PCR was small and a larger-scale study is further needed. Secondly, the biological function of identified DEGs wasn't investigated. Experiments in animal model and cell culture are further needed to explore underlying functional mechanism of breast cancer bone metastasis.
Author Contributions
YZ is in charge of data collection and data mining. WH is in charge of molecular biological experiments and revising manuscript. SZ is in charge of data analysis, molecular biological experiments, and manuscript draft.
Funding
This study was supported by CAMS Innovation Fund for Medical Sciences (CIFMS, no. 2016-I2M-3-011) and the 13th Five-Year Plan Major Scientific and Technological Projects-Major New Drug Creation (No.2018ZX09711001).
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
References
1. Adams S, Munz B. RIP4 is a target of multiple signal transduction pathways in keratinocytes: implications for epidermal differentiation and cutaneous wound repair. Exp Cell Res. (2010) 316:126–37. doi: 10.1016/j.yexcr.2009.10.006
2. Zheng R, Zeng H, Zhang S, Chen T, Chen W. National estimates of cancer prevalence in China, 2011. Cancer Lett. (2016) 370:33–8. doi: 10.1016/j.canlet.2015.10.003
3. Graham N, Qian BZ. Mesenchymal stromal cells: emerging roles in bone metastasis. Int J Mol Sci. (2018) 19:E1121. doi: 10.3390/ijms19041121
4. Coleman RE. The role of bisphosphonates in breast cancer. Breast. (2004) 13:19–28. doi: 10.1016/j.breast.2004.09.004
5. Braun S, Vogl FD, Naume B, Janni W, Osborne MP, Coombes RC, et al. A pooled analysis of bone marrow micrometastasis in breast cancer. N Engl J Med. (2005) 353:793–802. doi: 10.1056/NEJMoa050434
6. Mishra A, Shiozawa Y, Pienta KJ, Taichman RS. Homing of cancer cells to the bone. Cancer Microenviron. (2011) 4:221. doi: 10.1007/s12307-011-0083-6
7. Roodman GD. Mechanisms of bone metastasis. N Engl J Med. (2004) 350:1655–64. doi: 10.1056/NEJMra030831
8. Colleoni M, O'Neill A, Goldhirsch A, Gelber RD, Bonetti M, Thurlimann B, et al. Identifying breast cancer patients at high risk for bone metastases. J Clin Oncol. (2000) 18:3925–35. doi: 10.1200/JCO.2000.18.23.3925
9. Wei B, Wang J, Bourne P, Yang Q, Hicks D, Bu H, et al. Bone metastasis is strongly associated with estrogen receptor-positive/progesterone receptor-negative breast carcinomas. Hum Pathol. (2008) 39:1809–15. doi: 10.1016/j.humpath.2008.05.010
10. Mundy GR. Mechanisms of bone metastasis. Cancer. (1997) 80:1546–56. doi: 10.1002/(SICI)1097-0142(19971015)80:8+<1546::AID-CNCR4>3.0.CO;2-I
11. Zubeda S, Kaipa PR, Shaik NA, Mohiuddin MK, Vaidya S, Pavani B, et al. Her-2/neu status: a neglected marker of prognostication and management of breast cancer patients in India. Asian Pacific J Cancer Prevent. (2013) 14:2231–5. doi: 10.7314/APJCP.2013.14.4.2231
12. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Stat Soc. (1995) 57:289–300.
13. Reiner-Benaim A. FDR control by the BH procedure for two-sided correlated tests with implications to gene expression data analysis. Biometr J Biometr Zeitschrift. (2007) 49:107–26. doi: 10.1002/bimj.200510313
14. DeSantis CE, Fedewa SA, Goding Sauer A, Kramer JL, Smith RA, Jemal A. Breast cancer statistics, 2015: convergence of incidence rates between black and white women. Cancer J Clinic. (2016) 66:31–42. doi: 10.3322/caac.21320
15. Weigelt B, Peterse JL, van 't Veer LJ. Breast cancer metastasis: markers and models. Nat Rev Cancer. (2005) 5:591–602. doi: 10.1038/nrc1670
16. Awolaran O, Brooks SA, Lavender V. Breast cancer osteomimicry and its role in bone specific metastasis; an integrative, systematic review of preclinical evidence. Breast. (2016) 30:156–71. doi: 10.1016/j.breast.2016.09.017
17. Mercatali L, La Manna F, Groenewoud A, Casadei R, Recine F, Miserocchi G, et al. Development of a patient-derived xenograft (PDX) of breast cancer bone metastasis in a zebrafish model. Int J Mol Sci. (2016) 17:E1375. doi: 10.3390/ijms17081375
18. Li S, Zou D, Li C, Meng H, Sui W, Feng S, et al. Targeting stem cell niche can protect hematopoietic stem cells from chemotherapy and G-CSF treatment. Stem Cell Res Ther. (2015) 6:175. doi: 10.1186/s13287-015-0164-4
19. Javed A, Bae JS, Afzal F, Gutierrez S, Pratap J, Zaidi SK, et al. Structural coupling of Smad and Runx2 for execution of the BMP2 osteogenic signal. J Biol Chem. (2008) 283:8412–22. doi: 10.1074/jbc.M705578200
20. Nannuru KC, Futakuchi M, Varney ML, Vincent TM, Marcusson EG, Singh RK. Matrix metalloproteinase. (MMP)-13 regulates mammary tumor-induced osteolysis by activating MMP9 and transforming growth factor-beta signaling at the tumor-bone interface. Cancer Res. (2010) 70:3494–04. doi: 10.1158/0008-5472.CAN-09-3251
21. Brouwers B, Fumagalli D, Brohee S, Hatse S, Govaere O, Floris G, et al. The footprint of the ageing stroma in older patients with breast cancer. Breast Cancer Res. (2017) 19:78. doi: 10.1186/s13058-017-0871-0
22. Wang J, Du Q, Li C. Bioinformatics analysis of gene expression profiles to identify causal genes in luminal B2 breast cancer. Oncol Lett. (2017) 14:7880–8. doi: 10.3892/ol.2017.7256
23. Bellahcene A, Kroll M, Liebens F, Castronovo V. Bone sialoprotein expression in primary human breast cancer is associated with bone metastases development. J Bone Mineral Res. (1996) 11:665–70. doi: 10.1002/jbmr.5650110514
24. Bellahcene A, Menard S, Bufalino R, Moreau L, Castronovo V. Expression of bone sialoprotein in primary human breast cancer is associated with poor survival. Int J Cancer. (1996) 69:350–3. doi: 10.1002/(SICI)1097-0215(19960822)69:4<350::AID-IJC19>3.0.CO;2-9
25. Ohshiba T, Miyaura C, Inada M, Ito A. Role of RANKL-induced osteoclast formation and MMP-dependent matrix degradation in bone destruction by breast cancer metastasis. Br J Cancer. (2003) 88:1318–26. doi: 10.1038/sj.bjc.6600858
26. Cackowski FC, Anderson JL, Patrene KD, Choksi RJ, Shapiro SD, Windle JJ, et al. Osteoclasts are important for bone angiogenesis. Blood 115:140. doi: 10.1182/blood-2009-08-237628
27. Kohrmann A, Kammerer U, Kapp M, Dietl J, Anacker J. Expression of matrix metalloproteinases (MMPs) in primary human breast cancer and breast cancer cell lines: new findings and review of the literature. BMC Cancer. (2009) 9:188. doi: 10.1186/1471-2407-9-188
28. Lynch CC, Hikosaka A, Acuff HB, Martin MD, Kawai N, Singh RK, et al. MMP-7 promotes prostate cancer-induced osteolysis via the solubilization of RANKL. Cancer Cell. (2005) 7:485–96. doi: 10.1016/j.ccr.2005.04.013
29. Pellikainen JM, Ropponen KM, Kataja VV, Kellokoski JK, Eskelinen MJ, Kosma VM. Expression of matrix metalloproteinase. (MMP)-2 and MMP-9 in breast cancer with a special reference to activator protein-2, HER2, and prognosis. Clin Cancer Res. (2004) 10:7621–8. doi: 10.1158/1078-0432.CCR-04-1061
30. Paek AR, Mun JY, Hong KM, Lee J, Hong DW, You HJ. Zinc finger protein 143 expression is closely related to tumor malignancy via regulating cell motility in breast cancer. BMB Rep. (2017) 50:621–7. doi: 10.5483/BMBRep.2017.50.12.177
31. Zhang B, Cao X, Liu Y, Cao W, Zhang F, Zhang S, et al. Tumor-derived matrix metalloproteinase-13. (MMP-13). Correlates with poor prognoses of invasive breast cancer. BMC Cancer. (2008) 8:83. doi: 10.1186/1471-2407-8-83
32. Hanahan D, Weinberg RA. Hallmarks of cancer: the next generation. Cell. (2011) 144:646–74. doi: 10.1016/j.cell.2011.02.013
33. Erin N, Tanriover G, Curry A, Akman M, Duymus O, Gorczynski R. CD200fc enhances anti-tumoral immune response and inhibits visceral metastasis of breast carcinoma. Oncotarget. (2018) 9:19147–58. doi: 10.18632/oncotarget.24931
34. Choi H, Lee RH, Bazhanov N, Oh JY, Prockop DJ. Anti-inflammatory protein TSG-6 secreted by activated MSCs attenuates zymosan-induced mouse peritonitis by decreasing TLR2/NF-kappaB signaling in resident macrophages. Blood. (2011) 118:330–8. doi: 10.1182/blood-2010-12-327353
35. Lee RH, Pulin AA, Seo MJ, Kota DJ, Ylostalo J, Larson BL, et al. Intravenous hMSCs improve myocardial infarction in mice because cells embolized in lung are activated to secrete the anti-inflammatory protein TSG-6. Cell Stem Cell. (2009) 5:54–63. doi: 10.1016/j.stem.2009.05.003
36. Wong HS, Chang CM, Liu X, Huang WC, Chang WC. Characterization of cytokinome landscape for clinical responses in human cancers. Oncoimmunology. (2016) 5:e1214789. doi: 10.1080/2162402X.2016.1214789
37. Cerignoli F, Guo X, Cardinali B, Rinaldi C, Casaletto J, Frati L, et al. retSDR1, a short-chain retinol dehydrogenase/reductase, is retinoic acid-inducible and frequently deleted in human neuroblastoma cell lines. Cancer Res. (2002) 62:1196–204.
38. Haeseleer F, Huang J, Lebioda L, Saari JC, Palczewski K. Molecular characterization of a novel short-chain dehydrogenase/reductase that reduces all-trans-retinal. J Biol Chem. 273:21790. doi: 10.1074/jbc.273.34.21790
39. Lee CH, Kuo WH, Lin CC, Oyang YJ, Huang HC, Juan HF. MicroRNA-regulated protein-protein interaction networks and their functions in breast cancer. Int J Mol Sci. (2013) 14:11560–606. doi: 10.3390/ijms140611560
40. Schulten HJ, Bangash M, Karim S, Dallol A, Hussein D, Merdad A, et al. Comprehensive molecular biomarker identification in breast cancer brain metastases. J Transl Med. (2017) 15:269. doi: 10.1186/s12967-017-1370-x
41. Celik B, Yalcin AD, Genc GE, Bulut T, Kuloglu Genc S, Gumuslu S, et al. CXCL8, IL-1beta and sCD200 are pro-inflammatory cytokines and their levels increase in the circulation of breast carcinoma patients. Biomed Rep. (2016) 5:259–63. doi: 10.3892/br.2016.709
42. Gorczynski RM, Clark DA, Erin N, Khatri I. Role of CD200 expression in regulation of metastasis of EMT6 tumor cells in mice. Breast Cancer Res Treat. (2011) 130:49–60. doi: 10.1007/s10549-010-1259-3
43. Erin N, Podnos A, Tanriover G, Duymus O, Cote E, Khatri I, et al. Bidirectional effect of CD200 on breast cancer development and metastasis, with ultimate outcome determined by tumor aggressiveness and a cancer-induced inflammatory response. Oncogene. (2015) 34:3860–70. doi: 10.1038/onc.2014.317
44. Qiu F, Chen YR, Liu X, Chu CY, Shen LJ, Xu J, et al. Arginine starvation impairs mitochondrial respiratory function in ASS1-deficient breast cancer cells. Sci Signal. (2014) 7:ra31. doi: 10.1126/scisignal.2004761
45. Xiao B, Chen L, Ke Y, Hang J, Cao L, Zhang R, et al. Identification of methylation sites and signature genes with prognostic value for luminal breast cancer. BMC Cancer. (2018) 18:405. doi: 10.1186/s12885-018-4314-9
46. Sandholm J, Selander KS. Toll-like receptor 9 in breast cancer. Front Immunol. (2014) 5:330. doi: 10.3389/fimmu.2014.00330
47. Huang X, McGann JC, Liu BY, Hannoush RN, Lill JR, Pham V, et al. Phosphorylation of dishevelled by protein kinase RIPK4 regulates Wnt signaling. Science. (2013) 339:1441–5. doi: 10.1126/science.1232253
48. Liu DQ, Li FF, Zhang JB, Zhou TJ, Xue WQ, Zheng XH, et al. Increased RIPK4 expression is associated with progression and poor prognosis in cervical squamous cell carcinoma patients. Sci Rep. (2015) 5:11955. doi: 10.1038/srep11955
49. Body JJ, Mancini I. Bisphosphonates for cancer patients: why, how, and when? Support Care Cancer. (2002) 10:399–407. doi: 10.1007/s005200100292
50. Ji J, Liu X, Sundquist K, Sundquist J. Survival of cancer in patients with rheumatoid arthritis: a follow-up study in Sweden of patients hospitalized with rheumatoid arthritis 1 year before diagnosis of cancer. Rheumatology. (2011) 50:1513–8. doi: 10.1093/rheumatology/ker143
51. Joung YH, Darvin P, Kang DY, Sp N, Byun HJ, Lee CH, et al. Methylsulfonylmethane Inhibits RANKL-Induced Osteoclastogenesis in BMMs by Suppressing NF-kappaB and STAT3 Activities. PLoS ONE. (2016) 11:e0159891. doi: 10.1371/journal.pone.0159891
52. Sotiriou C, Wirapati P, Loi S, Harris A, Fox S, Smeds J, et al. Gene expression profiling in breast cancer: understanding the molecular basis of histologic grade to improve prognosis. J Natl Cancer Inst. (2006) 98:262–72. doi: 10.1093/jnci/djj052
53. Edwards CM, Mundy GR. Eph receptors and ephrin signaling pathways: a role in bone homeostasis. Int J Med Sci. (2008) 5:263–72. doi: 10.7150/ijms.5.263
Keywords: breast cancer, bone metastasis, diagnose, mRNA expression profile, protein-protein interaction, transcription factors-target gene
Citation: Zhang Y, He W and Zhang S (2019) Seeking for Correlative Genes and Signaling Pathways With Bone Metastasis From Breast Cancer by Integrated Analysis. Front. Oncol. 9:138. doi: 10.3389/fonc.2019.00138
Received: 25 September 2018; Accepted: 18 February 2019;
Published: 13 March 2019.
Edited by:
Tao Liu, University of New South Wales, AustraliaReviewed by:
Prashant Trikha, Nationwide Children's Hospital, United StatesSabarish Ramachandran, Texas Tech University Health Sciences Center, United States
Copyright © 2019 Zhang, He and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Sen Zhang, zhangs@imm.ac.cn
†These authors have contributed equally to this work