Revealing the role of VEGFA in clear cell sarcoma of the kidney by protein-protein interaction network and significant pathway analysis

Introduction

Clear cell sarcoma of the kidney (CCSK) is a rare, malignant childhood tumor occurring most frequently in the second and third years of life (1). CCSK has a wide histological spectrum and may often resemble other pediatric renal tumors, resulting in considerable diagnostic challenges. The clinical presentation of CCSK includes a large tumor size (mean diameter, 11.3 cm), a mucoid texture, foci of necrosis and prominent cyst formation (2). Typically, the therapeutic drugs used to treat CCSK include vincristine, doxorubicin and dactinomycin (3). The diagnosis of CCSK is important for its prevention, management and treatment. However, only 68% of CCSKs were correctly identified by institutional pathologists during the last cooperative group clinical protocol (4). Therefore, there is an urgent requirement to identify novel biomarkers for CCSK diagnosis.

Over the previous decades, genetic studies regarding CCSK have been limited and the underlying mechanism of the disease consequently remains unclear. Dhamne et al (5) reported that cyclin D1 (CCND1) may be a central molecule in the pathogenesis of CCSK, primarily regulated by nuclear factor κB (NF-κB) cells. Forkhead box D1 and cbp/p300-interacting transactivator, with Glu/Asp-rich carboxy-terminal domain 1, are highly expressed in CCSK, with oxidative-stress-responsive kinase 1, an early embryonic marker, also expressed at high levels in CCSK and at greater levels than observed in normal kidneys or Wilms' tumors (WT) (6). Little et al (7) identified the dysregulation of epidermal growth factor receptor (EGFR) at multiple levels within CCSK. Boo et al (8) detected the expression of vascular endothelial growth factor (VEGF) and angiopoietin 2 in CCSK tumor cells, indicating that angiogenesis in CCSK may involve these proteins. Furthermore, Short et al (9) reported that CCSK is typically hypermethylated at the CpG of thrombospondin-1, which is a suspected angiogenic factor, suggesting that angiogenesis serves a role in CCSK. Therefore, the specific role of angiogenesis in CCSK development requires full investigation, with the aim of novel targets being identified.

In the present study, the mechanism of CCSK was investigated using pathway analysis and protein-protein interaction (PPI) networks, whilst the role of VEGFA in CCSK development was predicted through gene ontology (GO) and transcription factor binding site (TFBS) analyses. The results demonstrated that VEGFA and its associated genes serve important roles during CCSK development and progression.

Materials and methods

Affymetrix chip data

The Gene Expression Omnibus (GEO) database located in the National Center for Biotechnology Information is currently the largest fully public gene expression resource, including 214,268 samples and 4,500 platforms (10). The chip dataset GSE2712, which includes 14 CCSK, 15 WT and 3 fetal kidney (FK) samples, was downloaded from GEO (4). In the present study, the gene expression patterns of 14 CCSKs were compared to 3 FK samples through the use of oligonucleotide arrays. The platform used was the GPL96 (HG-U133A) GeneChip® Human Genome 2.0 U133A Array (Affymetrix, Inc., Santa Clara, CA, USA).

Data preprocessing and the identification of differentially-expressed genes (DEGs)

The Affy R package is a package of functions and classes for the analysis of oligonucleotide arrays (11). Data preprocessing was performed using robust multi-array average in the Affy R package, including background-adjusted, normalized and log-transformed values (12). The original CEL file data was then transformed to probe-level. Following this, the probe-level data was converted to gene symbols by the perl procedure (13). When several probes corresponded to one gene symbol, the gene expression value was calculated as the mean of each probe. DEGs between the 14 CCSKs and 3 FK samples were identified by Significance Analysis of Microarrays using the R package, with |log fold-change| >1.5 and false discovery rate <0.05 (14,15).

Functional analysis

In the present study, the Database for Annotation, Visualization and Integrated Discovery (DAVID) was applied to conduct Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses of the DEGs. The DAVID bioinformatics resource consists of an integrated biological knowledge base and analytical tools aimed at systematically extracting biological meaning from large gene or protein lists (16). KEGG is a knowledge base that allows for the systematic analysis of gene functions (17). For the present study, the significant KEGG pathways, in which the DEGs were enriched, were screened out with P<0.05 and the most notable pathways were selected for further analysis.

PPI analysis

PPI analysis was conducted on the DEGs that were enriched in the most significant KEGG pathways through the use of the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING). The genes with PPI ≥30 were described as hot genes. The STRING database (http://string-db.org/) provides global protein interactions and associations (18). The PPI network was then visualized using Cytoscape software (19).

GO analysis

The up- and downregulated genes that interacted with hot genes were run through GO analysis by DAVID. GO terms are significantly overrepresented in a set of genes from three aspects, consisting of the cellular component, molecular function (MF) and biological process (BP) (20). The significant GO terms with P<0.05 were selected for further analysis.

TFBS analysis

TFBS analysis implements objects representing specific profile matrices, binding sites and sets thereof, pattern generators and pattern database interfaces (21). In the present study, the whole genome rVista (http://genome.lbl.gov/cgi-bin/WGRVistaInput5.pl?cfg_dir=gp_r4099_169) was used to select TFBS for hot genes (22). TFBS analysis was performed from 5,000 bp to the transcription start site and the significant TFBS was selected with P<0.01. Finally, the result was visualized by Cytoscape software.

Results

Identified DEGs and significant pathways

Based on the aforementioned analysis, a total of 2,681 DEGs between 14 CCSK and 3 FK samples were identified, including 543 upregulated and 2,138 downregulated genes. Furthermore, a total of 33 overexpressed genes, including wingless-type MMTV integration site family, member 5A and fibroblast growth factor receptor 1, were demonstrated to be significantly enriched in the cancer-associated pathways (P=1.35×10−6). Additionally, a total of 50 downregulated genes, including VEGFA and p21-activated kinase 4, were significantly enriched in focal adhesion (P=1.03×10−5) (Table I).

Table I. Most significant pathways in which the differentially-expressed genes were enriched (P<0.05).

Table I.

Most significant pathways in which the differentially-expressed genes were enriched (P<0.05).

Term	Count	P-value	Genes	Gene regulation
hsa05200: Pathways in cancer	33	1.35×10−6	WNT5A, FGFR1, FGF18, WNT5B, FGF9, PDGFA, FOXO1, FGF13, ZBTB16, KIT, GLI2, GLI1, AKT1, IGF1R, WNT4, CASP3, CASP8, RUNX1, CEBPA, EGFR, AR, BMP2, TCF7, CTBP2, EPAS1, CDK6, CDKN1A, CCND1, CDKN1B, PDGFRA, PTCH1, WNT11, PIAS1	Upregulated
hsa04510: Focal adhesion	50	1.03×10−5	PAK4, ILK, PDGFC, AKT3, ACTN4, ROCK1, ROCK2, FLNB, VASP, PRKCB, PPP1CA, LAMC3, VEGFA MAPK3, RELN, LAMC1, COL1A1, PARVB, PARVA, IBSP, CAV2, ERBB2, TNC, ITGB1, PXN, MYL9, LAMB2, DOCK1, ITGAV, SOS2, COL6A3, THBS1, SPP1, COL4A4, VAV3, MET, ITGA2, MYLPF, MYL12B, ITGA3, MYL12A, ITGA4, CAPN2, COL4A6, KDR, LAMA1, LAMA5, ITGA8, ITGA7, MYLK	Downregulated

[i] Count, number of genes enriched in the pathway.

PPI network

Combined with the STRING database, PPI analysis was performed on the 33 upregulated and 50 downregulated genes enriched in the two pathways. The PPI network is presented in Fig. 1, with five hot genes screened out, including three overexpressed genes [VEGFA, integrin αV (ITGαV) and integrin β1 (ITGβ1)] and two downregulated genes [v-akt murine thymoma viral oncogene homolog 1 (AKT1) and EGFR].

Figure 1. PPI network of the differentially-expressed genes enriched in the two significant pathways. Diamonds indicate genes associated with vascular endothelial growth factor A (VEGFA); circles indicate genes not associated with VEGFA; red indicates upregulated genes; green indicates downregulated genes; and larger diamonds indicate hot genes with >30 related genes.

GO analysis for VEGFA-associated genes

As presented in Table II, the VEGFA-associated upregulated genes [cyclin D1 and cyclin-dependent kinase inhibitor 1B (CDKN1B)] were significantly enriched in three MF terms, including cyclin-dependent protein kinase regulator activity (P=0.0137), protein kinase regulator activity (P=0.0408) and kinase regulator activity (P=0.0408). Furthermore, downregulated genes (receptor tyrosine-protein kinase erbB-2 (ERBB2), mesenchymal-epithelial transition (MET) tyrosine kinase receptor and kinase insert domain receptor (KDR)] markedly disturbed two MF terms and three BP terms. The most significant MF and BP terms were transmembrane receptor protein tyrosine kinase activity (P=3.39×10−4) and transmembrane receptor protein tyrosine kinase signaling pathway (P=7.93×10−4).

Table II. Significant GO terms for vascular endothelial growth factor A-related genes (P<0.05).

Table II.

Significant GO terms for vascular endothelial growth factor A-related genes (P<0.05).

Category	Term	P-value	Genes	Regulation
GOTERM_MF_FAT	GO:0016538 cyclin-dependent protein kinase regulator activity	1.37×10−2	CCND1, CDKN1B	Upregulated
GOTERM_MF_FAT	GO:0019887 protein kinase regulator activity	4.08×10−2	CCND1, CDKN1B	Upregulated
GOTERM_MF_FAT	GO:0019207 kinase regulator activity	4.08×10−2	CCND1, CDKN1B	Upregulated
GOTERM_MF_FAT	GO:0004714 transmembrane receptor protein tyrosine kinase activity	3.39×10−4	ERBB2, MET, KDR	Downregulated
GOTERM_MF_FAT	GO:0004713 protein tyrosine kinase activity	1.48×10−3	ERBB2, MET, KDR	Downregulated
GOTERM_BP_FAT	GO:0007169 transmembrane receptor protein tyrosine kinase signaling pathway	7.93×10−4	ERBB2, MET, KDR	Downregulated
GOTERM_BP_FAT	GO:0007167 enzyme-linked receptor protein signaling pathway	1.71×10−3	ERBB2, MET, KDR	Downregulated
GOTERM_BP_FAT	GO:0006928 cell motion	5.01×10−3	ERBB2, VASP, KDR	Downregulated

[i] MF, molecular function; BP, biological process; GO, gene ontology.

Results of TFBS analysis

TFBS analysis was performed on the four hot genes of VEGFA, ITGβ1, EGFR and AKTI. The results are presented in Fig. 2. VEGFA was demonstrated to be regulated by six TFs, including ETS domain-containing transcription factor (ERF), erythromycin resistance methylase (ERM), GA binding protein subunit α (GABPα), norepinephrine transporter (Net), NF-κB and Sp2 transcription factor. Similarly, the upregulated gene, AKT1, was also regulated by the same six TFs.

Figure 2. Transcription factor binding site analysis for hot genes. Diamonds indicate overexpressed TF names; red circles indicate upregulated genes; green circles indicate downregulated genes and blue lines indicate interactions between transcription factors and hot genes.

Discussion

CCSK was originally described in 1970, and often develops in children <3 years (23). CCSK has presented significant diagnostic and clinical challenges, thus numerous studies have aimed to investigate the mechanism of CCSK and also to identify markers for the diagnosis and treatment of the disease. In the present study, the underlying pathogenesis of CCSK was predicted through the use of bioinformatic approaches, with novel markers singled out for use in therapeutics.

A total of 2,681 DEGs, between 14 CCSK and 3 FK samples, were identified. In these DEGs, 2,138 genes were identified to be downregulated in CCSK samples, including VEGFA, whilst 543 were overexpressed. Such results suggest that these dysregulations may induce the development of CCSK, particularly in the downregulated genes. Furthermore, it was demonstrated that 33 upregulated genes were significantly enriched in cancer-associated pathways (P=1.35×10−6), whilst 50 downregulated genes, including VEGFA, were enriched in focal adhesion pathways (P=1.03×10−5). These findings indicate that the dysregulation of these genes may induce CCSK by disturbing these two pathways. Previous studies have demonstrated that these pathways serve important roles in cancer development and progression (24,25), which is consistent with the results from the current study.

In addition, a PPI network was constructed and hot genes were identified. By using DEGs enriched in the two most significant pathways and the STRING database, the PPI network was established and five hot genes were singled out. In the network, numerous genes interacted with VEGFA, which included the two upregulated hot genes AKT1 and EGFR, and the two downregulated genes ITGβ1 and ITGαV. VEGFA is a member of the platelet-derived growth factor/VEGF family and encodes a protein that is typically identified as a disulfide-linked homodimer; it is also an angiogenic factor that is required for tumor growth and metastasis (26). The protein it encodes is a glycosylated mitogen, specifically acting on endothelial cells, with its effects including increased vascular permeability, the promotion of cell migration, the induction of angiogenesis, vasculogenesis and endothelial cell growth, and the inhibition of apoptosis. EGFR, a growth factor receptor tyrosine kinase, has been detected as a common component of multiple types of cancer (27). AKT1 has been widely studied and is demonstrated to be frequently overexpressed in human cancer, with the ability to accelerate cell proliferation and suppress programmed cell death (28). Hunt et al (29) reported that the aberrant activation and expression of ITGβ1 has been identified in a number of malignancies. ITαV is associated with numerous types of cancer, including colon, gastric and breast cancer (30). Based on previous studies, these five hot genes have been widely investigated and detected in multiple tumors. Notably, in the present study, the dysregulation of these genes was also identified in CCSK, suggesting that the five hot genes serve important roles in CCSK development.

Furthermore, the functions of VEGFA-associated genes were analyzed by GO. The two upregulated genes, CCND1 and CDKN1B, were involved in kinase regulation, whilst the downregulated genes, ERBB2, MET and KDR, were significantly enriched in protein tyrosine kinase regulation. CCND1 is involved in several human malignancies, including breast, colorectal and non-small cell lung cancer (31). Variant genotypes of CDKN1B, encoding p27 KIP1, are associated with an increased risk of breast cancer, prostate carcinoma and human hepatocellular carcinoma (32). Overexpression of ERBB2, a receptor-like tyrosine kinase, has been detected in several types of human carcinomas (33). The MET tyrosine kinase receptor promotes tissue remodeling, which subsequently induces organ homeostasis, cancer metastasis and wound repair (34). Jang et al (35) demonstrated that VEGF and KDR serve important roles in angiogenesis, contributing to the development and progression of numerous tumors. The results of the present study demonstrated that these genes were associated with CCSK development, which is consistent with previous studies.

Through TFBS analysis, six TFs that regulated VEGFA, including ERF, ERM, GABPα, Net, NF-κB and SP2, were identified. These six TFs have the ability to regulate VEGFA expression, indicating that the conditions of these TFs may induce the upregulation, downregulation or normal expression of VEGFA, with the downregulated VEGFA believed to subsequently induce CCSK. Boo et al (8) reported the expression of VEGF in CCSK tumor cells, indicating that angiogenesis in CCSK may involve VEGF. However, the associations between VEGF and TFs in CCSK have not yet been investigated. Therefore, to the best of our knowledge, the present study predicts that these TFs are important in CCSK development, diagnosis and treatment for the first time.

In conclusion, the current study investigated the underlying mechanism of CCSK and identified novel diagnostic markers for the disease. Downregulated VEGFA and its interacting genes may initiate CCSK development by interrupting two different pathways, kinase regulation and protein tyrosine kinase metabolism. Additionally, VEGFA, ITGB1, ATK1 and EGFR would be the diagnostic markers for CCSK, while ERF, ERM, GABPα, Net, NF-κB and SP2 may function as targets for CCSK therapy. The present study focused on an original insight regarding the pathogenesis of CCSK and identified novel markers, as well as targets for the disease. However, due to the disadvantages and limitations of bioinformatic approaches, further investigation is required to verify such results.

References