Systematic module approach identifies altered genes and pathways in four types of ovarian cancer.

The present study aimed to identify altered genes and pathways associated with four histotypes of ovarian cancer, according to the systematic tracking of dysregulated modules of reweighted protein‑protein interaction (PPI) networks. Firstly, the PPI network and gene expression data were initially integrated to infer and reweight normal ovarian and four types of ovarian cancer (endometrioid, serous, mucinous and clear cell carcinoma) PPI networks based on Spearman's correlation coefficient. Secondly, modules in the PPI network were mined using a clique‑merging algorithm and the differential modules were identified through maximum weight bipartite matching. Finally, the gene compositions in the altered modules were analyzed, and pathway functional enrichment analyses for disrupted module genes were performed. In five conditional‑specific networks, universal alterations in gene correlations were revealed, which leads to the differential correlation density among disrupted module pairs. The analyses revealed 28, 133, 139 and 33 altered modules in endometrioid, serous, mucinous and clear cell carcinoma, respectively. Gene composition analyses of the disrupted modules revealed five common genes (mitogen‑activated protein kinase 1, phosphoinositide 3‑kinase‑encoding catalytic 110‑KDα, AKT serine/threonine kinase 1, cyclin D1 and tumor protein P53) across the four subtypes of ovarian cancer. In addition, pathway enrichment analysis confirmed one common pathway (pathways in cancer), in the four histotypes. This systematic module approach successfully identified altered genes and pathways in the four types of ovarian cancer. The extensive differences of gene correlations result in dysfunctional modules, and the coordinated disruption of these modules contributes to the development and progression of ovarian cancer.

Introduction

Among women, ovarian cancer has a high mortality rate and is the fifth leading cause of cancer-associated mortality, behind cancer of the lung and bronchus, and colorectal, breast and pancreatic cancer (1). The average age of ovarian cancer onset occurs later in reproductive life (2). The disease can advance rapidly with transcoelomic spreading from the ovary to other organs and peritoneal surfaces, and with ascites accumulation (3). Of cases of primary ovarian cancer, ~90% are epithelial carcinoma from the ovarian surface epithelium (4,5). Ovarian epithelial cancer predominantly contains four histotypes of epithelial tumor, including endometrioid, serous, mucinous and clear cell carcinoma. The serous type is the predominant form in women (3). It is well known that early-stage ovarian cancer (stage I/II) is difficult to diagnose, as it is frequently asymptomatic. Therefore, the majority of patients suffering from ovarian cancer are in advanced stages (III and IV) at the time of the initial diagnosis (6,7). To date, the treatment of ovarian cancer is primarily via platinum-based chemotherapy, debulking surgery and radiotherapy, however, the five-year survival rate has only improved only marginally in the last 40 years, remaining <40% (8). Thus, it is critical that effective and sensitive diagnostic biomarkers are examined, which can be applied in the early stage of ovarian cancer and improve survival rates of patients.

Ovarian carcinogenesis is caused by accumulated genetic or genomic alterations (9). DNA-microarray technology enables examination of the expression of thousands of genes simultaneously in tumor samples. Data-analysis software, a high-throughput technology, has made it possible to distinguish gene expression profiling between normal and cancer samples, and thus identify differentially expressed genes during cancer development and progression (10). Gene expression profiling can provide information for the mining of novel biomarkers. Based on oligonulceotide arrays, 275 genes have been predicted with increased/decreased expression in ovarian cancer (11). Several characteristic biomarkers involved in ovarian cancer have been determined, including E-cadherin (12), carbohydrate antigen-125 (13), cytochrome P450 1B1 (14), cyclooxygenase 1 (15), AKT serine/threonine kinase 2 (16), BRCA (17) and the human epidermal growth factor receptor family (16,18). The dual detection of hepatocyte nuclear factor-1β and napsin A have been reported to be sensitive markers for diagnosing ovarian clear cell carcinoma, which may also be useful for distinguishing ovarian clear cell carcinoma from endometrioid, serous carcinoma and metastatic Krukenberg tumors (19). However, the compensatory mechanisms in four histotypes of ovarian cancer remain to be fully elucidated, and the underlying diagnostic and therapeutic targets require further investigation.

It is known that protein complexes are key molecular entities and they integrate multiple gene products to perform cellular functions (20). Based on advances in high-throughput analysis technologies, substantial protein-protein interaction (PPI) data has been excavated, therefore, it is possible to investigate protein functions systematically (21). However, due to the technological limitations and dynamic nature of protein interaction maps, the protein interaction data produced by high-throughput experiments often possess high false positive and false negative rates, which lead to difficulties in predicting protein complexes accurately (22). Therefore, a systematic method is required to track gene and module behavior across diseases conditions in a controlled manner (23).

The present study aimed to further elucidate the mechanisms of four histotypes of ovarian cancer, therefore, the disrupted modules from reweighted PPI networks were tracked to systematically identify dysfunctional genes and pathways in samples of the four histotypes of ovarian cancer. Initially, based on Spearman's correlation coefficient (SCC) of gene interactions, normal and disease-specific PPI networks were inferred. Subsequently, the clique-merging algorithm was used to examine modules in the re-weighted PPI networks, and the modules obtained in cancer were compared with those in the normal condition to determine altered modules. Finally, the associated functional pathways of the different histotypes were identified, based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) database.

Materials and methods

Affymetrix chip data

The EMBL-EBI (http://www.ebi.ac.uk/) database provides freely available data from life science experiments, performs basic investigations in computational biology and offers a training program to extensive users. The E-GEOD-6008 dataset, which included 99 individual ovarian tumor samples (37 endometrioid, 41 serous, 13 mucinous and eight clear cell carcinomas) and four normal ovarian samples, was downloaded from EMBL-EBI (24). The RNA expression in each sample was analyzed using an Affymetrix GeneChip Human Genome HG-U133A array.

Data preprocessing

The Affy package (v1.48.0) in R (bioconductor.org/packages/release/bioc/html/affy.html) was used for the analysis of oligonucleotide arrays to delete undesirable values (25). The data preprocessing was performed using the robust multichip average method in the affy package, comprising background correction, normalization, perfect match/mismatch matching and expression value aggregate calculation (25). Subsequently, the raw data in the CEL files were converted to probe-level data, and the probe-level data were transformed to gene symbols. The FeatureFilter function was applied to discard probes that did not correspond to any gene symbol. A final total of 12,493 genes were obtained.

PPI network construction

As is already known, proteins rarely exert their functions individually, however, they are important in a variety of biological process in the form of large protein functional groups (26). Therefore, in the present study, the PPIs of the 12,493 genes were analyzed using the online Search Tool for the Retrieval of Interacting Genes (STRING) tool (string-db.org). Cytoscape software (v3.3.0; www.cytoscape.org), a biological graph visualization tool (27), was used to construct the PPI networks. All 1,048,576 interactions datasets were downloaded from the STRING database to construct the PPI networks. Following the elimination of self-loops, a complicated PPI network was constructed, which comprised 9,273 nodes and 58,617 interactions with a combine-score ≥0.75.

PPI network re-weighting

By obtaining intersection elements of the 12,493 genes in E-GEOD-6008 and the 9,273 nodes in the PPI network, a sub-network of 7,264 nodes and 45,286 interactions was obtained. The weights of interactions reflect their reliabilities, and low absolute scores of interactions may indicate false positives (28). In the present study, the SCC, which describes the association between two variables, was used to evaluate the strength of the association between two paired proteins in the PPI networks. The SCC value ranged from −1 to +1. The sign of the SCC indicates the direction of association between X (the independent variable) and Y (the dependent variable). If Y increases when X increases, the SCC is positive. If Y decreases when X increases, the SCC is negative. A coefficient of −1 indicates that there is a perfect inverse association between X and Y. A coefficient of +1 demonstrates that there is a perfect positive association between X and Y. An SCC of 0 indicates that there is no tendency for Y to either increase or decrease when X increases. The SCC increases in magnitude as X and Y become closer to being perfect monotone functions of each other. When X and Y are perfectly monotonically associated, the SCC is 1.

Spearman's rank formula was used to calculate the coefficient of two paired proteins X and Y in the PPI network. The formula was as follows:

where ‘R’ is the coefficient, ‘d’ is the difference between the ranks of corresponding values X and Y, ‘sum(d2)’ is the total of the ‘d2’ column, and ‘n’ is the number of observations. In the present study, the SCC of a gene-gene interaction was defined as the weight value of the interaction.

Module identification

Similar to the method described by Liu et al (28), the module-identification algorithm was performed in three steps based on clique-merging. Firstly, all of the maximal cliques from the weighted PPI networks of the normal sample and four histotypes were selected out, respectively. The maximal cliques were enumerated utilizing a fast depth-first method with a pruning-based algorithm, described by Tomita et al (29). Subsequently, a score was assigned to each clique, and the clique score (C) was referred to as its weighted density:

where w (u,v) represents the interaction weight between u and v, based on the fast depth-first method.

In the third step, the cliques were arrayed in descending order based on their weighted density and the highly overlapped cliques were removed in order to reduce the size of the result. The highly overlapped cliques were merged to construct larger, dense sub-graphs. The term ‘inter-connectivity’ was utilized to confirm whether the two overlapped cliques be merged together or not. The inter-connectivity score (C1,C2) between the non-overlapping proteins of C1 and C2 was computed according to the following formula:

The obtained clique scores were ranked in descending order and denoted as {C1,C2,…,Ck}. For every maximal clique Ci, if there existed another maximal clique Cj, and Cj possessed a lower score than Ci and |Ci∩Cj|/|Cj|≥to (a predefined overlap-threshold), the weighted inter-connecting score was calculated for the distinct nodes between the two cliques. Provided that Cj existed, the interconnectivity score (Ci Cj) was used as a standard to determine whether to remove Cj or merge Cj with Ci. If the inter-score (Ci, Cj) was higher than or equal to the predefined merge-threshold tm, Cj and Ci were merged to obtain a module; if not, Cj was removed. In the present study, the overlap-threshold was 0.5 and the merge-threshold was 0.25.

Differential module identification

Random statistical analysis was performed for the obtained modules in the four subtypes of ovarian cancer. A p-value cutoff of 0.01 can reveal numerous false positive results and requires another factor. False discovery rate (FDR), one of the most widely used multiple testing criterions for controlling errors of false discoveries, was utilized to adjust the P-value obtained via random statistical analysis. The FDR was first defined by Benjamini et al (30) as the expected proportion of the number of falsely rejected hypotheses among the total number of rejected hypotheses. In the present study, modules with an adjusted P<0.01 based on the FDR measure were considered to be disrupted modules.

Pathway enrichment analysis of genes in disrupted modules

The Database for Annotation, Visualization and Integrated Discovery (DAVID; david.ncifcrf.gov) is a comprehensive functional annotation software program, and is used for integrative and systematic analysis of large gene groups (31). In the present study, KEGG (www.genome.jp/kegg) pathway enrichment analysis was performed using DAVID with the threshold of FDR-adjusted P<0.001 for genes from the altered modules of endometrioid, serous, mucinous and clear cell carcinomas samples, respectively. For the enriched pathways, the appearance frequency of every gene was counted. The higher a frequency of a gene, the higher its level of involvement in pathways and the higher its importance.

Results

Disruptions in the PPI networks of four types of ovarian cancer

A total of 12,493 genes were obtained from normal ovarian sample and four ovarian cancer (endometrioid, serous, mucinous and clear cell carcinomas) samples using a data preprocessing procedure. The normal ovarian and four ovarian cancer PPI networks reflected equal numbers of nodes (7,264) and interactions (45,286). Subsequently, re-weighted PPI networks of the normal ovarian sample and the four stages of disease were examined using the SCC algorithm. In the normal ovarian and four ovarian cancer networks, the numbers of interactions and average scores (weights) were approximately equal; the 45,286 interactions had average scores of 0.083 (normal), 0.090 (endometrioid carcinoma), 0.085 (serous carcinoma), 0.087 (mucinous carcinoma) and 0.070 (clear cell carcinoma). The correlationwise frequency distributions were different across the normal ovarian and four ovarian cancer networks (Fig. 1). When the expression correlations varied between −1.0 and −0.8, −0.6 and −0.4, and 0.6 and 1.0, the number of interactions in the normal ovarian network was higher, compared with those in the four types of ovarian cancer. When the expression correlations varied between −0.4 and 0.6, the number of interactions in the normal network was almost lower, compared with those in the four types of ovarian cancer. In addition, the scores of the total 26,651 interactions in the four ovarian cancer networks were lower, compared with that in the normal network, whereas the total numbers (18,635) of interactions were higher in the disease conditions, compared with that in the normal condition.

Figure 1.

Correlation distribution of interactions in normal ovarian tissue and four types (endometrioid, serous, mucinous and clear cell) of ovarian cancer-specific networks based on Spearman's correlation coefficient.

Disruptions in the four ovarian cancer modules

The disrupted or altered modules from the normal and four ovarian cancer PPI sub-networks were identified based on the clique-merging algorithm. With the node threshold >5, a total of 951 modules were obtained under the five conditions. Comparative analysis for normal and disease modules was then performed to further elucidate the disruptions from a module perspective. Notably, as shown in Table I, the total number of modules (951) and average module size (31.83) were the same across the five conditions, which results from the same interactions. Furthermore, the average weighted density of mucinous carcinoma was marginally higher, compared with that of the other three cancer subtypes. The associations between the numbers of modules and weighted correlation density of the modules are shown in Fig. 2. No significant difference was found between the distribution of modules in the normal and disease conditions at the level of the overall correlation distribution based on the Kolmogorov-Smirnov test (P>0.05).

Table I.

Properties of the normal ovarian, and endometrioid, serous, mucinous and clear cell carcinoma modules.

			Correlation
Module set	Number of modules	Mean module size	Maximum	Average	Minimum
Normal	951	31.83	0.543	0.067	−0.021
Endometrioid	951	31.83	0.442	0.058	−0.055
Serous	951	31.83	0.693	0.058	−0.158
Nucinous	951	31.83	0.551	0.072	−0.088
Clear cell	951	31.83	0.396	0.049	−0.167

Figure 2.

Correlation in distribution of modules in normal and cancerous (endometrioid, serous, mucinous and clear cell carcinoma) ovarian tissues based on module correlation density.

Identification of differential modules

In the present study, a total of 28, 133, 139 and 33 differential modules (FDR-adjusted P<0.01) were identified in the endometrioid, serous, mucinous and clear cell carcinoma, respectively. Extracting genes from the differential modules in the four types of ovarian cancer resulted in 533, 491, 591 and 408 genes, respectively.

Pathway enrichment analysis of genes in differential modules

DAVID-based KEGG functional pathway analysis was performed for the genes involved in the differential modules. Based on the FDR-adjusted P<0.001, a total of 12, 16, 18 and nine significantly enriched pathways were identified in endometrioid, serous, mucinous, and clear cell carcinoma, respectively. As shown in Fig. 3, pathways in cancer was a common pathway across the four subtypes of ovarian cancer, which may be vital in the development of ovarian cancer. The cell cycle, DNA replication and ribosome pathways were also associated with the subtypes of ovarian cancer, with the exception of endometrioid cancer. In addition, unlike the endometrioid, serous and mucinous carcinoma, clear cell ovarian cancer was associated with several metabolic pathways, including drug metabolism and retinol metabolism. Overall, the cancer-associated pathways were important in the development of ovarian cancer.

Figure 3.

Significant pathways of genes involved in the differential modules in the four types of ovarian cancer (endometrioid, serous, mucinous and clear cell carcinoma) based on P<0.001. Dark green color indicates a present pathway; light green color indicates an absent pathway.

By identifying the intersection of the common genes across the four ovarian cancer subtypes and the most frequent genes in the disrupted pathways, five key genes were obtained: Mitogen-activated protein kinase 1 (MAPK1), phosphoinositide 3-kinase-encoding catalytic α (PIK3CA), AKT serine/threonine kinase 1 (AKT1), cyclin D1 (CCND1) and tumor protein P53 (TP53), which may perform an essential function in the pathogenesis of the four subtypes of ovarian cancer (Table II).

Table II.

Total 15 genes appearing most frequently in the disrupted pathways of the four types of ovarian carcinoma.

Endometrioid		Serous		Mucinous		Clear cell
Gene	Frequency	Gene	Frequency	Gene	Frequency	Gene	Frequency
MAPK1	20	MAPK1	24	MAPK1	24	MAPK1	20
PIK3CA	16	EGFR	16	PIK3CD	20	PIK3CA	18
EGFR	15	TP53	15	PIK3CA	20	AKT1	17
PIK3R3	15	CCND1	14	PIK3R2	19	CCND1	13
PIK3R2	15	AKT1	13	AKT1	18	TP53	12
AKT1	14	PRKCG	12	GRB2	16	PRKCG	11
TP53	12	CDK4	11	CCND1	16	PLCG1	11
CCND1	12	EGF	11	SOS1	16	PLCG2	10
PRKCG	10	E2F1	10	EGFR	15	RELA	9
PLCB3	10	E2F3	10	TP53	15	PLCB3	9
PLCB1	10	RB1	10	MYC	12	NFKB1	9
PLCB2	10	MYC	10	CDK4	11	CDK4	9
PTEN	9	PLCB3	9	PLCG1	11	PLCB1	9
ERBB2	8	PTEN	9	E2F1	10	ADCY9	9
MET	8	IGF1R	9	E2F2	10	PLCB2	9

MAPK1, mitogen-activated protein kinase 1; PIK3CA, phosphoinositide 3-kinase-encoding catalytic α; AKT1, AKT serine/threonine kinase 1; CCND1, cyclin D1; TP53, tumor protein P53; PLCB1, phospholipase C β1; PRKCG, protein kinase Cγ; GRB2, growth factor receptor-bound protein 2; CDK4, cyclin-dependent kinase 4; PLCG, phospholipase C γ; EGF, epidermal growth factor; SOS1, son of sevenless homolog 1; EGFR, EGF receptor; RB1, retinoblastoma 1; NFKB1, nuclear factor κB1; ERBB2, Erb-B2 receptor tyrosine kinase 2; PTEN, phosphatase and tensin homolog; ADCY9, adenylyl cyclase 9; IGF1R, insulin-like growth factor 1.

Discussion

The aim of the present study was to identify dysregulated genes and pathways in four histotypes of ovarian cancer via systematically tracking the dysregulated modules of reweighted PPI networks. The reweighted PPI networks of the normal and four ovarian cancer histotypes were obtained based on the SCC, and the modules in the PPI networks were identified. By comparing the modules of the normal and four ovarian cancer histotypes, 28, 133, 139 and 33 disrupted modules were obtained for endometrioid, serous, mucinous and clear cell carcinoma, respectively. A total of five common genes (MAPK1, PIK3CA, AKT1, CCND1 and TP53) and one common pathway (pathways in cancer) across the four histotypes were examined based on gene composition and pathway enrichment analyses.

The pathways in cancer pathway covers several types of pathway involved in cancer. Chen et al (32) documented that PI3K/AKT/hypoxia-inducible factor-1α/CCND1 pathway is vital in follicle-stimulating hormone-driven ovarian cancer cell proliferation. Genistein suppresses the epithelial-mesenchymal transition and migration efficacies of ovarian cancer cells via the estrogen receptor pathway and downregulation of the transforming growth factor-β signaling pathway (33). The extracellular-signal-regulated kinase (ERK) and c-Jun N-terminal kinase (JNK) signaling pathways can be regulated by interleukin-33, and promote ovarian cancer growth and metastasis (34). In addition, activation of the mammalian target of rapamycin (mTOR) signaling pathway has been demonstrated to promote epithelial ovarian cancer metastasis (35).

MAPKs, a family of serine/threonine protein kinases, including p38 MAPK, ERK1/2, and stress-activated protein kinases (JNK). MAPKs are mediators of comprehensive cellular programs, including cell proliferation, cell differentiation and cell apoptosis, in response to distinct stimuli (36). Studies have found that MAPKs are transiently triggered during mitosis and MAPK activation is involved in the spindle assembly checkpoint (37). Consistently, the protein levels of MAPK1 are increased following demecolcine treatment (38). The proto-oncoprotein, Mos, a serine/threonine kinase, has been recognized as a potent activator of MAPK1 during oocyte maturation (39,40). The overexpression of microRNA-378a-3p or silencing of MAPK1 can reduce the expression level of MAPK1 and enhance adipogenesis (41).

The PIK3CA protein modulates various signals to restrain apoptosis and facilitate cell survival and proliferation in several types of cell (42,43). It has been demonstrated that oncogenic mutations and amplification of PIK3CA can activate the PI3 K/Akt signaling pathway to initiate human papillomavirus-induced tumorigenesis and other types of cancer (44–46). Akt is a serine/threonine protein kinase comprising Akt1, Akt2 and Akt3. Akt1 encodes the principal Akt isoform associated with apoptosis regulation (47). In oropharyngeal cancer, the high protein level of Akt can be an unfavorable prognostic biomarker for relapse-free survival rates in patients (48). The importance of the PI3K/AKT signaling pathway in ovarian cancer has been well documented. In general, this pathway has significant roles in gene transcription, protein synthesis and membrane trafficking, however, the abnormal triggering of this pathway leads to cancer initiation, progression and invasion (49,50). PIK3CA and AKT1 amplification are regarded as prognostic factors for ovarian cancer, and the PI3K/AKT/mTOR axis may become a target for drugs.

Another oncogene, CCND1 is a dominating driver of several types of human tumor, including squamous cell and breast cancer, myeloma and Bcell lymphoma (51,52). A previous study indicated that CCND1 was overexpressed in >50% of human breast cancer cases (51). In mice, mammary-targeted gene overexpression resulted in mammary tumorigenesis (53). It has been demonstrated that jumonji and AT-rich interaction domain containing 2 can significantly inhibit leukemia cell proliferation by downregulating the expression of CCND1. In addition, the overexpression of CCND1 is closely associated with low-grade ovarian cancer, which is in line with the suggestion that CCND1 is a downstream target for the active MAPK constitutively expressed in ovarian tumors (54,55).

TP53 is a critical transcriptional regulator, which is involved in cell cycle and cell apoptosis upon activation by oncogenes and DNA damage (56). The activated TP53 protein is combined with the regulatory region of target genes to initiate the cell cycle (57). TP53 mutations are frequently screened genetic alterations in ovarian cancer (58). Reles et al (59) reported that the TP53 alteration closely correlates with poor response to chemotherapy, early recurrence and shortened survival rates in patients with ovarian cancer. The high prevalence of TP53 mutations in tubal epithelial carcinoma shows that the TP53 mutations occur in early carcinogenesis. Thus, TP53 mutations are considered to be poor prognostic factors (60). However, TP53 has been confirmed as an effective blood-based biomarker for the detection of ovarian cancer (61).

In conclusion, the present study successfully identified disrupted modules, including the pathways in cancer module, and hub genes (MAPK1, PIK3CA, AKT1, CCND1 and TP53) in four types of ovarian cancer based on the integrated PPI network. It was inferred that these pathways and genes may be potential biological processes and markers for understanding the mechanism underlying ovarian cancer.