This article is mentioned in:

Introduction

Breast cancer is one of the most common malignancies in women worldwide and results in relatively high rates of morbidity and mortality (1–4). In China, the incidence rate of breast cancer has been increasing over the past 20 years (5–7). Depending on the molecular classification, therapies used to treat breast cancer, including surgery, chemotherapy, radiation therapy, hormone (endocrine) therapy and molecule-targeted therapy, vary in their survival rates (8–13). As the prognosis for breast cancer patients remains unsatisfying (14–18), there is an urgent need to identify a more effective therapy.

Several chemical compounds used in traditional Chinese medicine (TCM) have proven useful in some conventional chemotherapies (19,20). Cinobufotalin, a member of the bufadienolide family, is isolated from the skin parotoid glands of toads, such as Bufo gargarizans and Duttaphrynus melanostictus (21). The broad-spectrum antineoplastic activity and chemosensitization of bufadienolide has also been previously reported (22). Another study on cinobufotalin revealed that it may serve as a cardiotonic, diuretic and hemostatic agent (23). Previous reports have also focused on the influences of cinobufotalin on cancers such as hepatocellular carcinoma (HCC) and lung cancer (24,25), but very few studies have examined its mechanism in these malignancies and there are even fewer, if any, reports on its functions in breast cancer. The mechanism of cinobufotalin against breast cancer cells remain unknown.

In the present study, the GSE85871 microarray data set from the Gene Expression Omnibus (GEO) database was used in an optimized analysis to identify differentially expressed genes (DEGs) in MCF-7 breast cancer treated with cinobufotalin. Subsequently, the potential molecular mechanism of cinobufotalin in breast cancer was explored through gene annotation, pathway analysis and protein-protein interaction (PPI) analysis. Connectivity Map (CMAP) analysis was used to identify drugs that may exhibit similar curative properties as cinobufotalin. Based on the mining of a large database, the present study comprehensively revealed the roles of cinobufotalin and its potential molecular mechanism in breast cancer, and offered a possible avenue for breast cancer treatment.

Materials and methods

Data sources and analysis of DEG expression profiles

The expression data of the GSE85871 data set were obtained from the National Center for Biotechnology Information GEO database (26). The subject of this microarray was Homo sapiens, and its research type was expression profiling by array. The expression profile of this microarray was provided by GPL571 (HG-U133A_2; Affymetrix Human Genome U133A 2.0 Array). GSE85871 included the gene expression profiles of 102 TCM compounds used to treat MCF-7 cells triple positive for the presence of estrogen, progesterone and human epidermal growth factor receptor-2 (HER-2) receptors. Additionally, in this microarray, comparisons were made between the experimental groups treated with each drug and their respective dimethylsulphoxide (DMSO) controls; there were two expression values for each drug. The expression profiles of genes in MCF-7 cells treated with cinobufotalin (GSM2286314 MCF-7_Cinobufotalin_1 µM_rep1 and GSM2286315 MCF-7_Cinobufotalin_1 µM_rep2) and the profiles of respective controls [GSM2286316 MCF-7_vehicle (DMSO)_rep1 and GSM2286317 MCF-7_vehicle (DMSO)_rep2] were downloaded from the GSE85871 data set. Fold change (FC) was set as the threshold for the mean value of gene expression in the experimental groups and in the respective DMSO controls, as previously described (26); DEGs were identified as FC≥2 or FC≤0.5, and categorized as upregulated or downregulated, respectively.

Functions of DEGs and pathway enrichment analysis

The STRING online tool (https://string-db.org) was used for a PPI network analysis. The top 10 genes that possessed the most protein interactions were considered hub genes, and their genetic alterations were displayed using The Cancer Genome Atlas (TCGA; http://cancergenome.nih.gov) data calculated with the cBioPortal (27,28). The cBioPortal dataset included 3 different data: TCGA RNA-sequencing, Provisional; Genomic Identification of Significant Targets In Cancer (GISTIC) data set; and, a microarray included in the TCGA project. mRNA and protein expression levels and clinical significance of these hub genes were computed by Gene Expression Profiling Interactive Analysis (GEPIA) and The Human Protein Atlas (29–34); TCGA and The Genotype-Tissue Expression (GTEx) RNA-sequencing data were included in GEPIA. The expression levels of genes of interest were extracted from both TCGA and GTEx projects. Transcripts Per Kilobase Million (TPM) was used to calculate mRNA expression levels, and the relative expression levels of genes were presented as Log2(TPM + 1). Student's t-test was used to compare the differences in expression between cancer and non-cancerous tissues; P<0.05 was considered to indicate a statistically significant difference. Gene Ontology (GO; http://www.geneontology.org) and Kyoto Encyclopedia of Genes and Genomes (KEGG; http://www.genome.jp/kegg) pathway analyses were carried out for DEGs post-cinobufotalin treatment of breast cancer cells in GSE85871 through DAVID (https://david.ncifcrf.gov) (35–39). Cytoscape (https://cytoscape.org) and R software were used to display the top 10 terms (40–45).

CMAP matching of cinobufotalin with existing compounds

The CMAP database (https://clue.io/cmap) was used to further investigate the relationship between the effectiveness of cinobufotalin on treating breast cancer cells and to explore its pharmacology and mechanism. The enrichment value in CMAP represented a similarity in gene mapping. A value close to 1 indicated a positive relationship between two drug molecules, which suggested that cinobufotalin exhibited similar curative properties to the existing compounds. Conversely, a negative value indicated that cinobufotalin exhibited opposite functions to the existing compounds (46–48). The molecular structures of certain compounds were retrieved from PubChem (https://pubchem.ncbi.nlm.nih.gov/substance), including cinobufotalin (PubChem CID: 259776; Fig. 1).

Figure 1. Molecular structure of cinobufotalin. Cinobufotalin chemical structure obtained from PubChem; PubChem CID: 259776. (A) 2D molecular structure. (B) 3D molecular structure.

Results

Analysis of DEGs and pathway enrichment

A total of 1,237 DEGs were identified in cinobufotalin-treated MCF-7 cells compared with the DMSO-treated controls. Of these, 641 genes were upregulated and 596 were downregulated (Fig. 2). To determine the protein interactions of these 1,237 DEGs, a PPI network was constructed using STRING (Fig. 3). PPI network analysis revealed several hub genes, including: Albumin (ALB); SRC proto-oncogene non-receptor tyrosine kinase (SRC); glucagon; protein tyrosine phosphatase receptor type C (PTPRC); spectrin α non-erythrocytic 1; coagulation factor II, thrombin; FYN proto-oncogene, SRC family tyrosine kinase (FYN); cyclin-dependent kinase inhibitor 2A (CDKN2A); transferrin (TF) and insulin-like growth factor 1. The associations between these 10 genes are demonstrated in Fig. 4. As these genes could be the targets of cinobufotalin treatment in breast cancer cells, the clinical roles of these genes in breast cancer were next investigated using various approaches. Through cBioPortal data mining (Fig. 5), varied genetic alterations were observed in these genes in clinical breast cancer tissue samples, including amplification, a number of different mutations, and changes in mRNA and protein expression levels detected by multiple approaches. These genetic alterations were identified by three independent sources as provided by cBioPortal, including TCGA RNA-sequencing, Provisional (Fig. 5A); GISTIC data set (Fig. 5B) and a microarray included in the TCGA project (Fig. 5C). Based on the TCGA data, the amplification of ALB, SRC and PTPRC were the main genetic alteration events identified (Fig. 5A). In addition, the corresponding mRNA and protein expression levels of ALB, FYN and TF tended to decrease in clinical breast cancer tissue samples (Fig. 6). The mRNA expression levels of ALB, FYN and TF were predominantly downregulated based on RNA-sequencing data with 1,085 cases of breast cancer and 291 non-cancerous breast tissues (all P<0.05), and the protein expression levels of ALB, FYN and TF were also decreased in breast cancer compared with the non-cancerous control tissues. The mRNA and protein expression levels of SRC and CDKN2A were increased in breast cancer tissue samples, compared with the control tissue (Fig. 7). However, owing to the limited number of cases provided by The Human Protein Atlas, statistical analysis was not possible. The protein expression levels of the above genes need to be confirmed using a larger cohort. Notably, following treatment with cinobufotalin, their expressions were remarkably reduced, with CDKN2A expression dropping to 42.2% and SRC plunging to 7.03% (data not shown). These results indicated that cinobufotalin is more likely to target these genes to exert its antitumor influences.

Figure 2. DEG microarray analysis in MCF-7 cells treated with cinobufotalin. A total of 1,237 DEGs were identified between cinobufotalin-treated and DMSO-treated MCF-7 breast cancer cells, of which 641 genes were upregulated and 596 were downregulated using cut-off values set as FC≥2 or FC≤0.5, respectively. DEGs, differentially expressed genes; DMSO, dimethylsulphoxide.

Figure 3. PPI network. A PPI network was created for the identified differentially expressed genes between MCF-7 breast cancer cells treated with cinobufotalin. Network nodes represent proteins; edges represent protein-protein associations. PPI, protein-protein interaction.

Figure 4. Protein-protein interaction network of the 10 hub genes in cinobufotalin-treated MCF-7 cells. Network nodes represent proteins; edges represent protein-protein associations. ALB, albumin; CDKN2A, cyclin-dependent kinase inhibitor 2A; GCG, glucagon; F2, coagulation factor II, thrombin; FYN, FYN proto-oncogene, SRC family tyrosine kinase; IGF1, insulin-like growth factor 1; PTPRC, protein tyrosine phosphatase receptor type C; SPTAN1, spectrin α non-erythrocytic 1; SRC, SRC proto-oncogene non-receptor tyrosine kinase; TF, transferrin.

Figure 5. Genetic alterations of the 10 hub genes in breast cancer tissue samples identified in cinobufotalin-treated breast cancer tissues. The genetic alterations were from cBioPortal comprising 1,105 cases of breast cancer tissues (TCGA, Provisional). (A) Mutations and putative copy-number alterations from Genomic Identification of Significant Targets In Cancer from the TCGA project. (B) mRNA expression Z-Score (a microarray from TCGA project), protein/phosphoprotein level Z-Score, threshold=1.5. (C) mRNA expression Z-Score (RNA Seq V2 RSEM from TCGA project), protein/phosphoprotein level Z-Score (mass spectrometry by CPTAC), threshold=1.5. ALB, albumin; CDKN2A, cyclin-dependent kinase inhibitor 2A; GCG, glucagon; F2, coagulation factor II, thrombin; FYN, FYN proto-oncogene, SRC family tyrosine kinase; IGF1, insulin-like growth factor 1; PTPRC, protein tyrosine phosphatase receptor type C; SPTAN1, spectrin α non-erythrocytic 1; SRC, SRC proto-oncogene non-receptor tyrosine kinase; TCGA, The Cancer Genome Atlas; TF, transferrin.

Figure 6. mRNA and protein expressions of ALB, FYN and TF in breast cancer tissues. (A-I) The GEPIA database, which included TCGA and GTEx RNA-sequencing data, was used to examine (A-C) ALB, (D-F) FYN and (G-I) TF mRNA and protein expression in breast cancer tissues and non-cancerous breast tissues; TPM was used to calculate relative mRNA gene expression levels from 1,085 tumor and 291 non-tumor tissues, which are presented as Log2(TPM + 1) scale. (A) ALB mRNA expression levels from GEPIA. (B) Non-cancerous breast myoepithelial cells stained with HPA031024 anti-ALB antibody (medium staining). (C) Breast cancer tissue stained with HPA031024 anti-ALB antibody (low staining with light brown color in the cytoplasm). (D) FYN mRNA expression levels from GEPIA. (E) Non-cancerous breast myoepithelial cells stained with anti-FYN HPA023887 antibody (medium staining). (F) Breast cancer tissue stained with anti-FYN HPA023887 antibody (low staining with light brown color in the cytoplasm). (G) TF mRNA expression levels from GEPIA. (H) Non-cancerous breast myoepithelial cells stained with HPA001527 anti-TF antibody (medium staining). (I) Breast cancer tissue stained with HPA001527 anti-TF antibody (low staining with light brown color in the cytoplasm). *P<0.01. Immunohistochemistry images in B, C, E, F, H and I are presented at ×400 magnification. ALB, albumin; FYN, FYN proto-oncogene, SRC family tyrosine kinase; GEPIA, Gene Expression Profiling Interactive Analysis database; GTEx, Genotype-Tissue Expression; TCGA, The Cancer Genome Atlas; TF, transferrin; TPM, transcripts per kilobase million.

Figure 7. mRNA and protein expressions of SRC and CDKN2A in breast cancer tissues. (A-F) The GEPIA database, which included TCGA and GTEx RNA-sequencing data, was used to examine (A-C) SRC and (D-F) CDKN2A mRNA and protein expression in breast cancer tissues and non-cancerous breast tissues; TPM was used to calculate relative mRNA gene expression levels from 1,085 tumor and 291 non-tumor tissues, which are presented as Log2(TPM + 1) scale. (A) SRC mRNA expression levels from GEPIA. (B) Non-cancerous breast myoepithelial cells stained with CAB004023 anti-SRC antibody (not detected; brown staining was in the lumen of the ducts, not in the myoepithelial cells. (C) Breast cancer tissue stained with CAB004023 anti-SRC antibody (medium staining). (D) CDKN2A mRNA expression levels from GEPIA. (E) Non-cancerous breast myoepithelial cells stained with CAB018232 anti CDKN2A antibody (not detected). (F) Breast cancer tissue stained with CAB018232 anti CDKN2A antibody (high staining). *P<0.01. Immunohistochemistry images in B, C, E and F are presented at ×400 magnification. CDKN2A, cyclin-dependent kinase inhibitor 2A; GEPIA, Gene Expression Profiling Interactive Analysis database; GTEx, Genotype-Tissue Expression; SRC, SRC proto-oncogene non-receptor tyrosine kinase; TCGA, The Cancer Genome Atlas; TPM, transcripts per kilobase million.

To explore the functions of these 1,237 DEGs, GO and KEGG analyses were conducted using DAVID. In the GO analysis, the genes were divided into three groups: i) Biological process (BP), ii) cellular component (CC) and iii) molecular function (MF). In BP, the three most significantly enriched processes were ‘immune response’, ‘chemical synaptic transmission’ and ‘cellular defense response’ (Fig. 8; Table IA). In CC, the three most significant cellular components were ‘integral component of plasma membrane’, ‘extracellular space’ and ‘plasma membrane’ (Fig. 9; Table IB). In MF, the top three functions were ‘receptor binding’, ‘immunoglobulin receptor binding’ and ‘cytokine activity’ (Fig. 10; Table IC). KEGG pathway analysis confirmed that DEGs were remarkably enriched in pathways of ‘neuroactive ligand-receptor interaction’, ‘calcium signaling pathway’ and ‘steroid hormone biosynthesis’ (Fig. 11; Table ID). The top five GO terms for each category and KEGG pathways concerning cinobufotalin treatment on breast cancer cells in GSE85871 are displayed in Table I.

Figure 8. Biological process of differentially expressed genes in MCF-7 cells treated with cinobufotalin. Circles represent different GO terms; colors of circles indicate the significance of the corresponding terms, and the sizes of circles display the numbers of genes enriched in this term. GO, Gene Ontology.

Figure 9. Cellular component of differentially expressed genes in MCF-7 cells treated with cinobufotalin. The circles represent different GO terms. The colors of circles indicate the significance of the corresponding terms, and the sizes of circles display the numbers of genes enriched in this term. GO, Gene Ontology.

Figure 10. Molecular function of differentially expressed genes in MCF-7 cells treated with cinobufotalin. The circles represent different GO terms. The colors of circles indicate the significance of the corresponding terms, and the sizes of circles display the numbers of genes enriched in this term. GO, Gene Ontology.

Figure 11. Kyoto Encyclopedia of Genes and Genomes pathways of differentially expressed genes in MCF-7 cells treated with cinobufotalin. Circle color indicates the significance of the corresponding terms, and circle size represents the number of genes enriched in this term.

Table I. Significant terms of GO annotation and KEGG pathway enrichment analysis of cinobufotalin in GSE85871.

Table I.

Significant terms of GO annotation and KEGG pathway enrichment analysis of cinobufotalin in GSE85871.

A, Biological processa
ID	Term	Count (%)	P-value	Genes
GO:0006955	Immune response	65 (5.36)	1.15×10−11	IGLV1-44, IL16, IL19, TRGC2, TNFSF14, IL15, ZEB1, CXCL11, CFP and CD96
GO:0007268	Chemical synaptic transmission	44 (6.36)	1.63×10−10	SLC6A2, OPRK1, TACR1, NPY2R, DRD5, GRIK4, KCNA1, CACNB1, TAC1 and PMCHL1
GO:0006968	Cellular defense response	21 (7.36)	3.51×10−10	IL4, NCF1, LY96, IL1RL2, CD160, GAGE2D, GAGE2E, CXCR2, UMOD and KIR2DS3
GO:0007165	Signal transduction	121 (8.36)	6.20×10−09	GNA13, EDN3, ADCY8, IL19, TRGC2, RRAD, IQGAP2, TNFSF14, IL15 and CXCL11
GO:0050776	Regulation of immune response	34 (9.36)	9.12×10−09	IGLV1-44, KLRK1, ITGB2, CD96, SH2D1A, IGHV3-23, IGLV6-57, KLRF1, IGKC and KLRD1
B, Cellular componenta
ID	Term	Count (%)	P-value	Genes
GO:0005887	Integral component of plasma membrane	186 (15.33)	6.73×10−27	F2RL2, ADCY1, SLC22A13, SLC6A2, GRIK4, F2RL1, CSPG4, LPAR4 and TRGC2
GO:0005615	Extracellular space	165 (16.33)	1.47×10−20	EDN3, IL16, IL19, SNCA, TNFSF14, IL15, CXCL11, IGHM, BMP15 and TGFB2.
GO:0005886	Plasma membrane	365 (17.33)	1.20×10−19	ADCY1, CROCC, SCN3B, ADCY8, SLC9A3, SNCA, SLC9A2, LPAR4, AQP4 and PVRIG
GO:0005576	Extracellular region	181 (18.33)	1.34×10−18	F2RL2, EDN3, IL16, IGLV1-44, MASP2, IL19, SNCA, CSPG4, IL15 and MMRN1
GO:0009897	External side of plasma membrane	39 (19.33)	5.71×10−10	IGHG1, IGHG2, IGHG3, IGHG4, LDLR, FCER2, KLRK1, AQP4, HFE and IGHM
C, Molecular functiona
ID	Term	Count (%)	P-value	Genes
GO:0005102	Receptor binding	51 (4.20)	9.76×10−09	INSL3, MBL2, EDN3, ACOX1, PPY, PLXNC1, BLK, F2RL1, HFE and TAC1
GO:0034987	Immunoglobulin receptor binding	11 (5.20)	1.18×10−06	IGHG1, IGHG2, IGHG3, IGHG4, IGHD, IGHV3-23, IGHA1, IGHA2, IGKC, IGHM, IGLC1
GO:0005125	Cytokine activity	29 (6.20)	1.78×10−06	CSF2, IFNA21, IL16, IL19, KITLG, TNFSF14, IL15, BMP15, CCL4 and FLT3LG
GO:0003823	Antigen binding	21 (7.20)	2.38×10−06	IGHG1, IGHG2, IGHG3, IGHG4, IGLV1-44, HFE, CD40, ITGA4, IGHM and LILRA2
GO:0004872	Receptor activity	32 (8.20)	5.21×10−06	NRP2, LDLR, NPY2R, F2RL1, KLRK1, ITGB2, ITGB3, TLR8, LILRA2 and LILRA3
D, KEGG pathwaya
ID	Term	Count (%)	P-value	Genes
hsa04080	Neuroactive ligand receptor interaction	61 (5.03)	2.65×10−14	CSH1, F2RL2, MCHR1, CSH2, GRIK3, TACR1, GRIK4, F2RL1, GLRA2 and LPAR4
hsa04020	Calcium signaling pathway	30 (6.03)	6.86×10−05	ADCY1, ADCY8, TACR1, DRD5, BDKRB2, EDNRA, ATP2B2, EDNRB, HRH1 and PDE1C
hsa00140	Steroid hormone biosynthesis	15 (7.03)	7.84×10−05	CYP3A4, STS, CYP1B1, CYP11A1, CYP2E1, UGT1A1, UGT1A9, UGT2B17, UGT1A8 and AKR1C4
hsa00980	Metabolism of xenobi otics by cytochrome P450	17 (8.03)	1.03×10−04	CYP3A4, CYP1B1, CYP2F1, SULT2A1, CYP2B6, ADH1C, CYP2E1, UGT1A1, ALDH3A1 and UGT1A9
hsa04060	Cytokine-cytokine receptor interaction	35 (9.03)	1.12×10−04	IFNA21, CSF2, KITLG, TNFSF14, CXCR2, IL15, CCL5, CXCL11, CCL4 and TGFB2

a Only the top five terms or pathways are provided as examples. GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Results of CMAP analysis

A gene expression profile comparison was conducted using the aforementioned 1,237 DEGs post-cinobufotalin treatment and those genes related to the known drugs from the whole CMAP database and 9 compounds yielded scores >0.962 (Table II), including trichostatin A (appearing twice), esculetin, metixene, niclosamide, 15-delta prostaglandin J2, pimethixene, acetohexamide, allantoin and pregnenolone. The 1,237 DEGs generated after cinobufotalin treatment in MCF-7 cells were also compared with the genes following treatment with other drugs from the CMAP project in MCF-7 and PC3 cell lines. A total of 11 compounds were chosen with P<0.05. The results revealed that the compounds exhibiting similar roles to cinobufotalin included miconazole, salbutamol, dexibuprofen, ciprofloxacin, nialamide, scopolamine N-oxide and cinnarizine, whereas those exhibiting an opposite effect included triamterene, iloprost, BCB000040 and BCB000038 (Table III; Fig. 12). Among these drugs, most of them were generated from the same cell line, MCF-7, which enhanced the power of drug prediction, e.g. miconazole, salbutamol, iloprost, dexibuprofen, ciprofloxacin, nialamide, scopolamine N-oxide and cinnarizine may have comparable target genes with cinobufotalin in breast cancer cells.

Figure 12. Chemical structures of 11 compounds acquired from the CMAP database. (A) Miconazole. (B) Salbutamol. (C) Dexibuprofen. (D) Ciprofloxacin. (E) Nialamide. (F) Scopolamine N-oxide. (G) Cinnarizine. (H) Triamterene. (I) Iloprost. (J) BCB000040. (K) BCB000038. CMAP, Connectivity Map database.

Table II. Total rank of the compounds from CMAP analysis compared with cinobufotalin in MCF-7 breast cancer cells.

Table II.

Total rank of the compounds from CMAP analysis compared with cinobufotalin in MCF-7 breast cancer cells.

Rank	Batch from CMAP	Compound name	Dose	Connectivity score
1	627	Trichostatin A	100 nM	1.00
2	772	Esculetin	22 µM	1.00
3	655	Metixene	12 µM	0.99
4	632	Niclosamide	12 µM	0.99
5	655	Trichostatin A	100 nM	0.99
6	765	15-delta prostaglandin J2	10 µM	0.98
7	772	Pimethixene	10 µM	0.97
8	772	Acetohexamide	12 µM	0.97
9	627	Allantoin	25 µM	0.97
10	656	Pregnenolone	13 µM	0.96

[i] CMAP, Connectivity Map database. Trichostatin A appeared twice.

Table III. Rank of the compounds from CMAP analysis compared with cinobufotalin in a single cell line.

Table III.

Rank of the compounds from CMAP analysis compared with cinobufotalin in a single cell line.

Rank	Compound name	Cell line	Enrichment	P-value
1	Miconazole	MCF-7	0.966	0.002
2	Triamterene	PC3	−0.964	0.003
3	Salbutamol	MCF-7	0.953	0.004
4	Iloprost	MCF-7	−0.937	0.008
5	BCB000040	PC3	−0.929	0.010
6	Dexibuprofen	MCF-7	0.920	0.013
7	BCB000038	PC3	−0.908	0.017
8	Ciprofloxacin	MCF-7	0.906	0.018
9	Nialamide	MCF-7	0.904	0.019
10	Scopolamine N-oxide	MCF-7	0.888	0.026
11	Cinnarizine	MCF-7	0.872	0.033

[i] CMAP, Connectivity Map database.

Discussion

Currently, to the best of the authors' knowledge, there have been no reports on the antitumor mechanisms of cinobufotalin in breast cancer cells through large data mining analyses. Microarrays and RNA-sequencing have facilitated research on functions and mechanisms of TCM (49–52). The present study was conducted by combining microarray analysis and RNA-sequencing data in breast cancer tissues. For the potential target gene of cinobufotalin, several genes were selected for confirmation and demonstrated that ALB, FYN, TF, SRC and CDKN2A may serve pivotal roles in the onset and development of breast cancer. These genes also were affected by cinobufotalin in treated MCF-7 cells, which may shed light on the potential mechanism of cinobufotalin on breast cancer cells. In the present study, the pathways of ‘neuroactive ligand-receptor interaction’ and ‘calcium signaling’ appeared to be significant pathways for cinobufotalin in MCF-7 breast cancer cells. In addition, connectivity mapping demonstrated that cinobufotalin had similar molecular mechanisms to drugs such as miconazole as they target consistent genes in breast cancer cells, which may provide a theoretical foundation for research on the anticancer mechanism of cinobufotalin in breast cancer cells.

The anticancer ability of cinobufotalin has been previously documented in hepatoblastoma (24) and lung cancer cells (25), based on in vitro models. In HepG2 hepatoblastoma cells, cinobufotalin was reported to inactivate Akt-S6K1 signaling, and in A549, H460 and HTB-58 lung cancer cells, cinobufotalin mainly induced Cyclophilin D-dependent non-apoptotic death. Data from PubChem also revealed that cinobufotalin exhibited effects on other cancer cells. For instance, cinobufotalin cytotoxicity against human Bel7402 cells, which have been identified as being derived from Hela cells, was detected by MTT assay (BioAssay AID: 343717) and the activity value (IC50) was 1.21 mM. Another BioAssay (AID: 1221865) indicated an activity value (IC 50) of 8.62 mM when cytotoxicity against human Bel7402 cells was assessed after treatment of 72 h. Interestingly, a phase I clinical trial sponsored by Shanghai University of TCM is now at recruitment stage and will use cinobufotalin injection as intervention to treat malignant tumor of small intestine metastatic to liver (https://clinicaltrials.gov; ClinicalTrials.gov Identifier: NCT03189992). However, no previous study has examined the effects and mechanism of cinobufotalin on breast cancer cells. From the expression data provided by the GSE85871 data set in MCF-7 cells following cinobufotalin treatment (26), the present study identified 1,237 DEGs, and subsequently conducted further analysis of the core genes disclosed by PPI. Additional analysis demonstrated that some of these core genes, to some extent, may influence the onset and development of breast cancer through their abnormal expression and genetic alteration. According to the data in TCGA and The Human Protein Atlas, the mRNA and protein expression levels of SRC and CDKN2A were increased in breast cancer tissues compared with non-cancerous tissues. Previous studies have also reported that increased SRC and CDKN2 expression levels correlated with the onset, metastasis and prognosis of breast cancer (53–59). Therefore, SRC and CDKN2 may be the most important hub genes in the biological function of cinobufotalin on breast cancer MCF-7 cells, as in this study, cinobufotalin was observed to inhibit the overexpression of SRC and CDKN2A (data not shown). The potential targeting of these genes suggested that cinobufotalin may have anticancer potential.

In addition, potential mechanisms of cinobufotalin in MCF-7 breast cancer cells were elucidated from the prospective signaling pathways. The 1,237 DEGs were annotated to examine how cinobufotalin functioned on breast cancer cells. Notably, the GO term results of the 1,237 DEGs were mainly linked to immunity, including ‘immune response’, ‘regulation of immune response’, ‘innate immune response’ and ‘adaptive immune response’. It was reported previously that cinobufotalin activated the nuclear factor-κB pathway and decreased expression levels of brain-derived neurotrophic factor to induce neuroinflammation in rats (60). Nonetheless, no direct evidence has revealed that cinobufotalin was associated with number of immune cells, organism immunity or tumor immunity. On a molecular level, the present study hypothesized that in treating breast cancer, cinobufotalin may exert its antitumor influences by activating the immune response; however, additional studies are required for verification.

To further interpret the potential mechanism of cinobufotalin in MCF-7 breast cancer cells, KEGG pathway enrichment was performed on the identified DEGs, which revealed that several pathways were connected, not only to immunity after cinobufotalin was applied to treat breast cancer samples but also to other pathways. ‘Neuroactive ligand-receptor interaction’, which contains numerous G protein-coupled receptors, continued to be the most significantly enriched pathway. As one of the most common pathways of malignancies, the pathway of ‘Neuroactive ligand-receptor interaction’ ranks fifth in genes with mutation in the central nervous system (61). A previous study reported that low expression of cannabinoid receptor-1 (CNR1), part of the neuroactive ligand-receptor interaction pathway, indicated that breast cancer patients may benefit from chemotherapy (62). Notably, this CNR1 was among the 1,237 DEGs following cinobufotalin treatment in breast cancer cells. In addition, this pathway was reported to serve a vital role in tumorigenesis and chemotherapy for breast cancer (62,63), which suggested that cinobufotalin may exhibit its curative potential for breast cancer by influencing the neuroactive ligand-receptor interaction pathway. The co-treatment of cinobufotalin to chemotherapeutics may induce a synergistic effect by suppressing the neuroactive ligand-receptor interaction pathway.

The second most enriched KEGG pathway was the ‘calcium signaling pathway’, the mechanism of which is rather complicated in breast cancer. Previous studies have reported that the calcium signaling pathway interacts with other pathways to contribute to the onset and progression of breast cancer (64–68). In invasive ductal carcinoma the calcium signaling pathway was reported to interact with pathways in cancer, including the pathways of glyoxylate and dicarboxylate metabolism, basal transcription factors, tyrosine metabolism, FcγR-mediated phagocytosis, metabolism of xenobiotics by cytochrome P450 and phagosome (69). Wnt5a in the noncanonical Wnt pathway was considered a possible anti-oncogene in breast cancer, as it was demonstrated to serve an essential role via the calcium signaling pathway (70). In addition, the calcium signaling pathway was strongly linked with the proliferation, migration and invasion of breast cancer cells (71). Based on these results, cinobufotalin may affect breast cancer cells by inhibiting the calcium signaling pathway.

TCMs differ from other chemical compounds owing to their natural ingredients and low toxicity (72–75). CMAP analysis revealed that cinobufotalin and miconazole shared similar mechanisms, as the varied genes post-miconazole treatment were comparable to those following cinobufotalin treatment in the same breast cancer cell line MCF-7. It was previously demonstrated that miconazole was able to induce apoptosis in bladder cancer cells through the death receptor 5-dependent and mitochondrial-mediated pathways (76); thus, cinobufotalin may also cause the death of cancer cells (24,25). Notably, miconazole activates the release of phospholipase C-dependent Ca2+ from the endoplasmic reticulum by influencing the elevation of calcium ions, thus inducing ZR-75-1 breast cancer cell apoptosis (77). Based on these previous studies and the present CMAP analysis, it was hypothesized that cinobufotalin may induce the apoptosis of cancer cells in breast cancer MCF-7 cells via the calcium signaling pathway, thus resembling the curative properties of miconazole. However, this hypothesis requires additional experiments for confirmation.

Several limitations exist in the present study. First, the cellular model MCF-7 only represents a specific subtype of breast cancer. It is probable that cinobufotalin may serve different functions through various genes and pathways in distinct subtypes of breast cancer. Therefore, additional cell lines should be examined in future studies. Second, the current findings are based on in silico analyses and verification with in vitro and in vivo experiments are needed, including the biological effect and potential molecular mechanism.

In conclusion, cinobufotalin is likely to act as an effective compound to treat this subtype of breast cancer, triple positive for estrogen, progesterone and HER-2 receptors, and its mechanism may correlate to various pathways. In addition, cinobufotalin may have anticancer functions in MCF-7 cells similar to those of miconazole.

Acknowledgements

The authors are grateful for the public microarray and RNA-sequencing data.

Funding

The present study was supported by The Guangxi Provincial Health Bureau Scientific Research Project (grant no. Z2013424).

Availability of data and materials

All data generated or analyzed during this study are included in this published article.

Authors' contributions

JL, MHR, YWD, RQH, JCZ, JM and GC designed the study. JL, MHR, XJL, DDX, LJZ, HQ, CXF and GC conceived and performed the experiments. JL, MHR, XJL, DDX, LJZ, HQ, CXF, QL and HY performed the statistical analysis and data interpretation. JL, MHR and RQH wrote the manuscript, HY, JCZ and GC corrected the paper. All authors approved the final version of the manuscript for publication.

Ethics approval and consent to participate

Not applicable.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

References