Identification of potential hub genes and therapeutic drugs in ovarian cancer via bioinformatics analysis

Abstract

This study aimed to identify hub genes and drugs for ovarian cancer (OC) using bioinformatics analysis. Gene expression data and clinical information were downloaded from public databases. Differentially expressed genes (DEGs) were identified, and ovarian disease-related genes were extracted from DISEASES database. Disease-associated DEGs were further screened. Functional enrichment, protein–protein interaction (PPI) analysis, drug prediction, and survival analysis were performed. Key gene-drug interactions were also predicted. Finally, hub genes expression in an OC cell line (SKOV3) was verified using real-time quantitative reverse transcription polymerase chain reaction (RT-qPCR). Totally 354 overlapping DEGs were identified. PPI network revealed genes such as CCL5 and CXCR4 with higher degrees, which were mainly involved in chemokine-related functions. CCND1, CCNE1, CXCL10, ERBB4, LPAR3, and SST were correlated with prognosis of OC patients. Among them, CXCL10 was considered as an independent prognostic factor. Moreover, 63 possible drugs were predicted and the sirolimus-DCN pair was screened via a literature search. RT-qPCR confirmed that the levels of PDGFRB, DCN, PDGFRA, LYN, and CENPE were consistent with the bioinformatics analysis results. Our findings provide new insights into the underlying mechanisms of OC pathogenesis. The identified hub genes and drugs may improve individualized OC diagnosis and therapy.

KEYWORDS:

• Ovarian cancer
• differentially expressed genes
• functional enrichment analysis
• survival analysis
• drug prediction

Introduction

Ovarian cancer (OC) has the highest mortality rate among female genital tract cancers (Torre et al. 2018). Although the therapeutic effects of surgery and chemotherapy have improved owing to advancements in medical technology, clinical therapeutic results are still not ideal (Brackmann et al. 2018). Early diagnosis of OC is anticipated to substantially improve survival rates compared to cases of late diagnosis (Chatterjee et al. 2017). Therefore, it is necessary to identify the molecular mechanisms of OC development and metastasis.

Molecular genetic studies have shown that mutations in tumor suppressor genes, such as BRCA1, are responsible for the majority of hereditary OCs (Holschneider and Berek 2000). A recent study showed that certain differentially expressed genes (DEGs), such as TRAF4, contribute to metastasis and poor prognosis (Javdan et al. 2018). In a bioinformatics study, DEGs, including CCND1 and CCNB2, were found to be involved in several pathways related to OC development (Fu and Wang 2013), such as the cell cycle, drug metabolism, and drug resistance (Richardson and Kaye 2005). A previous gene expression network analysis revealed that securitine inhibits the expression of its target gene, STAT3, which contributes to OC progression (Chen et al. 2017). Despite such studies, the complex mechanisms associated with OC progression are not yet fully understood.

The Cancer Genome Atlas (TCGA) contains data describing the molecular characteristics of over 20,000 primary cancers and is widely used in cancer research (Wang et al. 2016). Most previous studies have analyzed data from a single database, which does not reflect the complex regulatory mechanisms of OC. Therefore, multiple databases must be integrated to identify potential biomarkers and therapeutic targets for OC. In this study, gene expression data of OC tumors and normal samples were downloaded from the University of California Santa Cruz Xena (UCSC-Xena) platform. Furthermore, ovarian disease-related genes were screened using the DISEASES database. Next, these genes were intersected with DEGs, and overlapping DEGs were selected for further analysis. Subsequently, a series of bioinformatics analyses, including functional enrichment analysis, protein–protein interaction (PPI) analysis, survival analysis, and drug prediction, were performed to identify key genes and drugs. The expression levels of hub genes in human OC cell lines were determined through quantitative reverse transcription real-time polymerase chain reaction (qRT-PCR). The objective of this study was to investigate the detailed molecular mechanisms underlying OC progression and explore potential drugs for OC therapy.

Material and methods

Microarray data

RNA sequencing data were downloaded from the ‘UCSC Toil RNAseq Recompute’ Data Hub (https://toil.xenahubs.net) through the UCSC-Xena platform (https://xena.ucsc.edu/) (Goldman et al. 2020). RNA-seq expression data of 19,131 patients were downloaded from https://xenabrowser.net/datapages/?data set=TcgaTargetGtex_rsem_gene_tpm&host=https%3A%2F%2Ftoil.xxenahub.net&removeHub=https%3A%2F%2Fxena.treehouse.gi.ucsc.edu%3A443. Their corresponding clinical information was collected from https://xenabrowser.net/datapages/?dataset=TcgaTar-getGTEX_phenotype.txt&host=https%3A%2F%2Ft-oil.xexenahu.net&removeHub=https%3A%2F%2Fxe-na.treehouse.gi.ucsc.edu%3A443. Subsequently, the TCGA-OV samples with detailed_category set to ‘Ovarian Serous Cystadenocarcinoma’ were extracted and 427 cases with tumors were obtained. The GTEx samples with detailed_category set to ‘Ovary’ were extracted and 88 normal samples were obtained. The clinicopathological characteristics of these patients are presented in Supplementary Table 1.

Based on the gene annotation information (hg38) from the GENCODE database (https://www.gencode-genes.org/) (Harrow et al. 2012), 60,498 Ensembl IDs were converted to gene symbols. Genes with low expression levels were filtered out, and genes with expression levels > 0 in over 1/3 of samples were retained. Finally, a total of 16,746 genes were identified.

Identification of DEGs between OC and normal samples

Linear regression and empirical Bayesian methods (Quadrianto and Buntine 2016; Koop 2017) in the R limma package (version 3.10.3, http://www.biocon-ductor.org/packages/2.9/bioc/html/limma.html) (Smyth 2015) were used to analyze the DEGs between the tumor and normal groups. P-values were corrected using the Benjamini-Hochberg method. An adjusted P-value < 0.05 and |log2 fold change (FC)| > 2 were selected as the thresholds for DEG screening. Volcano plots were constructed using ggplot2 (version 3.2.1).

Screening of disease-related genes

Genes related to ovarian diseases including Brenner tumors, OC, and ovarian carcinosarcoma were screened from the DISEASES database (https://diseases.jensenlab.org/Search) (Pletscher-Frankild et al. 2015) using ‘ovarian’ as the keyword. Disease-related genes were integrated with the DEGs, and overlapping DEGs were selected for further analysis.

Functional enrichment analysis

Gene Ontology-Biological Process (GO-BP) (Ashburner et al. 2000) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway (Kanehisa and Goto 2000) enrichment analyses of the overlapping DEGs were performed using the Metascape online tool (http://metascape.org) (Zhou et al. 2019). Only terms with P-value < 0.01, a minimum count of 3, and an enrichment factor > 1.5 were considered significant. Clustering analysis was conducted according to the similarity of genes enriched in each term (similarity > 0.3). The most statistically significant terms (those with minimum P-values) were selected and used to construct the cluster. The top 20 clusters based on P-values were visualized using a histogram.

To further explore the relationships among terms, an interaction network of the top 20 clusters was constructed using Cytoscape (version 3.4.0, http://chianti.ucsd.edu/cytoscape-3.4.0/) (Shannon et al. 2003). Terms were included in the network only if they had the minimum P-value in their cluster and similarity > 0.3.

PPI network and module construction

PPI pairs of identified DEGs were identified using three databases in Metascape software (Zhou et al. 2019), namely BioGrid (Chatr-Aryamontri et al. 2017), InWeb_IM (Li et al. 2017), and OmniPath (Türei et al. 2016), with default parameters (Min Network Size: 3; Max Network Size: 500). After obtaining pairs, PPI networks were constructed using Cytoscape software (version 3.4.0) (Shannon et al. 2003).

Significantly enriched modules from the PPI network were identified using Molecular Complex Detection (MCODE; version 1.5.1) (Bader and Hogue 2003). GO-BP and KEGG analyses of genes in the modules were performed using Metascape (Zhou et al. 2019).

Drug-gene interaction and drug-related pathway analysis

Drug-gene pairs were predicted using the Drug-Gene Interaction database (DGIdb version 2.0; http://www.dgidb.org/search_interactions) (Cotto et al. 2018), and the drug-gene network was established using Cytoscape software (Shannon et al. 2003).

Additionally, the online database DrugPath (http://www.cuilab.cn/drugpath) (Shah et al. 2014) was used to search for KEGG pathways affected by the drugs. The drug-pathway network was constructed using Cytoscape software (Shannon et al. 2003). The three-dimensional (3D) structures and PubChem IDs of the predicted drugs were downloaded from PubChem (https://pubchem.ncbi.nlm.nih.gov/).

Survival analysis of potential genes

Nodes in both the PPI modules and drug-gene network were selected as potential genes. The expression levels of the genes were analyzed using Gene Expression Profiling Interactive Analysis 2 (GEPIA2) (http://gepia2.cancer-pku.cn/#index) (Tang et al. 2019). The prognostic value of each gene was assessed using the Survival Analysis package with the quartile grouping method, and the samples were divided into different groups based on gene expression. Briefly, the top and bottom 25% were defined as high– and low–expression groups, respectively. Differences in overall survival (OS) and disease-free survival (DFS) between the two groups were compared using Kaplan – Meier survival curve analysis (Bland and Altman 1998), and statistical significance was set at P-value < 0.05.

Univariate and multivariate Cox regression analyses were performed to assess whether the OS-related genes were independent prognostic factors. The following clinicopathological variables were considered: age, clinical stage, neoplasm histological grade, postoperative_rx_tx, and response to primary therapy. Statistical significance was set at P-value < 0.05.

Literature research and protein expression verification

Among the candidate genes, key genes were further screened from the National Center for Biotechnology Information (NCBI) PubMed database (https://pub-med.ncbi.nlm.nih.gov/) using Biopython (Cock et al. 2009) with the keyword ‘gene symbol + ovarian cancer’. Genes reported in more than three and less than 10 articles were considered hub genes. The Human Protein Atlas (HPA; https://www.proteinatlas.org/) (Uhlén et al. 2015) was used to verify the protein expression levels of key genes.

Investigation of key drug-gene interactions

The key genes and their corresponding target drugs were extracted. Briefly, the drugs were identified from the NCBI PubMed database using Biopython (Cock et al. 2009) with ‘drug + ovarian cancer’ as keywords. To narrow down the results, drugs reported in more than three and less than 10 studies were regarded as key drugs. Finally, the obtained drugs and their target genes were considered key drug-gene interactions.

Cell culture

The human normal ovarian cell line IOSE80 and human OC cell line SKOV3 were purchased from the Cell Bank of the Chinese Academy of Sciences (Shanghai, China). The cells were cultured in 1640 medium (Cuisi Bio, Shanghai, China) containing 10% fetal bovine serum (Gibco, Shanghai, China) in a 37°C incubator with 5% CO₂.

qRT-PCR

Total RNA was extracted from SKOV3 and IOSE80 cells using TRIzol reagent (TaKaRa Bio, Beijing, China) and RNA purification and integrity were assessed using an Agilent Bioanalyzer (Santa Clara, CA, USA). RNAs were then reverse-transcribed into complementary DNA (cDNA) using a Reverse Transcription kit (TOYOBO, Shanghai, China), and qRT-PCR was performed using SYBR Green PCR mix (Thermo Fisher Scientific, Waltham, USA) on an ABI 7900HT FAST Real-Time PCR System (Applied Biosystems, Foster City, CA, USA). (Glyceraldehyde-3-phosphate dehydrogenase) GADPH was used as an internal control and the relative mRNA expression levels were calculated using the 2^−ΔΔCT method. The following primers were used: platelet-derived growth factor receptor, beta polypeptide (PDGFRB, forward (F): 5’-AGCACCTTCGTTCTGACCTG-3’; reverse (R): 5’-TATTCTCCCGTGTCTAGCCCA-3’), decorin (DCN, F: 5’-ATGAAGGCCACTATCATCC-TCC-3’; R: 5’-GTCGCGGTCATCAGGAACTT-3’), platelet-derived growth factor receptor, alpha polypeptide (PDGFRA, F: 5’-TGGCAGTACCCCATGTCT-GAA-3’; R: 5’-CCAAGACCGTCACAAAAAGGC-3’), aurora kinase B (AURKB, F: 5’-CAGTGGGACACCC-GACATC-3’; R: 5’-GTACACGTTTCCAAACTT-GCC-3’), v-yes-1 Yamaguchi sarcoma viral related oncogene homolog (LYN, F: 5’-GCTTTTGGCACCA-GGAAATAGC-3’; R: 5’-GCTTTTGGCACCAGGAA-ATAGC-3’), centromere protein E (CENPE, F: 5’-GATTCTGCCATACAAGGCTACAA-3’; R: 5’-GAT-TCTGCCATACAAGGCTACAA-3’), and GAPDH (F: 5’-GGAGCGAGATCCCTCCAAAAT-3’; R: 5’-GGC-TGTTGTCATACTTCTCATGG-3’).

Statistical analysis

The qRT-PCR data were analyzed using GraphPad Prism v.8.0.2 (GraphPad Software, San Diego, CA) and reported as mean ± standard deviation. A two-tailed t-test was performed to compare the values between the normal and tumor groups. Statistical significance was set at P-value < 0.05.

Results

Identification of DEGs

A total of 3987 DEGs (2031 upregulated and 1956 downregulated) were identified between the tumor and normal groups (Figure 1(A)). Furthermore, 1442 ovarian disease-related genes were obtained from the DISEASES database. Subsequently, 354 overlapping DEGs were selected for further analyses (Figure 1(B)).

Figure 1. Identification of ovarian cancer (OC)-related differentially expressed genes (DEGs). (A) Volcano plots of DEGs between tumor and normal samples. Blue represents downregulated DEGs, red represents upregulated DEGs. (B) Venn diagram of the DEGs and ovarian disease-related genes.

Functional enrichment analysis of DEGs

DEGs were significantly enriched in 319 GO-BP terms and 20 KEGG pathways. The top 20 enriched categories are shown in Figure 2(A). Briefly, genes were mainly enriched in GO-BP terms such as negative regulation of cell cycle, regulation of protein kinase activity, and negative regulation of cell proliferation. DEGs were mainly enriched in pathways such as cell cycle and pathways in cancer. Furthermore, a network based on the top 20 enriched GO terms and KEGG pathways was constructed to observe the interactions among terms (Figure 2(B)). The terms were divided into different clusters, indicating potential close connections between different biological functions during OC pathogenesis.

Figure 2. Functional enrichment analysis of DEGs. (A) Top 20 enriched Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. (B) Network of enriched GO terms and KEGG pathways. Different node colors indicate different clusters, and the lines indicate similarities among genes annotated with each term.

PPI network and module analysis

A PPI network composed of 1134 edges and 292 nodes was constructed (Figure 3). Ten modules were selected from the PPI network based on MCODE analysis (Figure S1), and functional enrichment analysis was performed for the genes in these modules. Several genes, including chemokine (C-X-C motif) ligand 10 (CXCL10), chemokine (C–C motif) ligand 5 (CCL5) and chemokine (C-X-C motif) receptor 4 (CXCR4), were mainly enriched in GO-BP terms such as chemokine-mediated signaling, cellular response to chemokines, and response to chemokines. The genes in the modules were mainly enriched in KEGG pathways, such as cell cycle, pathways in cancer, and proteoglycans in cancer. The top three enriched functions and pathways are listed in Table 1.

Figure 3. Protein-protein interaction (PPI) network. Different colors indicate the different modules identified using Molecular Complex Detection (MCODE). Triangles and inverted triangles indicate the upregulated and downregulated genes, respectively. Node size indicates the degree of connectivity, and gray lines represent interactions between proteins.

Table 1. The GO-BP function and KEGG pathway of DEGs in 10 modules extracted from PPI network.

Modules	GO	Description	Count	Log10(P)
Module 1	GO:0070098	Chemokine-mediated signaling pathway	13	−28.7
Module 1	GO:1990869	Cellular response to chemokine	13	−28.1
Module 1	GO:1990868	Response to chemokine	13	−28.1
Module 2	GO:0051301	Cell division	15	−24.2
Module 2	GO:0000070	Mitotic sister chromatid segregation	11	−21.2
Module 2	GO:0140014	Mitotic nuclear division	12	−20.9
Module 3	GO:2001022	Positive regulation of response to DNA damage stimulus	4	−7.3
Module 3	GO:2001020	Regulation of response to DNA damage stimulus	4	−5.9
Module 3	GO:0051054	Positive regulation of DNA metabolic process	4	−5.8
Module 4	GO:0051321	Meiotic cell cycle	5	−8.2
Module 4	hsa04110	Cell cycle	4	−7.3
Module 4	hsa05200	Pathways in cancer	5	−7.2
Module 5	GO:0043966	Histone H3 acetylation	3	−6
Module 5	GO:0007098	Centrosome cycle	3	−5.1
Module 5	GO:0031023	Microtubule organizing center organization	3	−5
Module 6	GO:0038128	ERBB2 signaling pathway	3	−7.1
Module 6	GO:0097191	Extrinsic apoptotic signaling pathway	4	−6.6
Module 6	GO:0038034	Signal transduction in absence of ligand	3	−6.1
Module 7	GO:0090068	Positive regulation of cell cycle process	5	−8.8
Module 7	GO:0045787	Positive regulation of cell cycle	5	−8.2
Module 7	hsa04110	Cell cycle	4	−8
Module 8	GO:0018108	Peptidyl-tyrosine phosphorylation	5	−8.3
Module 8	GO:0018212	Peptidyl-tyrosine modification	5	−8.3
Module 8	GO:0030335	Positive regulation of cell migration	5	−7.5
Module 9	GO:0010389	Regulation of G2/M transition of mitotic cell cycle	3	−6.3
Module 9	GO:1902749	Regulation of cell cycle G2/M phase transition	3	−6.2
Module 9	GO:0000086	G2/M transition of mitotic cell cycle	3	−6
Module 10	hsa05205	Proteoglycans in cancer	3	−6.2
Module 10	GO:0030198	Extracellular matrix organization	3	−5.4
Module 10	GO:0043062	Extracellular structure organization	3	−5.3

Notes: GO-BP, Gene Ontology-biological process; KEGG, Kyoto Encyclopedia of Genes and Genomes; DEGs, differentially expressed genes; PPI, protein-protein interaction network.

Drug target prediction and drug pathway analysis

Based on the 354 DEGs, the online DGIdb database was used to predict the potential drugs. A total of 202 drug-gene interactions were identified, involving 70 DEGs and 63 drugs. The detailed results are shown in Figure 4(A). To further understand the pathways affected by the identified drugs, the online database, DrugPath, was used to search for KEGG pathways. As shown in Figure 4(B), the drug pathway network included 69 KEGG pathways (including the chemokine signaling pathway) and 16 drugs (including sirolimus). The 3D structures and corresponding PubChem IDs of the 16 drugs are shown in Figure S2.

Figure 4. Gene-drug interaction network (A) and drug-pathway network (B). (A) Red triangles indicate upregulated genes; green inverted triangles indicate downregulated genes; gray rectangles indicate predicted drugs. (B) Yellow squares indicate KEGG pathways and white squares indicate predicted drugs.

Survival analysis of potential genes

In total, 27 potential genes were identified after overlaying the genes in the PPI modules and the drug-target network. GEPIA2 was used to analyze the expression levels of 27 genes, which were consistent with the RNA-seq results (Figure S3). Furthermore, survival analysis showed that six genes (cyclin D1 (CCND1), cyclin E1 (CCNE1), CXCL10, v-erb-b2 avian erythroblastic leukemia viral oncogene homolog 4 (ERBB4), lysophosphatidic acid receptor 3 (LPAR3), and somatostatin (SST)) were significantly correlated with OS, and two genes (ERBB4 and LPAR3) were significantly associated with DFS (Figure 5).

Figure 5. Survival analysis of key genes in OC using Kaplan – Meier analysis. (A) CCND1, (B) CCNE1, (C) CXCL10, (D) ERBB4, (E) LPAR3, (F) SST, (G) ERBB4, and (H) LPAR3. Patients were split into high and low expression groups based on the median expression of genes.

Moreover, OS-related genes and clinical variables, including age, clinical stage, neoplasm histologic grade, postoperative_rx_tx, and primary therapy response, were analyzed using univariate and multivariate Cox regression analyses. In the univariate analysis, age, clinical stage, postoperative_rx_tx, primary therapy response, and four genes (SST, CXCL10, ERBB4, and LPAR3) were associated with clinical outcomes (P < 0.05; Figure S4). In the multivariate analysis, only age, primary therapy response, and CXCL10 levels were considered independent prognostic factors (P < 0.05).

Literature review and protein level verification

Following a literature review, 11 key genes meeting the screening criteria were identified. The protein levels of the genes were verified using the HPA database. The results revealed that the protein expression levels of PDGFRB, DCN, and PDGFRA were downregulated in OC tissues, whereas those of AURKB, LYN, and CENPE were upregulated (Figure 6).

Figure 6. Protein expression of key genes verified using the Human Protein Atlas (HPA) database. (A) AURKB, (B) CENPE, (C) DCN, (D) LYN, (E) PDGFRA, and (F) PDGFRB.

Identification of key drug-gene interactions

Among the six gene-targeted drugs, only sirolimus met the literature search criteria. By combining the results of the identified key genes and drugs, sirolimus-DCN was identified as a key drug-target pair.

Experimental verification of six key mRNAs

To verify the results of the bioinformatics analysis, qRT-PCR was performed on the six key nodes. Excluding AURKB, the expression of the remaining five genes differed significantly between normal and cancer cells (P < 0.05; Figure 7). Briefly, PDGFRB, DCN, and PDGFRA mRNA were downregulated, whereas LYN and CENPE were upregulated in the OC cell line SKOV3, consistent with our bioinformatics analysis results.

Figure 7. Validation of key genes using quantitative reverse transcription polymerase chain reaction (qRT-PCR). Data are presented as mean ± standard deviation. Statistical significance was determined using two-tailed t-test; * P < 0.05; *** P < 0.001.

Discussion

Although the expression of certain genes may contribute to OC progression, the details of the underlying molecular mechanisms remain unclear. A total of 354 DEGs were identified between tumor and normal samples. PPI analysis showed that genes such as CCL5 and CXCR4 with higher degrees, which were mainly enriched in chemokine-associated functions. Moreover, 63 drugs were predicted using DGIdb. Survival analysis showed that CCND1, CCNE1, CXCL10, ERBB4, LPAR3, and SST were correlated with OS, and ERBB4 and LPAR3 were associated with DFS. Among them, CXCL10 expression was considered an independent prognostic factor. Moreover, the protein and mRNA levels of PDGFRB, DCN, LYN, PDGFRA, and CENPE were consistent with the results of the bioinformatics analysis. Finally, sirolimus-DCN was selected as the key drug-gene pair.

CCL5 is a chemokine-encoding gene that promotes cancer cell migration under several pathological conditions (Liou et al. 2012), and has been widely used as a clinical prognostic factor to predict cancer progression (Yaalhahoshen et al. 2006). Notably, CCL5 expression is closely associated with OC invasion and migration (Long et al. 2012), and OC development has been reported to be influenced by the CCL5-triggered epithelial–mesenchymal transition (Long et al. 2015). CXCR4 is a prognostic marker in various types of cancers, including OC (Jiang et al. 2006), and the overexpression of CXCR4, which plays an important role in cancer-associated signaling, promotes metastasis (Obermajer et al. 2011). In an animal model, CXCR4 inhibition prolonged survival in mice with OC (Righi et al. 2011). Although CCL5 and CXCR4 have been detected in a large number of solid tumors (Weitzenfeld and Benbaruch 2014), the detailed molecular mechanisms whereby these genes influence OC progression remain unclear.

The biological functions of DEGs in OC are realized via certain mechanisms and signaling pathways (Zhou et al. 2016). A previous study showed that chemokine (C–C motif) ligand 2 (CCL2)/chemokine (C–C motif) receptor 2 (CCR2) influences the survival and motility of breast cancer cells via the chemokine signaling pathway (Fang et al. 2012). Studies on the chemokine signaling system in the OC microenvironment have contributed to the development of anti-angiogenic therapies (Agarwal et al. 2010). In the present study, we found that CCL5 and CXCR4 were enriched in chemokine-associated biological functions. Therefore, we speculated that CCL5 and CXCR4 may participate in OC progression via chemokine-associated functions.

Drug resistance and effectiveness affect the success of chemotherapy (Norouzi-Barough et al. 2018). Although various drugs, such as cisplatin, are widely used in the clinical treatment of OC, drug resistance prevents their widespread clinical use (Samuel et al. 2018). Thus, we screened for potential drugs for use in OC treatment. Sirolimus (also called rapamycin) is a class I mammalian target of rapamycin (mTOR) inhibitor that exhibits antitumor activity and has been approved for breast cancer treatment (Pintoleite et al. 2009; Acevedogadea et al. 2015). The efficacy and safety of sirolimus for the treatment of OC have been demonstrated in previous studies (Jiang and Feng 2004; Stippel et al. 2005). Additionally, rapamycin can block growth differentiation factor 15 (GDF-15)-mediated invasion of human ovarian tumor tissues by inhibiting mTOR-matrix metalloproteinases (MMP) signaling (Griner et al. 2013). In this study, we found that DCN was a target for sirolimus. DCN encodes a member of the leucine-rich proteoglycan family and is involved in the regulation of various important cellular functions, such as matrix assembly, cell proliferation, and migration (Halari et al. 2021). Binding of this protein to several cell surface receptors can mediate its tumor-suppressive effects, including the inhibition of angiogenesis and tumorigenesis. Recent studies have shown that DCN is downregulated in OC tissues and can inhibit growth of OC cells (Nash et al. 2002). DCN is also associated with the upregulation of TNFα and downregulation of VEGF, suggesting that it may play a key role in OC angiogenesis and tumor progression (Grisaru et al. 2007). These studies have revealed that DCN may be a tumor suppressor and a potential biomarker of OC. However, the relationship between DCN and sirolimus in OC has not been reported. Therefore, whether sirolimus exerts a therapeutic effect on OC by regulating the expression of DCN remains to be confirmed.

CXCL10 has been widely used as a prognostic indicator of the progression of various diseases (Jiang et al. 2010; Bronger et al. 2016). A previous study showed that CXCL10 participates in the development of ovarian tumor cells (Rainczuk et al. 2014). In an animal model, CXCL10 altered the OC tumor immune microenvironment (Au et al. 2017) and could be used to predict survival in patients with advanced serous OC (Bronger et al. 2016), which is consistent with our findings. Thus, we speculate that CXCL10 may be a novel prognostic marker for OC.

The drug-gene interactions explored in the present study may provide insights into novel targets and therapeutic drugs associated with the development and metastasis of OC; however, this study has several limitations. Although there are different types of OC, the TCGA database does not provide a detailed OC classification; therefore, OC type was not considered in this analysis. We plan to include patients with different types of OC in future studies. Moreover, drug-gene interactions, such as sirolimus-DCN, were obtained via bioinformatics analyses, and the results need to be verified via additional experiments.

Conclusion

In summary, key genes and potential drugs associated with OC were identified using bioinformatics analyses. These key genes, such as DCN, may serve as reliable molecular biomarkers for OC diagnosis, and predicted drugs, such as sirolimus, may serve as potential drugs for OC treatment. Further research is required to explore the biological functions of these genes and the mechanisms underlying OC pathogenesis.

Author contributions

Conceptualization and design of the research: JS and XZ; acquisition of data: ZS, JC, and DW; analysis and interpretation of data: ZS, JC, and DW; statistical analysis: ZS, JC, and DW; funding acquisition: JS; drafting of the manuscript: XZ; revision of the manuscript for important intellectual content: JS. All authors have read and approved the final version of the manuscript.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Data availability statement

The data were derived from the following resources available in the public domain: the UCSC Xena database (https://xena.ucsc.edu/), GENCODE database (https://www.gencodegenes.org/), and DISEASES (https://diseases.jensenlab.org/Search). The qRT-PCR data in this study are available in 4TU. ResearchData under DOI: https://doi.org/10.4121/22267030 at https://doi.org/10.4121/22267030.

Additional information

Funding

This work was supported by the Doctor Initiation Research Fund of Liaoning Province (grant number: 2020-BS-034).