This article is mentioned in:

Introduction

Breast cancer is a type of cancer that originates from breast tissue and most commonly from the inner lining of the milk ducts or from the lobules that supply the ducts (1). Breast cancer affects ~1.2 million women worldwide and accounts for ~50,000 mortalities every year (2). Despite major advances in surgical and nonsurgical management of the disease, breast cancer metastasis remains a significant clinical challenge affecting numerous of patients (3). The prognosis and survival rates for breast cancer are highly variable, and depend on the cancer type, treatment strategy, stage of the disease and geographical location of the patient (4).

Microarray technology, which may be used to simultaneously interrogate 10,000–40,000 genes, has provided new insight into the molecular classification of different cancer types (5). Many genes involved in the regulation of cell cycle, invasion, metastasis and angiogenesis have been indicated to be prognostic biomarkers based on microarray analyses. Furthermore, numerous researchers have proposed that the phenotypic diversity of breast tumors may be accompanied by a corresponding diversity in gene expression patterns (6). Therefore, the systematic investigation of gene expression patterns in human breast tumors may aid in understanding the pathogenesis of this disease (7).

The present study aimed to investigate the differences between breast cancer cells and normal cells, and examine the possible underlying mechanisms of breast cancer. Biological microarray analysis was used to analyze the expression profile of breast cancer and normal cells, and identify the differentially expressed genes (DEGs). In addition, altered metabolic pathways in breast cancer were identified using a bioinformatics approach, while the target sites of potential transcription factors and miRNAs were screened. Furthermore, the current study aimed to improve the understanding of the occurrence and development of breast cancer and facilitate the discovery of potential novel biomarkers for its treatment.

Materials and methods

Gene expression microarray

In order to investigate the alterations in breast cancer cells compared with normal cells, DEGs were screened at the gene level and the possible mechanisms were examined. The gene expression series, GSE9574 (8), was downloaded from the Gene Expression Omnibus (GEO; http://www.ncbi.nlm.nih.gov/geo/) of the National Center for Biotechnology Information (Bethesda, MD, USA), based on the GPL96 [HG-U133] platform data (Affymetrix Human Genome U133 Array), and included 14 breast cancer and 15 normal samples.

Identification of DEGs

Microarray data were analyzed using R software v.2.13.0 (9) and further processed using Geoquery (10) and Limma (11) packages. Geoquery is used to rapidly obtain gene expression profiles from the GEO database (10), while the Limma package is the most popular method for the analysis of DEGs (11,12). The preprocessed expression data were obtained using Geoquery, and subjected to a log2 transformation. The breast cancer and normal samples were compared using Limma in order to identify the DEGs between the two tissue types. Gene P-values were determined by R software using the Student's t-test and P<0.001 in a Bayesian model was considered to indicate a DEG (12).

Gene Ontology (GO) analysis of DEGs

In order to assess the changes in DEGs occurring at the cellular level and the functional clustering of DEGs, the GO Enrichment Analysis Software Toolkit (GOEAST) (13) in the GO database (14) was used. Hypergeometric algorithms were used for statistical analysis, and terms associated with biological process and molecular function were enriched.

Bio-pathway analysis of DEGs

In order to detect the changes in DEGs at the molecular level, all the metabolic and nonmetabolic pathways were obtained from the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. Pathway enrichment of DEGs was assessed using the Gene Set Analysis Toolkit v2 (15,16), and the number of genes was counted for each term. A gene number of ≥2 and P<0.05 were used as the cut-off values.

Protein-protein interaction (PPI) network construction

PPI data were integrated and verified using data obtained from the following databases (all accessed on November 11, 2012): Human Protein Reference Database (http://www.hprd.org/); Biological General Repository for Interaction Datasets (http://thebiogrid.org/); Biomolecular Object Network Databank (http://bond.unleashedinformatics.com/); Database of Interacting Proteins (http://dip.doe-mbi.ucla.edu/); IntAct (http://www.ebi.ac.uk/intact/); Molecular INTeraction database (http://mint.bio.uniroma2.it/mint/welcome.do); and Reactome (http://www.reactome.org/). Subsequently, a PPI network was constructed based on the identified DEGs. A hypergeometric algorithm was used, and P<0.05 was considered to indicate statistically significant differences.

Screening the target sites of potential transcription factors and miRNAs

Based on the gene annotation data from the Molecular Signatures Database (http://www.broadinstitute.org/gsea/msigdb/index.jsp; accessed November 11, 2012), the abundance of the gene sets were analyzed. In addition, hypergeometric algorithms and false discovery rate (FDR) correction were performed using the Benjamini & Hochberg method (17). FDR<0.05 was selected as the statistical significance threshold to indicate target sites of potential transcription factors and miRNAs.

Results

Identification of DEGs in breast cancer

Using P<0.001 as the statistical significance threshold, a total of 123 probes were identified, which presented altered expression levels in breast cancer when compared with normal tissues, and involved 106 DEGs (Table I).

Table I. Differentially expressed genes in breast cancer tissue as compared with normal tissue.

Table I.

Differentially expressed genes in breast cancer tissue as compared with normal tissue.

Gene symbol	P-value
PTP4A1	7.48×10−9
IER2	1.10×10−8
FOSB	1.23×10−8
NR4A3	4.97×10−8
ATF3	1.36×10−7
BTG2	1.61×10−7
NOL12	2.01×10−7
FOS	2.21×10−7
TACSTD2	5.92×10−7
H3F3B	9.46×10−7
JUN	1.00×10−6
EIF1	1.17×10−6
NR4A2	1.72×10−6
DUSP1	5.29×10−6
EIF5	7.45×10−6
TGFB2	1.13×10−5
APOH	1.31×10−5
ZFP36	2.30×10−5
GPR183	2.34×10−5
JUND	2.51×10−5
EXOC7	2.71×10−5
MCL1	2.73×10−5
REXO4	3.17×10−5
CLDN1	3.69×10−5
HEY2	4.63×10−5
KLF6	5.98×10−5
CYLD	6.32×10−5
SH3BP2	6.76×10−5
SERPINE1	7.88×10−5
HIST2H2BE	8.21×10−5
PRR5-ARHGAP8 / ARHGAP8	8.36×10−5
EGR1	9.68×10−5
KLF4	1.09×10−4
IL5RA	1.24×10−4
PMAIP1	1.33×10−4
FANCG	1.35×10−4
NOS1AP	1.65×10−4
C1orf50	1.69×10−4
MUM1	1.69×10−4
STK17B	1.72×10−4
ILF3	1.91×10−4
CLNS1A	2.01×10−4
BHLHE40	2.21×10−4
ATP2B2	2.30×10−4
TSC22D2	2.31×10−4
LAGE3	2.35×10−4
ZNHIT1	2.25×10−4
CXCL2	2.60×10−4
UNC119B	2.62×10−4
C9orf3	2.78×10−4
GBAP1	2.85×10−4
TTC38	2.86×10−4
IL8	2.86×10−4
SAFB2	2.98×10−4
SPATS2L	3.01×10−4
MST1P9	3.12×10−4
CCL2	3.16×10−4
AKR7A2	3.34×10−4
LSR	3.35×10−4
DNAJB4	3.38×10−4
KLF11	3.42×10−4
TP53TG1	3.60×10−4
GNAS	3.79×10−4
CD69	3.90×10−4
NR4A1	3.92×10−4
TIMM23 / TIMM23B	4.22×10−4
JUNB	4.30×10−4
FGF20	4.49×10−4
HIC2	4.51×10−4
GP6	4.53×10−4
SIK1	4.77×10−4
PGS1	4.77×10−4
SIDT2	4.89×10−4
IL33	5.04×10−4
FLRT2	5.05×10−4
TRAIP	5.10×10−4
PLXND1	5.21×10−4
YWHAZ	5.55×10−4
LINC00094	5.60×10−4
ZNF451	5.90×10−4
ZNF232	6.16×10−4
DUSP2	6.24×10−4
EBLN2	6.31×10−4
XPNPEP1	6.40×10−4
FAM46C	6.61×10−4
PDE4B	6.77×10−4
RBM14	6.79×10−4
RGS1	6.83×10−4
SLC33A1	7.28×10−4
SLC35E2	7.34×10−4
CNDP2	7.43×10−4
GTF2H1	7.64×10−4
CARKD	7.77×10−4
CCNL1	7.78×10−4
CDC42	7.84×10−4
TAF1D	8.08×10−4
CLGN	8.22×10−4
SGSM2	8.23×10−4
KIF1A	8.48×10−4
NRAS	8.87×10−4
KLF10	8.96×10−4
ORC5	9.17×10−4
CCDC170	9.41×10−4
C14orf105	9.54×10−4
MORC2	9.74×10−4
EIF4A1	9.99×10−4

GO cluster of DEGs

The molecular functions enriched in the identified DEGs included nucleic acid binding transcription factor activity, sequence-specific DNA binding transcription factor activity and double-stranded DNA binding (Fig. 1). In addition, the biological processes enriched are shown in Fig. 2, and include positive regulation of biological process, positive regulation of cellular process, cellular response to organic substance and positive regulation of transcription from RNA polymerase II promoter. The GO clustering results provide a preliminary description of the potential functions of the DEGs and their effects on cells.

Figure 1. Molecular function clustering of differentially expressed genes. Yellow boxes are significantly enriched Gene Ontology terms (false discovery rate <0.051). Brighter color indicates stronger statistical significance. GO, Gene Ontology.

Figure 2. Biological process clustering of differentially expressed genes. Colored boxes are significantly enriched GO terms (false discovery rate <0.051). A brighter color indicates stronger statistical significance. GO, Gene Ontology.

Bio-pathways altered in breast cancer

To further investigate changes of the biological pathways within cancer cells in detail, a KEGG pathway enrichment analysis of the identified DEGs was performed. KEGG clustering results indicated that a number of bio-pathways were altered in breast cancer cells, primarily signaling and disease-associated pathways (Table II). The alteration of the RNA transport pathway was consistent with the GO clustering results, suggesting that gene expression in breast cancer cells differs from that in normal cells. In addition, signaling pathways in the cell surface were altered, including the nucleotide-binding oligomerization domain (NOD)-like receptor signaling pathway, epithelial cell signaling in Helicobacter pylori infection, chemokine signaling pathway, T cell receptor signaling pathway and B cell receptor signaling pathway. The latter three are involved in immune reactions, while the NOD-like receptor signaling pathway is closely associated with cell differentiation and development. Alterations in other metabolic pathways were also indicated, including the mismatch repair pathway, proteasome, metabolic pathways, glutathione metabolism and glycolysis/gluconeogenesis. Notably, KEGG clustering results indicated that a number of other disease pathways may also be altered, including Chagas disease (American trypanosomiasis), osteoclast differentiation, rheumatoid arthritis, renal cell carcinoma, pathogenic Escherichia coli infection, malaria, colorectal cancer, pathways in cancer and leishmaniasis.

Table II. Bio-pathways altered in breast cancer.

Table II.

Bio-pathways altered in breast cancer.

KEGG pathway	P-value
MAPK signaling pathway	0.0006
Chagas disease (American trypanosomiasis)	0.0006
Osteoclast differentiation	0.0009
Rheumatoid arthritis	0.0048
Renal cell carcinoma	0.0202
NOD-like receptor signaling pathway	0.0432
GnRH signaling pathway	0.0432
Pathogenic Escherichia coli infection	0.0432
Malaria	0.0432
T cell receptor signaling pathway	0.0432
Chemokine signaling pathway	0.0441
Colorectal cancer	0.0445
Epithelial cell signaling in Helicobacter pylori infection	0.0445
Neurotrophin signaling pathway	0.0445
Pathways in cancer	0.0445
Leishmaniasis	0.0475
RNA transport	0.0475
B-cell receptor signaling pathway	0.0475

[i] KEGG, Kyoto Encyclopedia of Genes and Genomes; NOD, nucleotide-binding oligomerization domain.

Construction and analysis of PPI network

Four DEGs (JUND, JUNB, FOSB and ATF3) were found to exhibit close associations with other genes, via the proteins identified to construct the PPI network (Fig. 3). These results show that the four DEGs interact indirectly, indicating that they may exhibit a joint role in the pathogenesis of breast cancer.

Figure 3. Protein-protein interaction network of differentially expressed genes (green circles) in breast cancer.

Screening of target sites of potential transcription factors

The expression of transcription factors affects the regulation downstream gene expression. Exploring the target sites of possible transcription factors is of great importance. The target sites of potential transcription factors were explored, and the ten most significant sites are listed in Table III. As shown in Table III, FOSB, ATF3 and JUND shared the same binding sites [Homo sapiens (hsa)_V$ATF_01 and hsa_V$ATF3_Q6]. Therefore, the regulation of these binding sites may present a potential treatment for breast cancer by controlling the expression of FOSB, ATF3 and JUND.

Table III. Target sites of potential transcription factors.

Table III.

Target sites of potential transcription factors.

Target	P-value	Gene count	Gene symbols
hsa_V$CREB_Q2	3.51×10−7	20	GNAS, EGR1, CDC42, CYLD, FOSB, ATF3, FOS
hsa_V$CREB_02	1.18×10−6	16	C9orf3, FOSB, NR4A2, DUSP1, JUND
hsa_V$CREB_Q4	9.21×10−6	15	GNAS, CDC42, EIF1, NR4A2, ATF3, FOS
hsa_V$ATF_01	9.21×10−6	15	EIF1, FOSB, NR4A2, ATF3, SIK1, FOS, JUND
hsa_V$E4F1_Q6	1.10×10−5	16	GNAS, NR4A2, EIF1, FOSB, EIF4A1, HIC2
hsa_V$TAXCREB_01	1.10×10−5	10	CYLD, CDC42, FOSB, FOS, JUND, DUSP1
hsa_GTGACGY_V$E4F1_Q6	2.86×10−5	22	GNAS, PTP4A1, NR4A3, EGR1, JUNB, EIF1, JUND
hsa_V$CREB_Q2_01	4.56×10−5	16	JUN, CYLD, EGR1, FOSB, NR4A1, DUSP1
hsa_V$SRF_C	4.56×10−5	13	DUSP2, IER2, EGR1, JUNB, FOSB, KLF6, FOS
hsa_V$ATF3_Q6	1.00×10−4	13	JUN, CYLD, CDC42, FOSB, ATF3, JUND

[i] hsa, Homo sapiens; V$, separator between transcription factor and its binding sequence.

Screening of potential miRNAs

miRNAs regulate gene expression by controlling the stability of RNA. Therefore, identifying the regulatory miRNAs of specific DEGs may increase understanding with regard to the involvement of miRNAs in breast cancer. Potential regulatory miRNAs were analyzed based on the sequences of the DEGs, and the target sites with P-values of <0.05 are listed in Table IV. The regulatory miRNAs of two target sites, hsa_AGCACTT and hsa_ACTTTAT, were collected, which included miR-93, miR-302A, miR-302B, miR-302C, miR-373 and miR-520. These miRNAs may present potential therapeutic targets for breast cancer treatment.

Table IV. Potential miRNAs of differentially expressed genes with P<0.05.

Table IV.

Potential miRNAs of differentially expressed genes with P<0.05.

Target sequence	Potential miRNAs	P-value
hsa_AGCACTT	miR-93, miR-302A, miR-302B, miR-302C, miR-302D, miR-372, miR-373	0.0182
	miR-520E, miR-520A, miR-526B, miR-520B, miR-520C, miR-520D
hsa_ACTTTAT	miR-142-5P	0.0182

[i] miR or miRNA, microRNA; hsa, Homo sapiens.

Discussion

Breast cancer, similar to other cancer types, occurs due to an interaction between the environment and a defective gene (3). The incidence of breast cancer is increasing rapidly in the majority of Asian countries (18). Worldwide, breast cancer accounts for 22.9% of all cancer cases (excluding nonmelanoma skin cancers) in females. In 2008, breast cancer was the cause of 458,503 mortalities worldwide (13.7% of cancer-associated mortalities in females) (2). Therefore, the study of breast cancer has important significance for human health. In the present study, a total of 106 DEGs between breast cancer and normal cells were identified. With regard to molecular function, the identified DEGs were predominantly involved in DNA-binding and regulation of downstream gene expression, as well as in metabolic and synthetic pathways. In terms of bio-pathways, signaling-associated pathways and other disease-associated pathways were found to be altered. ATF3, JUND, FOSB and JUNB were the only four genes composing the PPI network of DEGs. Finally, the most significant target sites of potential transcription factors and miRNAs in breast cancer were identified.

Functionally, the differential expression of genes results in a variety of abnormal physiological processes, including alterations in gene expression regulation, RNA transcription and protein translation processes (19,20). In the present study, changes in the RNA transport pathway were consistent with the GO clustering results, indicating that gene expression in breast cancer cells differs from that in normal cells. In addition, by KEGG clustering, a number of bio-pathways were revealed to be altered in breast cancer. Among the altered cell surface signaling pathways, alteration of the NOD-like receptor signaling pathway may be an important mechanism involved in the abnormal differentiation of breast cancer cells. By contrast, changes associated with the immune response, including the T cell receptor, the chemokine and the B-cell receptor signaling pathways, may be involved in the evasion of autoimmunity of breast cancer cells. In addition, the MAPK signaling pathway, which was also enriched, transports extracellular signals into the intracellular environment, and may affect downstream gene expression in breast cancer cells (21). Alterations in the aforementioned pathways may result in changes to a number of other metabolic signaling pathways, including the GnRH and the neurotrophin signaling pathways. In addition, differences in the mismatch repair pathway indicated that the self-repair capacity may be affected. Changes in proteasome, metabolic pathways, glutathione metabolism and glycolysis/gluconeogenesis further confirmed the results of the GO clustering analysis, which indicated that the metabolic capability of breast cancer cells was altered compared with that of normal cells.

ATF3 is a member of the mammalian activation transcription factor/cyclic adenosine monophosphate responsive element-binding protein family, which is induced by numerous signals in cancer tissues (22). The ATF3 product forms the activator protein 1 complex by interacting with the protein products of JUND, FOSB and JUNB, thus regulating the expression of downstream genes in response to cytokines, growth factors and cell stress (23). Certain studies have demonstrated that it may be involved in the resistance to p53-dependent cellular senescence and apoptosis (24). Therefore, the existence of the PPI network formed by ATF3, JUND, FOSB and JUNB may result in the inability of the p53 pathway to kill breast cancer cells.

Of the screened transcription factors, the majority are associated with cylindromatosis, cell division cycle 42, FOSB, nuclear receptor subfamily 4 and JUND. Among these genes, JUND has been proposed to protect cells from p53-dependent senescence and apoptosis (25), while FOSB has been implicated as a regulator of cell proliferation, differentiation and transformation (26) These genes are key factors in the pathogenesis of cancer, and the identification of associated transcription factors may help to elucidate the underlying mechanisms of breast cancer.

In the current study, miR-93 was identified as a potential regulatory miRNA of the DEGs in breast cancer. It has been reported that miR-93 is the substrate of caspase-3 during apoptosis (27). In the present study, miR-93 was found to be a potential tumor suppressor using bioinformatics. In addition, miR-302A has been strongly linked with the hypoxia pathway and is upregulated in response to hypoxia-inducible factors (28). Furthermore, this miRNA is also overexpressed in breast tumor cells and was found to be a potential regulatory miRNA of certain DEGs in breast cancer in the current study. miR-302B, a frequently amplified miRNA, is associated with intrahepatic metastasis of hepatocellular carcinoma (29). A previous study identified that miR-302C suppresses the expression of the prostate-specific antigen and prostate cancer cell proliferation (30). The results of the present study suggested that the role of miR-302C in breast cancer was consistent with its role in prostate cancer. Furthermore, a previous study has demonstrated that miR-373 is involved in the development of normal and cancer cells (31). In breast cancer tissues, miR-373 (31) is downregulated compared with normal breast tissues, which is consistent with the results of the current study. miR-520 is a short RNA molecule (32), with only a limited number of studies investigating its role in cancer; however, the present study identified that it may be involved in breast cancer.

In conclusion, a large number of DEGs in breast cancer share the same transcription factors and miRNAs. The target sites of these molecules may be important in the regulation of DEG expression. Therefore, regulating the expression of these genes through the aforementioned target sites may contribute towards the development of novel treatments for breast cancer.

References