Identification of modules and hub genes associated with platinum-based chemotherapy resistance and treatment response in ovarian cancer by weighted gene co-expression network analysis.

High-grade serous ovarian carcinoma (HGSOC) is the most prevalent and malignant ovarian tumor.To identify co-expression modules and hub genes correlated with platinum-based chemotherapy resistant and sensitive HGSOC, we performed weighted gene co-expression network analysis (WGCNA) on microarray data of HGSOC with 12 resistant samples and 16 sensitive samples of GSE51373 dataset.A total of 5122 genes were included in WGCNA, and 16 modules were identified. Module-trait analysis identified that the module salmon (cor = 0.50), magenta (cor = 0.49), and black (cor = 0.45) were discovered associated with chemotherapy resistant, and the significance for these platinum-resistant modules were validated in the GSE63885 dataset. Given that the black module was validated to be the most related one, hub genes of this module, alcohol dehydrogenase 1B, cadherin 11, and vestigial like family member 3were revealed to be expressional related with platinum resistance, and could serve as prognostic markers for ovarian cancer.Our analysis might provide insight for molecular mechanisms of platinum-based chemotherapy resistance and treatment response in ovarian cancer.

Introduction

Ovarian cancer (OV) is one of the commonest gynecologic malignancies, which causes many deaths each year. Epithelial ovarian cancer (EOC) comprises more than 85% OV, and could be further classified into 5 categories, including clear cell, mucinous, endometrioid, low-grade, and high-grade serous ovarian carcinoma (HGSOC).[ HGSOC accounts for nearly 70% of EOC,[ and is always diagnosed at advanced stages, which makes HGSOC the most prevalent and malignant ovarian carcinoma. Conventionally, EOCs will be treated with cytoreductive surgery followed by chemotherapy, such as platinum and taxane. However, in contrast to the 4 other types of EOCs, which exhibit low degree of response to the mentioned therapy, patients with HGSOC initially respond well to the treatment, but subsequently recur and become chemoresistant at high proportion.[ And the outcomes of HGSOC are poor, with 5-year survival being less than 30%.[ The progression of resistance seriously restrained the efficacy of platinum-based chemotherapy. It is crucial to reveal the underlying mechanisms of platinum resistance and uncover novel therapeutic strategies to improve treatment outcomes of patients with HGSOC.

HGSOC is thought to be a complicated trait determined by a series of mutations occurring in thousands of cancer-associated genes.[ However, mechanistically these genes do not all contribute in the same way to HGSOC progression. Transcriptomic sequencing technology, such as microarray and RNA_seq, is an effective choice for detecting potential participants of complicated traits on a genome-wide scale in animals, fungus, plants, and microorganisms.[ Thousands of analyses about cancers,[ immune systems,[ and genetic diseases[ in model animals have been widely conducted using sequencing technologies. The integrations of genes’ spatio-temporal expression patterns and HGSOC related traits have helped to identify the potential genes and controlling machineries of ovarian carcinomas.[ The accumulation of these data provides us great resource to investigate the underlying molecular mechanisms, pathways, and networks of disease including HGSOC.

Gene co-expression networks based on sequencing data facilitate a global view of gene-gene co-expressed relationships.[ Genes involved in related biological processes (BPs) tend to be co-expressed and clustered as functional modules which can help to discover how the interplay between inter-connected genes accomplish specific biological functions.[ Weighted gene co-expression network analysis (WGCNA) is an efficient method widely used in systematic biological analysis, which explores co-expression modules in correlation with traits or phenotypes based on gene expression data.[ By depicting relevant genes in the samples, disclosing the highly correlated genes in the module, summarizing such clusters using the module eigengene or an intramodular hub gene, and relating modules to the external sample traits, WGCNA can discover the genes for unknown function with BPs and identify candidate biomarkers or therapeutic targets.[

In this study, we performed WGCNA based on data from microarray of chemotherapy resistant and sensitive HGSOC to identify critical co-expression modules and calculated the module-trait relationship, functionally annotated the trait related modules, and determined hub genes related with chemotherapy resistance and sensitivity, which might provide new insights for exploring the underlying molecular mechanisms of HGSOC.

Materials and methods

Data collection

The gene chip GSE51373[ (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE51373) and GSE63885[ (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE63885) of chemotherapy resistance and treatment response in OV with its clinical manifestation data were downloaded from the gene expression omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). The platform was GPL570 (Affymetrix Human Genome U133 Plus 2.0) for both datasets. In this study, the GSE51373 was used to establish the WGCNA co-expression network, and the GSE63885 served as the validation dataset. There are 12 samples from the resistant cohort, and 16 samples from the sensitive cohort of GSE51373. For the validation dataset of GSE63885, there are 34 samples and 41 samples for platinum-resistance and sensitive cohorts, respectively. The raw data had been processed, and the gene expression matrix provided by the website was directly used for the analysis. The annotation information for microarray analysis was downloaded at the platform GPL570. The ENTREZ gene ID was used for network construction. The probe mapped to multiple ENTREZ gene IDs were omitted, and the median of expression values of probes assigned to 1 gene ID was used to represent the expression of this gene. Microarray data were normalized to the mean or median of each array, as described in the original papers.[ More than 20,000 genes were entered into the BRB-ArrayTools (https://linus.nci.nih.gov/BRB-ArrayTools/) for filtering: genes with the top 25% variance of expression values among samples were used for further analysis.

WGCNA analysis of the filtered genes

The WGCNA package in R (Version 3.3.2) was used to construct a co-expression network for the filtered genes.[ After sample clustering, scale independence, and mean connectivity analysis of modules with different power values were performed to determine the soft threshold of module analysis. The power value was set from 1 to 20, and then the values of scale independence and mean connectivity were generated according to these power values. The power value was determined when the scale independence value was 0.9. A hierarchical clustering dendrogram of the TOM matrix was constructed by the average distance with a minimum size threshold of 30 and the merge cut height of 0.25 to classify the similar genes expression profiles into different gene modules. A cluster dendrogram among modules and an eigengene adjacency heatmap between modules were generated. The co-expression networks were visualized using Cytoscope 2.8.2. The betweenness of modules were calculated by igraph package in R.

Module-trait analysis based on patients’ varying sensitivity to first-line chemotherapy

Information on the platinum-based chemotherapy for treatment of serous EOC of 28 patients, was collected from the original study to identify significant co-expression modules related with the clinical characteristics (as trait). All 28 patients had advanced stage, HGSOC, and were treated with the same standard platinum-based chemotherapy. Twelve patient tumors demonstrating relative resistance to platinum chemotherapy corresponding to shorter progression-free survival (PFS) (<8 months) were compared to 16 tumors from platinum-sensitive patients (PFS >18 months). Module-trait relationships were calculated according to the correlation between modules and traits; modules that were significantly correlated with individual traits (P-value <.05, module size <500) were identified; and genes in significant modules were then exported for further analysis. The differentially expressed genes (DEGs) identified in GSE51373 between resistance to platinum samples and platinum-sensitive samples were screened by using limma package in R with the cutoff of P-value <.05.

Functional enrichment analysis of significant co-expression modules

The functional verified OV marker genes were referenced from the study published by Zhang.[ The package TopGO of R 3.3.2[ was used in enrichment analysis for gene ontology (GO) of the studied genes, including OV marker genes and the gene members of each module. The GO functional annotations file of Homo sapiens was downloaded from GO database (submission date: May, 2018, http://geneontology.org/gene-associations/goa_human_rna.gaf.gz). The DAVID database (https://david.ncifcrf.gov/) was used for gene information retrieving, gene ID conversion and gene function annotation. Pathways with a P-value <.05 were considered as significant pathways.

Validating key modules

The modulePreservation function of WGCNA package in R was used to calculate the preservation Zsummary of each module in the GSE63885 dataset to determine if a module was reproducible. The Experimental design for GSE63885 was similar with that of GSE51373, and the gene expressional microarray data of GSE63885 was retrieved from GEO dataset, which is consisted of 34 platinum-resistant and 41 platinum-sensitive surgical OV samples. When the Zsummary of the certain module was approaching 10 or even higher, and median rank was relatively low, the module was considered to be preserved.

Kaplan–Meier survival analysis

The prognostic value of the hub genes was assessed by relapse-free survival analysis. Log-rank tests were performed and Kaplan–Meier survival curves were plotted based on the online software Kaplan–Meier plotter which integrates several other independent data sets (http://kmplot.com/analysis/index.php?p=service&cancer=ovar). In these analyses, P < .05 was considered statistically significant.

Results

WGCNA co-expression network construction

By selecting genes with the threshold of the top 25% variance of expression values, a total of 5122 genes were included in WGCNA. By sample clustering, no outliers were observed in 28 samples, thus all samples were included in the analysis. Then the soft threshold was determined by scale independence and mean connectivity analysis of modules with different power values ranging from 1 to 20. In our study, power = 3 was set to guarantee high scale independence (near 0.9) and low mean connectivity (near 0) (Fig. 1A). The mergeCutHeight was set as 0.25, and a total of 16 modules were generated and displayed with different colors (Fig. 1B). Of the 16 modules, 1297 and 1115 genes were assigned to grey and turquoise module, blue and brow module included 475 and 447 genes, 51 genes were grouped to salmon module (Fig. 1B; Supplementary Table 1).

Figure 1

Construction of a weighted gene co-expression network. (A) Scale independence and mean connectivity analysis. (B) Cluster dendrogram among modules.

Module-trait relationship calculation

The chemotherapy (resistant or sensitive) result of samples was collected as a trait, and module-trait relationships were calculated according to the correlation between 16 modules and traits. The module trait relationship is shown in Figure 2A. The salmon (cor = 0.50, P-value = 7.00E-03), magenta (cor = 0.49, P-value = 8.00E-03), and black (cor = 0.45, P-value = 2.00E-02) module were discovered associated with chemotherapy resistant (Fig. 2A and Table 1), whereas brown (cor = 0.23) and green (cor = 0.22) module shew weakly related with chemotherapy sensitive with the P-value >.05 (Fig. 2A and Table 1). The eigengene adjacency heatmap is shown in Figure 2B, which indicates that the salmon, magenta and black module were adjacent and the brown, green and yellow module clustered.

Figure 2

Association between modules and platinum-resistance and sensitivity in ovarian cancer. (A) Module-trait relationships. The darker the module color, the more significant their relationship. (B) Eigengene adjacency heatmap between modules.

Table 1

Module-trait relation and gene function enrichment information of genes in resistant or sensitive related modules.

Function annotation of the trait related modules

The results of the GO analysis of the chemotherapy resistant related modules (salmon, magenta, and black) and potential chemotherapy sensitive related modules (brown and green) were shown in Table 1 and Supplementary Table 2. The BP GO terms enriched by OV marker genes were used as OV related GO terms. For the salmon module, we identified 120 GO terms, including OV related terms “positive regulation of phosphatidylinositol 3-kinase (PI3K) activity” (GO:0043552), “positive regulation of phospholipase C activity” (GO:0010863), and “somatic stem cell population maintenance” (GO:0035019). A total of 158 terms were annotated in resistant module magenta, including OV related terms “branching involved in ureteric bud morphogenesis” (GO:0001658), “palate development” (GO:0060021), and “positive regulation of bone mineralization” (GO:0030501). For the sensitive related green module, we identified 201 GO terms, including OV related “embryonic skeletal system morphogenesis” (GO:0048704), “type I interferon signaling pathway” (GO:0060337), and “anterior/posterior pattern specification” (GO:0009952).

Validation for preservation of the modules

To evaluate whether the platinum-resistance related modules, especially salmon, magenta, and black modules are robust and reproducible across datasets, the modulePreservation function was applied to validate the preservation of all modules in the GSE63885 dataset. The results showed that for the black module and magenta module, the Zsummary were more than 10 and the preservation median rank were below 5, and they were considered to be highly preserved in GSE63885. Whereas for salmon module with its Zsummary score of 5.9 and preservation median rank score of 6, was revealed to be relatively less preserved (Fig. 3A and B).

Figure 3

Preservation of modules. (A) The medianRank preservation statistics of the modules. Each module was represented by a point, labeled with the corresponding color and names. Low value of numbers on the y-axis indicate high preservation. (B) The Zsummary preservation statistics of the modules. The modules are labeled as in (A). The dashed blue and green lines indicate the thresholds. Values over 2 of Zsummary represent moderately preserved modules, and values over 10 provides strong evidence of module preservation.

Identification of OV related hub genes in chemotherapy resistant/sensitive modules

The edges signifying the correlations in the 3 chemotherapy resistant related modules and 2 sensitive modules were ranked by the weight value and the top 5% edge were used for screening OV related hub genes. The co-expression relations of OV marker genes and DEGs were obtained from the top 5% edges, and a total of 17/45 (nodes/edges), 47/156, 77/186, 69/335, and 45/63 were obtained from the salmon, magenta, black, brown, and green module (Fig. 4; Supplementary Fig. 1). The degree of the filtered edges was calculated, and genes with top degree values were identified as OV related hub genes in each module and shown in Table 2. For resistant module salmon, KIT proto-oncogene receptor tyrosine kinase (KIT) and tumor growth related fibroblast growth factor 2 (FGF2) co-expressed with 18 OV genes or DEGs, separately (Fig. 4A and B). Aldehyde dehydrogenase 1 family member A1 (ALDH1A1), aldehyde oxidase 1 (AOX1), and adenylate cyclase 8 (ADCY8) are the most important hub genes in the module of magenta (Fig. 4C and D). Alcohol dehydrogenase 1B (ADH1B), cancer blood-based biomarker calbindin 2, and a tumor suppressor caveolin 1 were screened as hubs in resistant module black (Fig. 4E and F). In chemotherapy sensitive module brown, biomarker for diagnosis of OV acid phosphatase, prostate (ACPP) and alkaline phosphatase, germ cell (ALPG) was filtered as hubs (Supplementary Fig. 1A and B). For green module, annexin A4 was revealed to be hub gene (Supplementary Fig. 1C and D).

Figure 4

The co-expression relations of OV marker genes and DEGs were obtained from the top 5% edges in chemotherapy resistant modules. (A) The co-expression network of the significant genes in the salmon module. It has 17 nodes and 45 edges. The triangle nodes indicate the OV marker genes and the blue shape indicate DEGs upregulated in resistant samples. (B) The heatmap of 15 DEGs in chemotherapy resistant and sensitive samples in salmon module. A red color indicates that the gene is highly expressed in resistant patients. (C) The co-expression network of the significant genes in the magenta module. (D) The heatmap of 39 DEGs in chemotherapy resistant and sensitive samples in magenta module. (E) The co-expression network of the significant genes in the black module. (F) The heatmap of 48 DEGs in chemotherapy resistant and sensitive samples in black module. DEGs = differently expressed genes, OV = ovarian cancer.

Table 2

Information of hub genes co-expressed with OV marker genes in resistant or sensitive related modules (ranked by degree).

ADH1B, cadherin 11, and vestigial like family member 3are differentially expression and their prognostic roles for OV were validated in TCGA database

Because the black module was validated by GSE63885 dataset to be the most correlated module with platinum resistance, we checked the expression and prognostic roles of the hub genes to support with clues for further analysis in the aspect of platinum resistance or sensitivity for OV cases. We found that the enriched hubs of ADH1B, cadherin 11 (CDH11), and vestigial like family member 3 (VGLL3) were representative for platinum resistant ones. First, we separated the OV patients into platinum resistant and sensitive cases and compared ADH1B, CDH11, and VGLL3 expression in 2 groups. We found all 3 genes got higher expression in platinum resistant patients compared with that of platinum sensitive ones (Fig. 5A, C, E). Then, to discover the prognostic roles of these genes, the overall survival data revealed that the prognosis of patient with higher expression of ADH1B, CDH11, and VGLL3 is poorer (Fig. 5B, D, F). Besides, several hub genes were also screened to be related with prognosis for patients with OV, including KIT and AQP5 for salmon and brown modules. KIT was revealed to be expressed with higher levels in platinum-resistant compared with platinum-sensitive patients, and higher expression of KIT correspond to poorer survival (Supplementary Fig. 2A and B). Whereas, AQP5 represents the hubs of modules that correlated with platinum-sensitivity. Although AQP5 showed no differences in expression between platinum resistant and sensitive patients, it is still meaningful for predicting patient survival after chemotherapy (Supplementary Fig. 2C and D).

Figure 5

Expressional changes in platinum-resistant ovarian cancer samples and prognostic roles for ADH1B, CDH11, and VGLL3. (A) The expression of chemotherapy resistant hub gene ADH1B in platinum-resistant and sensitive samples of GSE51373. (B) Survival plot of the significant hub gene ADH1B by Kaplan–Meier test. The data was extracted from the TCGA website. The Kaplan–Meier test P-value <.05. (C) The expression of chemotherapy resistant hub gene CDH11 in platinum-resistant and sensitive samples of GSE51373. (D) Survival plot of the significant hub gene CDH11 by Kaplan–Meier test. The data was extracted from the TCGA website. The Kaplan–Meier test P-value <.05. (E) The expression of chemotherapy resistant hub gene VGLL3 in platinum-resistant and sensitive samples of GSE51373. (F) Survival plot of the significant hub gene VGLL3 by Kaplan–Meier test. The data was extracted from the TCGA website. The Kaplan–Meier test P-value <.05. ADH1B = alcohol dehydrogenase 1B, CDH11= cadherin 11, VGLL3 = vestigial like family member 3.

Discussion

Approximately 80% of women with advanced OV, especially HGSOC will have tumor progression, or more commonly a recurrence, which is largely due to the emergence of chemotherapy resistance.[ Despite more than 10 years of array-based genomic and gene expressional profiling, no distinguishable signature of drug response has been identified for OVs.[ So, it is urgent and obligatory to understand the underlying mechanisms and discover the novel molecular targets of drug resistance for HGSOC relief. At the genomic level, researchers have found that gene breakage, which commonly inactivated the tumor suppressors RB1 (transcriptional corepressor 1), neurofibromin 1, RAD51B (RAD51 paralog B), and phosphatase and tensin homolog would contribute to acquired chemotherapy resistance.[ Also, a couple of molecular events, including “multiple independent reversions of germline BRCA1 (DNA repair associated) or BRCA2 (DNA repair associated) mutations in individual patients,” “loss of BRCA1 promoter methylation,” and “recurrent promoter fusion associated with overexpression of the drug efflux pump MDR1 (ATP binding cassette subfamily B member 1)” showed close relationship with acquired resistance.[ However, at the gene expression level, it remains for us to explore the marker genes, which might serve as therapeutic targets to surpass drug resistance. Therefore, in this study, we attempted to identify both chemotherapy resistant and sensitive genes of OV using WGCNA. By constructing a WGCNA co-expression network, we calculated the module-trait relationship, functionally annotated the trait related modules, and determined hub genes related with chemotherapy resistance and sensitivity. Finally, ADH1B, CDH11, VGLL3 for black module (Fig. 5A–F) and KIT for salmon module (Supplementary Fig. 2A and B) were found to be over-expressed in chemotherapy resistant samples. Besides, these 4 upregulated genes in platinum-resistance patients together with AQP5 (Supplementary Fig. 2C and D), a hub gene representative for platinum-sensitive responses, were verified to be prognostic biomarkers for OV patients by survival analysis in TCGA (The Cancer Genome Atlas) database.

For functional annotation analysis, “positive regulation of PI3K activity” was top enriched BPs in the salmon module, which indicates that the signaling pathway play important roles in platinum resistance. PI3Ks mediated pathway is a complex network, which coordinates signaling from a series of tyrosine kinase receptors,[ including epidermal growth factor receptor, vascular endothelial growth factor receptor, and fibroblast growth factor receptor to activate cellular growth and survival. PI3K proteins are composed of a catalytic p110 subunit (PIK3CA) and a regulatory subunit p85 (PIK3R) that mediates receptor binding.[ The direct binding to phosphotyrosine of such receptors leads to the allosteric activation of the catalytic subunit of PI3K and subsequently yield of the second messenger PIP3 from phosphatidylinositol 4,5-bisphosphate (PIP2).[ As a result, PI3K complex is plasma membrane localized and activates pyruvate dehydrogenase kinase 1 and AKT (serine/threonine kinase 1) proteins, which ultimately control crucial cellular survival and metabolic processes.[ The PI3K pathway has been identified as one of the most frequently altered pathway in OV, and the activation of PI3K pathway confer chemotherapy resistance via an extended cell survival.[ PIK3CA amplification was reported to associate with platinum-based chemotherapy on a large set of patients with primary OVs, and the presence of the driving amplification may predict resistance of chemotherapy for OV patients.[ Besides, in vitro analysis proved that overexpression of AKT converted platinum-sensitive OV cell line into platinum resistant cells, and downregulation of AKT reverse platinum resistance.[ Therefore, blockade of the PI3K pathway makes tumor cell vulnerable to platinum, and combinatorial chemotherapy in multiple human xenograft models.[

For OV and other complicated traits, like diseases and immunodeficiency in animals and abiotic stress response in plants, a series of physical and biochemical mechanisms were recruited by organisms to respond to the damages caused by genes’ mutations, such as signal transduction, tissue/organ tolerance, macromolecular compound synthesis, membrane structure conversion, and biochemistry homeostasis adjustment.[ Cancer genes could be divided into 2 categories based on their functional characters.[ The “caretakers,” which participate in fundamental functions that support genome stability and mutations in “caretakers” promote tumor via increasing the chances that mutations will hit some genes.[ The other is “gatekeepers,” which involve in cellular signaling and regulating processes[ and mutations promote tumor progression directly by changing cell differentiation. However, both “gatekeepers” and “caretakers” genes can possibly influence the convergence and/or divergence of complex traits.[ In the regulating network of drug resistance for HGSOC, the hub genes were speculated as “gatekeepers” and are related to influencing cooperation among cells (oncogenes) or to prevent the expansion of cheater cells (tumor suppressor genes).[

For hub gene analysis, which co-expressed with marker genes of OV, several genes were revealed to be of potential importance in platinum resistance and sensitivity listed in Table 2. For the black module of platinum resistance, 3 of the top-ranked hub genes including alcohol dehydrogenase 1B (class I), beta polypeptide (ADH1B) (Fig. 5A and B), CDH11 (Fig. 5C and D), and VGLL3 (Fig. 5E and F) were revealed to be dysregulated in platinum-resistant patients with OV. ADH1B encodes an essential enzyme in alcohol metabolism, which leads to variations in the amount of production of acetaldehyde between individuals.[ Consistent with our results of upregulation of ADH1B in chemoresistant patients, the mRNA expression of ADH1B was significantly increased when patients were treated with chemotherapies,[ which implied that ADH1B might closely related with drug resistance in OV. And it is also reported that high expression of ADH1B was associated with significantly higher risk of residual disease in HGSOC and correlated with promoted metastasis and proliferation of HGSOC.[CDH11 encodes CDH11, a mesenchymal cadherin, found to be upregulated in metastatic lesions and associated with advanced stage and nodal involvement of OV.[ Although there is no evidence for CDH11 participation in platinum-resistance for OV, it has been reported to mediate paclitaxel resistance in lung, skin, and liver cancers.[VGLL3 encodes a putative transcription co-factor, which was reported to be under-expressed in malignant ovarian tumors compared with benign samples, and might serve as a tumor suppressor gene.[ However, in our study, the gene was found to be upregulated in platinum-resistant patients with OV. So it is still worthy of studying the roles for VGLL3 in chemoresistance and illustrating its therapeutic potential in OV.

Gamma-aminobutyric acid type A receptor beta2 subunit (GABRB2) and DENN domain containing 2A (DENND2A) were identified as hub genes co-expressed with OV marker genes KIT and FGF2 in the salmon module. KIT gene encodes the human homolog of the proto-oncogene c-kit. C-kit protein is a type 3 transmembrane receptor for MGF (mast cell growth factor), its mutations are associated with diverse types of tumors, and is also expressed in OV cells.[ Consistent with our results, a clinical study has proved that c-kit positive cases were less responsive to chemotherapy than c-kit negative cases of OV,[ and prognosis of patients exhibiting c-kit expression is usually poorer, which is also validated by using TCGA database.[FGF2 gene encoded protein is a member of the fibroblast growth factor (FGF) family, which always bind to heparin and has angiogenic activity.[ Although, FGF2 has not been proved to be functional in platinum resistance, it has been pointed out that FGF2/PI3K/AKT autocrine loop is required for tumor protection upon bevacizumab treatment for OV by endothelial cells, which makes FGF2 a potential therapeutic target for bevacizumab resistant OV.[ For the other 2 hub genes GABRB2 and DENND2A, there is relatively few literatures about their relevance with cancer, or even platinum resistance, but our result may provide insight for functional study in OV chemotherapy resistance.

For magenta module, top 5 enriched hub genes, include ALDH1A1, FOXL2 (forkhead box L2), AOX1, ADCY8, and MUM1L1 (PWWP domain containing 3B). ALDH1A1 and FOXL2 are OV marker genes. Cancer stem cell is one of critical factors contributing to chemo-resistance, which is quiescent and not eliminated by chemotherapy.[ ALDH activity is revealed to be increased in cells grown as sphere versus monolayer cultures under differentiating conditions and in OV cells after treatment with platinum. Knockdown and inhibition of ALDH1A1 by its inhibitor CM37 could both reduce OV proliferation as sphere, expression of stemness markers, and delayed tumor initiation capacity in vivo.[ Mutation of FOXL2 (402C->G) is present in 97% of adult granulosa cell tumors ovary,[ and its role in chemo-resistance remains to be explored. The functions of hub gene encoded proteins, AOX1, ADCY8, and MUM1L1 have not been well demonstrated yet, and their potential impact on chemo-resistance is still worthy of study.

For platinum sensitivity analysis, the brown module contains 5 important hub genes with 4 reported marker genes of OV (ACPP, ALPPL2, AQP5, and BST2) and 1 suspected functional gene (apolipoprotein F [APOF]) in platinum response sensitivity. ACPP (prostatic acid phosphatase) is a glycoprotein, which gets high expression in epithelial cells of ovarian carcinoma.[ ALPPL2 (alkaline phosphatase, germ cell) protein is a membrane bound glycosylated enzyme, present in testis, thymus and certain germ cell tumors. It has been reported that increased level of serum alkaline phosphatase might serve as an early diagnosis of OV.[ AQP5 belongs to a family of transmembrane channel, and higher expression showed better overall survival in OV patients,[ which is consistent with our data explored in TCGA database that AQP5 correlated with better overall survival for OV cases. BST2 (bone marrow stromal cell antigen 2) known as a tetherin, is associated with lipid rafts and involved in the growth and development of B-cells.[ But the expression of BST2 was downregulated in cisplatin resistant OV cell lines, which predict a chemosensitive role of BST2 for OV patients.[ APOF protein forms complex with lipoproteins and may be functional in cholesterol transportation.[ Although, little is known about the linkage between APOF and OV, it is enriched with 4 marker genes in the platinum sensitive module and deserve deep research for its conceivable roles in OV.

Although we identified modules and hub genes platinum-based chemotherapy resistance and treatment response in OV through bioinformatics methods, we did not conduct experimental test for any of these selected genes, which was a limitation of this study. Functional analysis is necessary to further study the functions of these genes in platinum-based chemotherapy resistant and sensitive HGSOC regulation.

Conclusions

In this study, we performed WGCNA by using gene expressional data of platinum treated OV patients to identify both chemotherapy resistant and sensitive genes for OV. By constructing a WGCNA co-expression network, the module-trait relationship, functionally annotated the trait related modules were calculated, and hub genes in relation with chemotherapy resistance and sensitivity were determined. 16 modules with different colors were generated, of which, 51 genes were grouped to salmon module, 1297 and 1115 genes were assigned to grey and turquoise module, blue and brow module included 475 and 447 genes. “positive regulation of PI3K activity” (GO:0043552) may be important signaling pathway in OV platinum resistance. ADH1B, CDH11, VGLL3, and KIT were revealed to be expressional related with platinum resistance, and could serve as prognostic markers for OV. And the potential roles for other genes in OV for platinum resistance or sensitivity, including FGF2, GABRB2, DENND2A, AOX1, ALDH1A1, ADCY8, FOXL2, MUM1L1, ACPP, ALPPL2, AQP5, APOF, and BST2, remain worthy of exploration.

Author contributions

Conceptualization: Luoyan Zhang, Xuejie Zhang, Zhen Zhang.

Data curation: Luoyan Zhang, Xuejie Zhang, Zhen Zhang.

Formal analysis: Luoyan Zhang, Xuejie Zhang, Shoujin Fan.

Funding acquisition: Luoyan Zhang, Zhen Zhang.

Writing – original draft: Luoyan Zhang.

Writing – review and editing: Zhen Zhang.