Network Models of Prostate Cancer Immune Microenvironments Identify ROMO1 as Heterogeneity and Prognostic Marker

doi:10.21203/rs.3.rs-786726/v1

Download PDF

Research

Network Models of Prostate Cancer Immune Microenvironments Identify ROMO1 as Heterogeneity and Prognostic Marker

https://doi.org/10.21203/rs.3.rs-786726/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 07 Jan, 2022

Read the published version in Scientific Reports →

Version 1

posted

You are reading this latest preprint version

Background: Prostate cancer (PCa) is the fifth leading cause of death from cancer in men worldwide. Its treatment remains challenging due to the heterogeneity of the tumor, mainly because of the lack of effective and targeted prognostic markers at the system biology level.

Methods: First, samples of PCa were retrieved from UCSC Xena data, and valid samples were obtained by consistent clustering and principal component analysis; next, key genes were analyzed for prognosis of PCa using WGCNA, MEGENA, and LASSO Cox regression model analysis, while key genes were screened based on disease-free survival (DFS) significance. Finally, TIMER data were selected to explore the relationship between genes and tumor immune infiltration, and GSCAlite was used to explore the small-molecule targeted drugs that act with them.

Results: In this study, we used tumor subtype analysis and an energetic co-expression network algorithm of WGCNA and MEGENA to identify a signal dominated by the ROMO1 tag to predict PCa prognosis. Cox regression analysis of this signal was an independent influence, and the prognostic value of this biomarker was validated in the training set, the validated data itself, and external data, respectively. This biomarker correlates with tumor immune infiltration and has a high degree of infiltration, poor prognosis, and strong correlation with CD8+ T cells. Gene function annotation and other analyses also implied a potential molecular mechanism for ROMO1.

Conclusions: In conclusion, our study suggests ROMO1 as a portal key prognostic gene for the diagnosis and prognosis of PCa. Moveover, we identify nutlin-3, fluorouracil, and other drugs as possible key drugs for diagnosis and prognosis, which provide new insights into the diagnosis and treatment of prostate cancer.

General Biochemistry

Molecular Biology

Prostate cancer

WGCNA

MEGENA

Intratumoral heterogeneity

tumor immune microenvironment

Prostate cancer (PCa), an epithelial malignancy of the prostate, is the second most common cancer in the world and one of the most common malignancies in men [1, 2]. At present, PCa has become a global public problem threatening human health and life [3], with more than 1,275,000 new diagnoses and 350,000 deaths each year. The factors associated with the development of PCa are mainly involved in age, race, family heritage, geography, and diet [4–6]. According to the Epidemiology and End Results report (https://seer.cancer.gov//), the five-year survival rate for PCa is 97.5% [8], and it will account for 13.1% of all new cancers in 2021, according to the Epidemiology and End Results report. Although there are many types of research on PCa, its etiology and pathogenesis remain not fully understood. The current treatments for the PCa include radical surgery [7], external beam radiotherapy [8], brachytherapy [9], local treatment of experimental PCa [10], endocrine therapy [11], and chemotherapy [12], etc. Serum prostate-specific antigen (PSA) is now commonly used as the preferred marker for screening PCa, with low specificity [13]. The immune system provides a huge role in the prevention and treatment of tumors [14]. Drugs developed to target the immune system, such as cancer vaccines, can effectively target specific antigens expressed on the surface of cancer cells, thereby protecting normal cells from drug damage[15]. However, due to the heterogeneity within tumors and because cancer vaccines target only a limited number of specific antigens, cancer cells that do not express these antigens or a group of cancer cells that have mutations that result in altered surface antigens can evade immune control and thus create new tumor populations that can resist treatment with vaccines encoding the same tumor-associated antigens[16, 17]. Therefore, the development of new effective oncology drugs against heterogeneous populations of tumours has become an urgent task. However, immunotherapy for prostate cancer has not been efficacious for patients in the past. Therefore, comprehensive bioinformatics screening for new diagnostic and prognostic indicators and the proposal of new effective cancer drugs are of great importance for the treatment and diagnosis of PCa.

Intratumoral heterogeneity is a vital feature of tumorigenesis [18]. As cancer cells divide and proliferate, somatic mutations accumulate, some of which give cancer cells a more significant adaptive advantage [19–22]. Meanwhile, alterations in epigenetics with the tumor microenvironment may lead to different molecular subtypes that have varying tissue types, invasive capabilities, and degrees of differentiation [23–26]. The presence of heterogeneity makes tumor treatment very challenging. While the discovery of tumor driver genes and associated drugs targeting these genes offers the possibility of treating these tumors [27, 28]. Precision medicine for PCa needs to consider different molecular subtypes so that tumor heterogeneity can be achieved and patients can be offered personalized treatment options.

With the rapid development of next-generation sequencing technologies, bioinformatics analysis allows us to understand tumor characterizes at a multi-omics level [29]. Research in RNA sequencing (RNA-seq), Co-expression network analysis, immunohistochemistry, immuno-infiltration, and targeted drug studies, etc, are now widely used. RNA sequencing can study the transcriptional profiles of cell populations or the average expression levels of tissues [30]. The Cancer Genome Atlas (TCGA) [31], Gene Expression Profiling Interactive Analysis (GEPIA)[32], The Human Protein Atlas (HPA)[33], etc are open databases for the collection and storage of patient genomic data. This study explores the PCa RNA-seq data to discover new biomarkers and provide new insights for the diagnosis and treatment of PCa.

Data acquisition

Gene expression profiles and clinical data for PCa were downloaded using the UCSC Xena (http://xena.ucsc.edu/) database [34], containing 436 PCa tissue and 114 normal tissue. Patients with unclear survival time and survival status characteristics were excluded from the cohort, and all data were collected on 10 April 2021.

Determination of cluster count and membership

Based on the expression profile of PCa, ConensusClusterPlus (v. 1.52.0) R language software package has carried out unsupervised clustering classification for 550 samples [35]. The project selection criteria were selected 1000 iterations and 80% resampling rate Pearson correlation. Reshape2 (v. 1.4.4) and factoextra (v. 1.0.7) two R language software packages were used to reduce the dimensions of clusters of different clusters to verify the results, and the first two PC components of the principal component analysis were selected for subsequent analysis and ggplot2 (v. 3.3.3) R language software package was used to draw the graph to display the results.

Identification of DEGs based on consensus clustering and PCA

The Limma (v. 3.46.0) R language package was used to screen differentially expressed genes (DEGs) among noncancerous samples with different prostate carcinoma subtypes [36]. The DEGs are identified by significant cutoffs. Genes with |log2FC| (fold change) ≥ 1 and P-value ≤ 0.01 is considered as DEGs. The R language software packages of ggplot2 (v. 3.3.3) and RColorBrewer (v. 1.1.2) were used to draw the graph to display the results. The Venn diagram of the overlapping DEGs was outputted by the Draw Venn diagram website (http://bioinformatics.psb.ugent.be/webtools/Venn/), and the Upset plot was drawn by the UpSetR R language package.

Functional enrichment analysis of noncancerous samples and different tumor subtypes of DEGsTo elucidate the biological functions of the DEGs, we performed further functional annotation and pathway enrichment analysis using the Metascape (https://metascape.org) online tools [37]. P-value < 0.05 was considered as the cutoff value.

Gene co-expression network analysis to obtain key gene modules associated with tumorigenesis and malignant progression

Gene co-expression network analysis is widely used in identifying functional gene modules associated with human disease and has proven to be highly effective [38, 39]. The complexity of tumors is further expressed through intratumoral heterogeneity and the tumor microenvironment, which encompasses intratumoral cell differences and interactions and the expression of several biomarkers between normal and different tumor cells[40, 41]. In this study, we used Weighted Gene Co-expression Network Analysis (WGCNA) and Multiscale Embedded Gene Co-expression Network Analysis (MEGENA) to identify key genes associated with prostate carcinogenesis as well as tumor malignant progression WGCNA (v. 1.70-3) R language package was used to perform a WGCNA of all DEGs with clinical traits, mining modules, and key genes associated with clinical characteristics [42]:

The gene correlation matrix was transformed into a scale-free network.
Unsupervised hierarchical clustering was performed using a dynamic tree chopping algorithm, and the clustered tree branches were defined as modules.
The correlation between each module and clinical traits was calculated, and the modules with higher correlations were selected.

Subsequently, the hub genes in the selected modules were ranked by intra-module connectivity and correlation with module trait genes, and finally, candidate genes were identified.

MEGENA is an innovative co-expression network analysis method that offers unique advantages over WGCNA for efficiently constructing large-scale co-expression planar filter networks and preserving gene interactions [43]. The R software package MEGENA (v. 1.3.7) is used to perform MEGENA, which consists of the following steps: (1) constructing a planar filtered network (PFN); firstly, calculating correlation coefficients based on gene expression profiles, and then filtering and clustering gene pairs using a parallel filtering method to obtain a fast planar filter network; (2) multi-scale clustering analysis; from the initial PFN of the connected components, multi-scale clustering of each parent cluster can obtain more sub-modules, followed by hierarchical clustering results; (3) downstream analysis, using multiscale hub analysis (MHA) to identify important hubs based on the network topology; (4) Finally, the correlation between clustering results and clinical information was analyzed by cluster-trait association analysis (CTA). To test whether selected modules in WGCN and multiscale hub genes in MEGCN were highly associated with tumorigenesis and malignant progression, we performed GO and KEGG pathway analysis on all genes in the selected DAVID (https://david.ncifcrf.gov/) database [44] for hub modules. P ≤ 0.05 was considered as the cut-off criterion for identifying enrichment.

Validation and survival analysis of pivotal genes

In WGCNA, when hub modules were identified, we selected candidate hub genes by module connectivity (cor. gene module-Membership (MM) > 0.8) and clinical feature correlation (cor.geneTraitSignificance (GS) > 0.2), both by Pearson correlation absolute values were used to determine. At the same time, we chose the STRING database (https://www.string-db.org/) [45] to identify the connectivity of candidate genes. Next, Cytoscape's plug-in CytoHubba was used to identify top50 genes as candidate genes [46]. The gene modules from the WGCNA and MEGENA analyses are saved separately. The userListEnrichment R function was used to calculate the overlap of genes between the WGCNA and MEGENA modules. The WGCNA and MEGENA modules that overlapped significantly and were significantly associated with PCa were retained. In this study, we used the Gene Expression Profiling Interactive Analysis (GEPIA) dataset (http://gepia.cancer-pku.cn/) for prognostic analysis of PCa[47]. We analyzed the differential expression of key hub genes in PCa and their association with RFS using GEPIA. P-value and fold change were defined as 0.05 and group limit of 50% for two survival analyses. Also, normal distribution tests and differential analysis of hub genes in normal versus tumor tissue and normal versus different tumor subtypes were performed based on the Wilcoxon test. Finally, Human Protein Atlas Dataset (https://www.proteinatlas.org/) [48] was used to validate the immunohistochemistry (IHC) of the hub gene.

Development and validation of the hub gene prognostic model

To assess the prognostic value of hub genes, Cox regression analysis was used to evaluate the correlation between genes and survival status in a cohort of 494 PCas. Next, we chose the glmnet R language package for LASSO Cox regression modeling to narrow down the candidate genes and build prognostic models [49]. Ultimately, the screening retained three genes and their coefficients, and the penalty parameter λ was determined based on the minimum criterion. The algorithm obtained a more accurate model by constructing a penalty function. As a method for complex covariance data, variable selection can be achieved and parameter estimation to better address multiple covariances in regression analysis. The PCa expression data were centrally normalized using a scale function, and a risk score was calculated.

Patients with PCa OC were divided into low and high-risk subgroups based on median risk scores. Kaplan Meier analysis was performed to compare OS time assessment risk scores and overall survival between the two subgroups. To validate the accuracy and predictive power of the model, ROC curve analyses were performed using the survivor, survminer, and timeROC R packages for one year, three years, and five years, respectively, to calculate the area under The curve (AUC) and to compare the effect of the timeROC R package on classifier performance.

Table 1

Information on the clinical features of patients with prostate cancer
VARIABLE	TYPE
PCA (260)
Tumor	232	89%
Normol	28	11%
CCF
Cluster1	41	16%
Cluster2	68	26%
Cluster3	72	28%
Cluster4	51	20%
Cluster5	28	10%
Relaplse
NA	13	5%
NO	206	79%
YES	41	16%
Age
>=60	163	63%
< 60	97	37%
Gleason_score
6	17	6%
7	144	55%
8	28	11%
9	69	27%
10	2	1%
OS
0 < = OS < = 1000	142	55%
1000 < OS < = 2000	84	32%
2000 < OS < = 3000	24	9%
3000 < OS < = 4000	5	2%
4000 < OS < = 5000	3	1%
OS > = 5000	2	1%

Relationship between key genes and tumor-infiltrating immune cells

The Tumor Immunity Evaluation Resource (TIMER) (https://cistrome.shinyapps.io/- -timer/) is a comprehensive resource for detecting immune cell infiltration in tumor tissue using RNA-Seq expression profiling data and contains immune infiltration in different cancer types[50]. Hub gene expression was correlated with six types of immune infiltration (B cells, CD4 + T cells, CD8 + T cells, neutrophils, macrophages, and dendritic cells) were correlated and assessed using the Gene module. The SCNA module was also used to explore the correlation between somatic cell copy number changes and the abundance of immune infiltrates.

Exploring the drug sensitivity of the hub gene

The GSCAlite (http://bioinfo.life.hust.edu.cn/web/GSCALite/) database is a genomic cancer analysis platform[51]. The database can be used for genomic and immunogenomic analyses. It also enables researchers to combine clinical information and small molecule drugs to mine candidate biomarkers and valuable small molecule drugs for better experimental design and further clinical trials. The GSCA data contains 33 cancer types from TCGA and normal tissue data from GTEx, and over 750 small-molecule drugs from GDSC and CTRP for 10,000 genomic data. Spearman correlations represent gene expression associated with drugs. A positive correlation means that a gene with high expression is resistant to the drug and vice versa.

Combined consistency clustering and principal component analysis to obtain sample cohorts

The analytical process used in this study is shown in Figure S1. To follow up on the molecular heterogeneity of PCa, we performed an unsupervised consensus analysis on all samples. In this study, we chose k = 5, which allowed us to divide all samples into five groups. Among them, C1:59, C2:146, C3:129, C4:102, and C5:114. Moreover, we found that C1-C4 were tumor subgroups, and C5 were normal samples (Fig. 1A-C). Next, to verify the robustness of this classification, we performed another principal component analysis based on the results of consistent clustering and observed the subgroup differentiation. PC1 and PC2 were selected as the main components for the analysis, and 260 samples with significant differences were obtained after excluding outlier samples as well as those with insignificant differences. As shown in Fig. 1D, there was an individual crossover between tumor subgroups, indicating good differentiation between subgroups. The clustering results were as follows: C1:41, C2:68, C3:72, C4:51, and C5:28, reflecting the impact of differences between tumor subgroups and normal samples on transcriptional profiles.

Identification of differentially expressed genes and GO functional annotation in PCa

We performed differential gene expression analysis of the four tumor subgroups separately from normal samples. A total of 3521 differentially expressed genes were obtained. Of these, 581 up-regulated and 1432 down-regulated genes were obtained for C1 vs. C5, 231 up-regulated and 637 down-regulated genes were obtained for C2 vs. C5, 234 up-regulated and 684 down-regulated genes were obtained for C3 vs. C5, and 282 up-regulated and 2265 down-regulated genes were obtained for C4 vs. C5 (Fig. 2A-B). Comparisons between tumor subgroups revealed 842 DEGs specific to C1, 47 DEGs specific to C2, 29 DEGs specific to C3, and 1045 DEGs specific to C4. Three hundred ninety-six overlap genes were also identified, suggesting a degree of similarity in expression profiles between tumor subgroups but also a high degree of heterogeneity between tumor subgroups (Fig. 2C-D).

Functional enrichment analysis of four groups of differentially expressed genes

The researchers used the Metascape website to carry out a biological process analysis of the GO function of 3521 DEGs. As shown in Fig. 3A, there was a strong interaction between the four groups of differentially expressed genes. Also, DEGs were found to be mainly enriched in functions such as leukocyte migration, inflammatory response, regulation of MAPK cascade, circulatory system processes cell component morphogenesis, and regulation of cell adhesion. Moreover, Metascape enrichment analysis also showed that DEGs differed significantly between the four groups in the pathways regulating cell adhesion and leukocyte migration, demonstrating genetic and functional differences in tumor subtypes (Fig. 3B-E).

Co-expression network analysis reveals key modules of PCa network

We used the expression profiles and clinical information of 260 samples from HCC as input files for WGCNA by integrating them. We eliminated two outlier samples based on the clustering tree and then plotted the sample dendrograms and expression heat maps of the traits (Figure. 3A). After excluding the two outlier samples, the remaining 258 samples and 3521 differentially expressed genes were subjected to WGCNA analysis. We chose a soft threshold power of 8 (scale-free R^2 = 0.85) by calculation, ensuring that the scale-free network was reasonable (Figure. 3B). We built the network by choosing the soft threshold. And the minimum number of genes in the modules was required to be 50. Subsequently, the initial modules were divided using a dynamic tree, and the genes were divided into modules based on the similarity of the characteristic genes in the gene clustering map, and ten gene modules were finally identified (Fig. 3C). One color indicates a gene module, while grey gene modules are considered invalid genes that cannot be assigned to a module. Next, we created a heat map of the correlation between gene modules and clinical information (Fig. 4D). We found that the blue module was found to have the highest correlation coefficient with consensus clustering (Pearson cor = 0.71, a significant positive correlation). This suggests that the genes in this module correlate with tumor typing as well as malignancy; meanwhile, we found a large correlation coefficient between the green module and tumorigenesis (Pearson cor = 0.53, significant positive correlation). In contrast, the absolute values of the correlation coefficients between genes and age, recurrence status, and gleason_score were smaller in the other modules, suggesting that the genes are less relevant to other clinical information. Finally, we calculated the connectivity of the genes within the blue and green modules and selected the top50 most closely connected genes for subsequent analysis.

Meanwhile, the expression profiles of 3521 differentially expressed genes (DEGs) were applied to the multiscale embedded gene co-expression network (MEGENA) analysis. After parallelization, early termination, and pre-quality control steps, 10,333 gene pairs of PFNs were used to enter 1755 MCAs for multiscale cluster identification. A total of 5 different scales were identified by multi-scale cluster analysis, and 226 significant modular clusters were constructed (Figure S4), where the modular p-values were all less than 0.05.

Table 2

Functional annotation of GO_BP for top20 of four sets of DEGs.
NO.	Category	Term ID	Description	LogP	Ratio
1	GO_BP	GO:0003013	circulatory system process	-71.68	207/611
2	GO_BP	GO:0032989	cellular component morphogenesis	-65.57	226/766
3	GO_BP	GO:0061061	muscle structure development	-56.30	186/609
4	GO_BP	GO:0030029	actin filament-based process	-55.21	215/794
5	GO_BP	GO:0045229	external encapsulating structure organization	-52.90	143/398
6	GO_BP	GO:0034330	cell junction organization	-51.96	195/701
7	GO_BP	GO:0040017	positive regulation of locomotion	-48.52	171/587
8	GO_BP	GO:0048729	tissue morphogenesis	-46.92	177/638
9	GO_BP	GO:0070848	response to growth factor	-46.55	192/736
10	GO_BP	GO:0006954	inflammatory response	-45.06	196/778
11	GO_BP	GO:0043408	regulation of MAPK cascade	-44.56	185/712
12	GO_BP	GO:0009611	response to wounding	-41.77	172/658
13	GO_BP	GO:0050801	ion homeostasis	-40.81	190/787
14	GO_BP	GO:0030155	regulation of cell adhesion	-40.19	181/734
15	GO_BP	GO:0048598	embryonic morphogenesis	-39.35	154/569
16	GO_BP	GO:0034762	regulation of transmembrane transport	-38.50	152/565
17	GO_BP	GO:0001501	skeletal system development	-37.11	137/486
18	GO_BP	GO:0050900	leukocyte migration	-36.41	140/511
19	GO_BP	GO:0048871	multicellular organismal homeostasi	-35.28	139/516
20	GO_BP	GO:0060485	mesenchyme development	-34.80	99/288

Comparing the significant modules of WGCNA with the MEGENA modules, we found a degree of overlap between the key gene modules of WGCNA and those of MEGENA. We found a large degree of overlap between the blue module (Fig. 5A) and green (Fig. 5B) module of WGCNA and the C1_5 module (Fig. 6A) and C1_7 (Fig. 6B) module of MEGENA (Fig. 6A). We also annotated the genes of the four modules and found that the four groups of genes are mainly involved in the cell cycle, cell division, and transcriptional functions (Figure S5).WGCNA and MEGENA can complement each other and help us to preserve the maximum range of gene modules specific to PCa. We found that the expression profiles of this class of overlapping genes differed greatly between tumor and paraneoplastic tissues, and some of the genes also differed greatly between different subtypes of tumors (Fig. 7A-B). In addition, we found that NDUFB7, AURKAIP1, ROMO1, SCAND1, GADD45GIP1, and NDUFA13 tumors all showed high expression between the five subtypes of these genes (Fig. 7C). Next, we annotated the GO_BP function of 29 overlapping genes based on the DAVID database and found that these genes are mainly involved in cell division, cell cycle, mitosis, and positive regulation of ubiquitin-protein ligase activity (Fig. 7D).

Validation of key candidate genes

We first performed a disease-free survival analysis of key genes using the GEPIA database to explore the relationship between genes and phenotype and prognostic status. We selected median as a criterion to differentiate gene expression levels and found that TPX2, ROMO1, PLK1, UBE2C, and KIF4A were considered hub genes with the most significant p-values selected from the 29 core genes (Fig. 8A). We then perform a Wilcoxon test on the five hub genes screened, and we found significant differences (p < 0.05) in the expression of the five genes between tumor and para-cancer tissues (Fig. 8B). Meanwhile, we further analyzed the expression differences of the five hub genes among different tumor subtypes, and we found that these genes also differed significantly among tumor subtypes of PCa (Fig. 8C). Moreover, we found that ROMO1 was highly expressed in both tumor and normal tissues after immunohistochemistry (Fig. 8D).

Construction and validation of prognostic gene markers in TCGA

The data for 497 tumors were obtained from The Cancer Genome Atlas (TCGA) dataset, and corresponding clinical information was incorporated into the prognostic training model. The expression data of five genes, TPX2, ROMO1, PLK1, UBE2C, and KIF4A, were pooled, and LASSO Cox was applied to construct a regression model. The model identified ROMO1, PLK1, and KIF4A5 based on the optimal value of λ to construct an OS prognostic model for liver cancer patients (Fig. 8A, B). The risk scores, patient survival status distributions, and gene expression profiles associated with the three genetic features in the training dataset are shown in Fig. 8C. OS was significantly lower in the high-risk group than in the low-risk group (p = 0.0021) (Fig. 8D). Finally, ROC analysis was performed, and the area under the ROC curve (AUC) for overall survival (OS) at 1, 3, and 5 years was 0.789 and 0.72 and 0.626, respectively (Fig. 8E). In a subsequent study, we also did a LASSO regression analysis on ROMO1 and found that ROMO1 could also be used as an independent prognostic model and that its low-risk group had a better survival outcome than the high-risk group (Figure S3).

Alterations in the immune microenvironment may result from aberrant regulation of key genes in PCa

To explore the relationship between hub genes expression and tumor-infiltrating immune cells, we used the TIMER tool to analyze the correlation between key gene expression and the level of immune infiltration. Our results showed that ROMO1 expression was strongly correlated with immunity to PCa and was found to be negatively correlated with B cells, CD8 + T cells, macrophages, neutrophils, dendritic cells and positively correlated with CD4 + T cells; PLK1 and KLF4A expression was positively correlated with B cells, CD8 + T cells, macrophages, neutrophils, dendritic cells (Fig. 11A). In addition, we found that different mutant forms of ROMO1 correlated with immune infiltration of CD8 + T cells, macrophages, and neutrophils; different mutant forms of KLF4A correlated with immune levels of CD8 + T cells, macrophages, neutrophils, and dendritic cells; and we found that ROMO1, PLK1, and KLF4A were found to have higher levels of immune infiltration in dendritic cells (Fig. 11B). This also reveals the influence of key genes on the immune microenvironment of PCa.

Drug sensitivity and resistance of key genes

To investigate the relationship between the expression of hub genes and drug interactions, we used the GSCA database to screen for drugs that interacted with it. The GSCA database integrated drug sensitivity and gene expression profile data for cancer cell lines from both the GDSC and CTRP databases into GSCALite, and then, using Spearman correlation analysis, the expression levels of genes were correlated with the IC50. FDR < 0.05 was considered significant. During the research, we found that ROMO1 expression was significantly positively correlated with nutlin-3, fluorouracil, etc. KIF4A expression was significantly positively correlated with PD318088 and selumetinib, and negatively correlated with ML239; PLK1 expression was significantly negatively correlated with STF − 31 was significantly negatively correlated (Fig. 12). This reveals a diversity of expression levels of key genes and drug sensitivity.

In recent years, PCa has become the most common and deadly solid cancer and genitourinary tumor in men worldwide [52], with a diagnosis rate of 12% and a mortality rate of 9%, and its incidence rises with age [53]. The traditional clinical treatment options include resection, radiotherapy, chemotherapy, and endocrine therapy [54]. The tumor heterogeneity leads to limitations in conventional treatment and makes it difficult to effectively manage and risk assess patients [55, 56]. Recent studies have demonstrated the imperative of identifying new effective biomarkers and promising immune-related therapeutic targets that can be used to guide cancer treatment in clinical practice [57]. In this study, we found a high correlation between ROMO1 expression and clinical features such as the occurrence of PCa, subtypes classified by consistent clustering. Also, ROMO1 correlated with the level of immune cell infiltration and immune pathways in PCa.

Table 3

Functional annotation of GO_BP for candidate genes.
NO.	Term Category_BP	Count	PValue	Fold.Enrichment	FDR
1	GO:0006468:protein phosphorylation	6	0.00	8.18	0.14
2	GO:0006366:transcription from RNA polymerase II promoter	5	0.00	6.06	0.64
3	GO:1904668:positive regulation of ubiquitin protein ligase activity	2	0.01	138.21	0.62
4	GO:0030071:regulation of mitotic metaphase/anaphase transition	2	0.01	138.21	0.64
5	GO:0018105:peptidyl-serine phosphorylation	3	0.02	14.93	0.64
6	GO:0051301:cell division	4	0.02	7.11	0.64
7	GO:0000086:G2/M transition of mitotic cell cycle	3	0.02	13.62	0.64
8	GO:0006355:regulation of transcription, DNA-templated	7	0.02	2.89	0.69
9	GO:0045862:positive regulation of proteolysis	2	0.02	65.47	0.69
10	GO:0031572:G2 DNA damage checkpoint	2	0.03	62.19	0.69
11	GO:0045736:negative regulation of cyclin-dependent protein serine/threonine kinase activity	2	0.03	56.54	0.69
12	GO:0051439:regulation of ubiquitin-protein ligase activity involved in mitotic cell cycle	2	0.05	54.08	0.69
13	GO:0045892:negative regulation of transcription, DNA-templated	4	0.04	4.99	0.74
14	GO:0000281:mitotic cytokinesis	2	0.04	42.89	0.74

Firstly, we selected consistent clustering and principal component analysis to identify 260 valid samples containing one normal group and four tumor subtypes. Secondly, the limmaR package was selected to differentially analyze the four groups of tumor subtypes from the normal group, respectively. Finally, 3521 DEGs were identified and functionally annotated, and these differentially expressed genes were found to be mainly enriched in ah immune-related pathways. Next, WGCNA and MEGENA were applied to construct co-expression networks of differentially expressed genes to identify gene modules associated with tumorigenesis and heterogeneity, and 29 overlapping genes were identified by calculation. Then, the Lasso Cox regression model was constructed to identify ROMO1, PLK1, and KIF4A5 as the optimal core genes. Subsequently, when analyzing the relationship between the expression of hub genes and tumor-infiltrating immune cells, it was found that the expression of ROMO1, PLK1, and KIF4A were all associated with tumor-infiltrating immune cells, which in turn led to the alteration of the tumor microenvironment alterations, which in turn lead to tumor heterogeneity. Finally, we found that nutlin-3, fluorouracil, and others could be potential therapeutic agents for ROMO1, and PD318088 and selumetinib could be potential therapeutic agents for PLK1. In conclusion, our study provides a new theoretical basis for the diagnosis and treatment of PCa.

Polo-like kinase 1 (PLK1) is a member of a family of serine/threonine protein kinases that are widely found in eukaryotic cells [58]. Its specific functions are mainly in cell cycle processes, including controlling mitotic entry and the G2/M checkpoint, coordinating centrosomes and the cell cycle, regulating spindle assembly and chromosome segregation, performing multiple functions during mid-spindle and abscission, promoting DNA replication, participating in the cytoplasmic division and meiosis, and playing an important role in the initiation, maintenance, and completion of mitosis [59–61]. The pharmacological inhibition of PLK1 in triple-negative breast cancer has been reported to increase the anti-proliferative activity of drug-resistant cells, which in turn causes G2/M phase block and increasing the phosphorylation of cell cycle proteins inducing apoptosis.[62]. Also, in PCa, mitotic kinase polo-like kinase 1 (PLK1) is expressed at elevated levels and is associated with tumor grade [63]. In the present study, we also found a significant prognostic effect of PLK1, with expression in different subtypes of the prostate gland also differing significantly from normal samples, and found that PLK1 expression was also associated with tumor-infiltrating immune cells, further demonstrating the reliability of PLK1 as a biomarker for PCa.

Table 4

Drugs that interact with the hub gene.
NO.	SYMBOL	DRUE	COR	FDR
1	KIF4A	Compound 23 citrate	0.12	0.00
2	KIF4A	ML239	-0.15	0.00
3	KIF4A	PD318088	0.17	0.00
4	KIF4A	SID 26681509	-0.12	0.01
5	KIF4A	STF-31	-0.11	0.01
6	KIF4A	VAF-347	0.12	0.05
7	KIF4A	YM-155	-0.14	0.00
8	KIF4A	erlotinib	0.13	0.00
9	KIF4A	fluorouracil	0.12	0.00
10	KIF4A	lapatinib	0.14	0.00
11	KIF4A	necrosulfonamide	-0.11	0.02
12	KIF4A	niclosamide	-0.13	0.00
13	KIF4A	nutlin-3	0.13	0.00
14	KIF4A	selumetinib	0.17	0.00
15	KIF4A	serdemetan	0.11	0.01
16	KIF4A	trametinib	0.16	0.02
17	KIF4A	trifluoperazine	0.12	0.02
18	KIF4A	vandetanib	0.10	0.03
19	PLK1	Compound 23 citrate	0.03	0.47
20	PLK1	ML239	-0.08	0.05
21	PLK1	PD318088	0.06	0.23
22	PLK1	SID 26681509	-0.06	0.23
23	PLK1	STF-31	-0.17	0.00
24	PLK1	VAF-347	0.07	0.28
25	PLK1	YM-155	-0.09	0.044
26	PLK1	erlotinib	0.01	0.91
27	PLK1	fluorouracil	-0.08	0.04
28	PLK1	lapatinib	0.05	0.29
29	PLK1	necrosulfonamide	-0.04	0.41
30	PLK1	niclosamide	-0.10	0.016
31	PLK1	nutlin-3	0.05	0.25
32	PLK1	selumetinib	0.07	0.14
33	PLK1	serdemetan	0.01	0.81
34	PLK1	trametinib	0.13	0.07
35	PLK1	trifluoperazine	0.06	0.24
36	PLK1	vandetanib	0.00	0.99
37	ROMO1	Compound 23 citrate	0.12	0.00
38	ROMO1	ML239	0.02	0.71
39	ROMO1	PD318088	0.10	0.02
40	ROMO1	SID 26681509	0.02	0.63
41	ROMO1	STF-31	0.04	0.38
42	ROMO1	VAF-347	-0.09	0.13
43	ROMO1	YM-155	-0.02	0.78
44	ROMO1	erlotinib	0.12	0.00
45	ROMO1	fluorouracil	0.16	0.00
46	ROMO1	lapatinib	0.11	0.01
47	ROMO1	necrosulfonamide	0.11	0.02
48	ROMO1	niclosamide	-0.05	0.23
49	ROMO1	nutlin-3	0.18	2.01
50	ROMO1	selumetinib	0.10	0.03
51	ROMO1	serdemetan	0.12	0.00
52	ROMO1	trametinib	0.07	0.39
53	ROMO1	trifluoperazine	0.13	0.00
54	ROMO1	vandetanib	0.05	0.32

Kinesin superfamily protein 4A (KIF4A) is found in all eukaryotes and belongs to a family of KIFs that are highly conserved [64]. KIF4A has important roles in DNA repair, DNA replication, spindle organization, cytoplasmic division, and intracellular transport [65]. KIF4A has been previously reported to be aberrantly expressed in many cancers, revealing its function and role in different tumors [66–68]. KIF4A can promote PCa cell growth through AR and AR-V7-dependent signaling [66]. In the present study, KIF4A expression in different subtypes of the prostate was also significantly different from normal samples, and KIF4A expression was also found to be associated with tumor-infiltrating immune cells.

Reactive oxygen species (ROS) modulator 1 (ROMO1) is a membrane protein found in mitochondria that is important for regulating mitochondrial ROS production and redox sensing [69]. Romo1 is capable of triggering and exacerbating cancer through extracellular signal-regulated kinases (ERKs) and nuclear factor-kB (NF-kB)-induced reactive oxygen species (ROS) [70]. ROS can trigger and exacerbate cancer in a variety of malignancies [71, 72]. It can also influence cancer cell invasion by affecting the epithelial-mesenchymal transition (EMT) pathway [73, 74]. This protein affects signaling pathways and ROS homeostasis, affects the G2/M phase cell cycle, and leads to cell overgrowth through increased expression levels[75]. Studies have shown that ROMO1 is overexpressed in hepatocellular carcinoma, colorectal cancer, and glioma [76–78] but has not been reported in PCa. In the present study, we selected the Wilcoxon test and found that ROMO1 was highly expressed in tumor tissue and significantly different from normal tissue; we also found that the four identified tumor subtypes were significantly different. The expression of ROMO1 was also found to be associated with tumor-infiltrating immune cells, leading to changes in the tumor microenvironment and further increasing tumor heterogeneity; drug sensitivity analysis revealed that nutlin-3 and fluorouracil could be used as potential therapeutic agents for ROMO1.

In general, tumor pathogenesis involves many interacting signaling pathways, including tumor cell proliferation, cell immortalization, invasion, and migration, etc[79]. The complexity of cancer can be reflected through the tumor microenvironment, while protein interactions can further increase heterogeneity between tumors [80, 81]. In the present study, we have used weighted gene co-expression network analysis (WGCNA), which classifies the gene co-expression network of PCa into 10 highly correlated signature modules. The modules were then correlated with specific clinical features to identify genes that are key to tumorigenesis and transformation, to help identify potential mechanisms involved, and to explore candidate biomarkers. However, WGCNA has the limitation of not being able to coexist at different levels of clustering within a single network, thus not reflecting the multi-scale hierarchical nature of complex networks. Multi-scale embedded gene co-expression network analysis (MEGENA), on the other hand, allows the construction and analysis of large-scale planar filtered co-expression networks to the greatest extent possible [43]. Parallelization of embedded network construction and the identification of multiscale clustering structures are two key components of MEGENA, which is an essential complement to existing co-expression network analysis methods by identifying multiscale modular systems and co-expression networks with varying degrees of sparse and tight connectivity. Here, by using WGCNA in conjunction with MEGENA to construct a gene co-expression network for PCa, we identified more meaningful clusters of co-expressed genes and identified key biomarkers associated with prostate carcinogenesis in transformation.

In this study, we combined various bioinformatic analysis methods, especially the introduction of weighted gene co-expression network analysis and multi-scale chimeric network analysis, to reveal that ROMO1 may serve as a new key prognostic marker for PCa. However, the article still has some limitations. Firstly, the role of ROMO1 has not been validated experimentally in vivo and in vitro, and secondly, the fact that the number of cancer samples and normal samples are not identical to each other has led to some preference in our data. We believe that if we had direct access to a larger sample of clinical sequencing data and sample information, we would obtain better and more accurate results.

In summary, we have performed a systematic bioinformatics analysis of differentially expressed genes between normal samples and samples with different PCa subtypes and identified ROMO1, a key gene closely associated with PCa development and progression. This provides a more detailed molecular mechanism for the development and progression of PCa. We believe that ROMO1 is expected to serve as a potential biomarker and therapeutic We suggest that ROMO1 has potential as a potential biomarker and therapeutic target for PCa. We suggest that ROMO1 may be a potential biomarker and therapeutic target for PCa, providing new theoretical support for the diagnosis prognosis of PCa.

PCa: Prostate cancer; WGCNA: Weighted Gene Co-expression Network Analysis; MEGENA: Multiscale Embedded Gene Co-expression Network Analysis; DEG: Differentially expressed genes; PFN: Planar filtered network; MHA: Multiscale hub analysis; CTA: Cluster-trait association analysis; GS: Gene TraitSignificance; MM: Module-Membership; IHC: Immunohistochemistry.

Authors’ contributions

L.W., W.L., and X.S. conceived and designed the study. M.F., J.Y. and Y.M.. carried out experiments and statistical analysis. Y.Z., Y.M., J.L., R.L., X.Z., P.Z., L.Z., J.F., and Y.C., provided some critical revision; X.L., Y.W. and Z.L. wrote the manuscript. The author(s) read and approved the final manuscript

Funding

Not applicable.

Availability of data and materials

The datasets used and/or analysed during in the study are available from the corresponding author on reasonable request.

Ethics approval and consent to participate

Not applicable.

Acknowledgements

Not applicable.

Consent for publication

Not applicable.

Competing interests

All the authors declare that there are no conflicts of interest to disclose.

Author details

¹College of Life Sciences, Xinyang Normal University, Xinyang 464000, China ²College of Medicine, Xinyang Normal University, Xinyang 464000, China. ³Department of Computer Science, City University of Hong Kong, Hong Kong, China. ⁴Department of Recovery Medicine, People's Liberation Army 990 Hospital, Xinyang, China. ⁵Department of Urology，Central Hospital of Xinyang，Xinyang, China.

Ge R, Wang Z, Montironi R, Jiang Z, Cheng M, Santoni M, Huang K, Massari F, Lu X, Cimadamore A, Lopez-Beltran A, Cheng L (2020) Epigenetic modulations and lineage plasticity in advanced prostate cancer. Ann Oncol 31 (4):470-479. doi:10.1016/j.annonc.2020.02.002
Auchus RJ, Sharifi N (2020) Sex Hormones and Prostate Cancer. Annu Rev Med 71:33-45. doi:10.1146/annurev-med-051418-060357
Chen RC, Basak R, Meyer AM, Kuo TM, Carpenter WR, Agans RP, Broughman JR, Reeve BB, Nielsen ME, Usinger DS, Spearman KC, Walden S, Kaleel D, Anderson M, Sturmer T, Godley PA (2017) Association Between Choice of Radical Prostatectomy, External Beam Radiotherapy, Brachytherapy, or Active Surveillance and Patient-Reported Quality of Life Among Men With Localized Prostate Cancer. JAMA 317 (11):1141-1150. doi:10.1001/jama.2017.1652
Siegel DA, O'Neil ME, Richards TB, Dowling NF, Weir HK (2020) Prostate Cancer Incidence and Survival, by Stage and Race/Ethnicity - United States, 2001-2017. MMWR Morb Mortal Wkly Rep 69 (41):1473-1480. doi:10.15585/mmwr.mm6941a1
Bourke L, Smith D, Steed L, Hooper R, Carter A, Catto J, Albertsen PC, Tombal B, Payne HA, Rosario DJ (2016) Exercise for Men with Prostate Cancer: A Systematic Review and Meta-analysis. Eur Urol 69 (4):693-703. doi:10.1016/j.eururo.2015.10.047
Venkateswaran V, Klotz LH (2010) Diet and prostate cancer: mechanisms of action and implications for chemoprevention. Nat Rev Urol 7 (8):442-453. doi:10.1038/nrurol.2010.102
Costello AJ (2020) Considering the role of radical prostatectomy in 21st century prostate cancer care. Nat Rev Urol 17 (3):177-188. doi:10.1038/s41585-020-0287-y
Wallis CJD, Saskin R, Choo R, Herschorn S, Kodama RT, Satkunasivam R, Shah PS, Danjoux C, Nam RK (2016) Surgery Versus Radiotherapy for Clinically-localized Prostate Cancer: A Systematic Review and Meta-analysis. Eur Urol 70 (1):21-30. doi:10.1016/j.eururo.2015.11.010
Zaorsky NG, Davis BJ, Nguyen PL, Showalter TN, Hoskin PJ, Yoshioka Y, Morton GC, Horwitz EM (2017) The evolution of brachytherapy for prostate cancer. Nat Rev Urol 14 (7):415-439. doi:10.1038/nrurol.2017.76
Mottet N, Bellmunt J, Bolla M, Briers E, Cumberbatch MG, De Santis M, Fossati N, Gross T, Henry AM, Joniau S, Lam TB, Mason MD, Matveev VB, Moldovan PC, van den Bergh RCN, Van den Broeck T, van der Poel HG, van der Kwast TH, Rouviere O, Schoots IG, Wiegel T, Cornford P (2017) EAU-ESTRO-SIOG Guidelines on Prostate Cancer. Part 1: Screening, Diagnosis, and Local Treatment with Curative Intent. Eur Urol 71 (4):618-629. doi:10.1016/j.eururo.2016.08.003
Widmark A, Klepp O, Solberg A, Damber JE, Angelsen A, Fransson P, Lund JA, Tasdemir I, Hoyer M, Wiklund F, Fossa SD, Scandinavian Prostate Cancer Group S, Swedish Association for Urological O (2009) Endocrine treatment, with or without radiotherapy, in locally advanced prostate cancer (SPCG-7/SFUO-3): an open randomised phase III trial. Lancet 373 (9660):301-308. doi:10.1016/S0140-6736(08)61815-2
Aparicio AM, Harzstark AL, Corn PG, Wen S, Araujo JC, Tu SM, Pagliaro LC, Kim J, Millikan RE, Ryan C, Tannir NM, Zurita AJ, Mathew P, Arap W, Troncoso P, Thall PF, Logothetis CJ (2013) Platinum-based chemotherapy for variant castrate-resistant prostate cancer. Clin Cancer Res 19 (13):3621-3630. doi:10.1158/1078-0432.CCR-12-3791
Balk SP, Ko YJ, Bubley GJ (2003) Biology of prostate-specific antigen. J Clin Oncol 21 (2):383-391. doi:10.1200/JCO.2003.02.083
Muenst S, Laubli H, Soysal SD, Zippelius A, Tzankov A, Hoeller S (2016) The immune system and cancer evasion strategies: therapeutic concepts. J Intern Med 279 (6):541-562. doi:10.1111/joim.12470
Perez CR, De Palma M (2019) Engineering dendritic cell vaccines to improve cancer immunotherapy. Nat Commun 10 (1):5408. doi:10.1038/s41467-019-13368-y
McLaughlin J, Han G, Schalper KA, Carvajal-Hausdorf D, Pelekanou V, Rehman J, Velcheti V, Herbst R, LoRusso P, Rimm DL (2016) Quantitative Assessment of the Heterogeneity of PD-L1 Expression in Non-Small-Cell Lung Cancer. JAMA Oncol 2 (1):46-54. doi:10.1001/jamaoncol.2015.3638
Ebert LM, Yu W, Gargett T, Brown MP (2018) Logic-gated approaches to extend the utility of chimeric antigen receptor T-cell technology. Biochem Soc Trans 46 (2):391-401. doi:10.1042/BST20170178
Dagogo-Jack I, Shaw AT (2018) Tumour heterogeneity and resistance to cancer therapies. Nat Rev Clin Oncol 15 (2):81-94. doi:10.1038/nrclinonc.2017.166
Cacchione S, Biroccio A, Rizzo A (2019) Emerging roles of telomeric chromatin alterations in cancer. J Exp Clin Cancer Res 38 (1):21. doi:10.1186/s13046-019-1030-5
Kuppers R, Stevenson FK (2018) Critical influences on the pathogenesis of follicular lymphoma. Blood 131 (21):2297-2306. doi:10.1182/blood-2017-11-764365
Martincorena I, Campbell PJ (2015) Somatic mutation in cancer and normal cells. Science 349 (6255):1483-1489. doi:10.1126/science.aab4082
Sahin U, Tureci O (2018) Personalized vaccines for cancer immunotherapy. Science 359 (6382):1355-1360. doi:10.1126/science.aar7112
Lee E, Yong RL, Paddison P, Zhu J (2018) Comparison of glioblastoma (GBM) molecular classification methods. Semin Cancer Biol 53:201-211. doi:10.1016/j.semcancer.2018.07.006
Wang Q, Hu B, Hu X, Kim H, Squatrito M, Scarpace L, deCarvalho AC, Lyu S, Li P, Li Y, Barthel F, Cho HJ, Lin YH, Satani N, Martinez-Ledesma E, Zheng S, Chang E, Sauve CG, Olar A, Lan ZD, Finocchiaro G, Phillips JJ, Berger MS, Gabrusiewicz KR, Wang G, Eskilsson E, Hu J, Mikkelsen T, DePinho RA, Muller F, Heimberger AB, Sulman EP, Nam DH, Verhaak RGW (2017) Tumor Evolution of Glioma-Intrinsic Gene Expression Subtypes Associates with Immunological Changes in the Microenvironment. Cancer Cell 32 (1):42-56 e46. doi:10.1016/j.ccell.2017.06.003
Jackstadt R, van Hooff SR, Leach JD, Cortes-Lavaud X, Lohuis JO, Ridgway RA, Wouters VM, Roper J, Kendall TJ, Roxburgh CS, Horgan PG, Nixon C, Nourse C, Gunzer M, Clark W, Hedley A, Yilmaz OH, Rashid M, Bailey P, Biankin AV, Campbell AD, Adams DJ, Barry ST, Steele CW, Medema JP, Sansom OJ (2019) Epithelial NOTCH Signaling Rewires the Tumor Microenvironment of Colorectal Cancer to Drive Poor-Prognosis Subtypes and Metastasis. Cancer Cell 36 (3):319-336 e317. doi:10.1016/j.ccell.2019.08.003
Ulz P, Perakis S, Zhou Q, Moser T, Belic J, Lazzeri I, Wolfler A, Zebisch A, Gerger A, Pristauz G, Petru E, White B, Roberts CES, John JS, Schimek MG, Geigl JB, Bauernhofer T, Sill H, Bock C, Heitzer E, Speicher MR (2019) Inference of transcription factor binding from cell-free DNA enables tumor subtype prediction and early detection. Nat Commun 10 (1):4666. doi:10.1038/s41467-019-12714-4
Kris MG, Johnson BE, Berry LD, Kwiatkowski DJ, Iafrate AJ, Wistuba, II, Varella-Garcia M, Franklin WA, Aronson SL, Su PF, Shyr Y, Camidge DR, Sequist LV, Glisson BS, Khuri FR, Garon EB, Pao W, Rudin C, Schiller J, Haura EB, Socinski M, Shirai K, Chen H, Giaccone G, Ladanyi M, Kugler K, Minna JD, Bunn PA (2014) Using multiplexed assays of oncogenic drivers in lung cancers to select targeted drugs. JAMA 311 (19):1998-2006. doi:10.1001/jama.2014.3741
Chow RD, Guzman CD, Wang G, Schmidt F, Youngblood MW, Ye L, Errami Y, Dong MB, Martinez MA, Zhang S, Renauer P, Bilguvar K, Gunel M, Sharp PA, Zhang F, Platt RJ, Chen S (2017) AAV-mediated direct in vivo CRISPR screen identifies functional suppressors in glioblastoma. Nat Neurosci 20 (10):1329-1341. doi:10.1038/nn.4620
Liu Y, Yang M, Deng Y, Su G, Enninful A, Guo CC, Tebaldi T, Zhang D, Kim D, Bai Z, Norris E, Pan A, Li J, Xiao Y, Halene S, Fan R (2020) High-Spatial-Resolution Multi-Omics Sequencing via Deterministic Barcoding in Tissue. Cell 183 (6):1665-1681 e1618. doi:10.1016/j.cell.2020.10.026
Grindberg RV, Yee-Greenbaum JL, McConnell MJ, Novotny M, O'Shaughnessy AL, Lambert GM, Arauzo-Bravo MJ, Lee J, Fishman M, Robbins GE, Lin X, Venepally P, Badger JH, Galbraith DW, Gage FH, Lasken RS (2013) RNA-sequencing from single nuclei. Proc Natl Acad Sci U S A 110 (49):19802-19807. doi:10.1073/pnas.1319700110
Hutter C, Zenklusen JC (2018) The Cancer Genome Atlas: Creating Lasting Value beyond Its Data. Cell 173 (2):283-285. doi:10.1016/j.cell.2018.03.042
Tang Z, Li C, Kang B, Gao G, Li C, Zhang Z (2017) GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res 45 (W1):W98-W102. doi:10.1093/nar/gkx247
Thul PJ, Lindskog C (2018) The human protein atlas: A spatial map of the human proteome. Protein Sci 27 (1):233-244. doi:10.1002/pro.3307
Goldman MJ, Craft B, Hastie M, Repecka K, McDade F, Kamath A, Banerjee A, Luo Y, Rogers D, Brooks AN, Zhu J, Haussler D (2020) Visualizing and interpreting cancer genomics data via the Xena platform. Nat Biotechnol 38 (6):675-678. doi:10.1038/s41587-020-0546-8
Wilkerson MD, Hayes DN (2010) ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics 26 (12):1572-1573. doi:10.1093/bioinformatics/btq170
Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK (2015) limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res 43 (7):e47. doi:10.1093/nar/gkv007
Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK (2019) Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun 10 (1):1523. doi:10.1038/s41467-019-09234-6
Van Dam S, Vosa U, van der Graaf A, Franke L, de Magalhaes JP (2018) Gene co-expression analysis for functional classification and gene-disease predictions. Brief Bioinform 19 (4):575-592. doi:10.1093/bib/bbw139
Kang HJ, Kawasawa YI, Cheng F, Zhu Y, Xu X, Li M, Sousa AM, Pletikos M, Meyer KA, Sedmak G, Guennel T, Shin Y, Johnson MB, Krsnik Z, Mayer S, Fertuzinhos S, Umlauf S, Lisgo SN, Vortmeyer A, Weinberger DR, Mane S, Hyde TM, Huttner A, Reimers M, Kleinman JE, Sestan N (2011) Spatio-temporal transcriptome of the human brain. Nature 478 (7370):483-489. doi:10.1038/nature10523
Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang TH, Porta-Pardo E, Gao GF, Plaisier CL, Eddy JA, Ziv E, Culhane AC, Paull EO, Sivakumar IKA, Gentles AJ, Malhotra R, Farshidfar F, Colaprico A, Parker JS, Mose LE, Vo NS, Liu J, Liu Y, Rader J, Dhankani V, Reynolds SM, Bowlby R, Califano A, Cherniack AD, Anastassiou D, Bedognetti D, Mokrab Y, Newman AM, Rao A, Chen K, Krasnitz A, Hu H, Malta TM, Noushmehr H, Pedamallu CS, Bullman S, Ojesina AI, Lamb A, Zhou W, Shen H, Choueiri TK, Weinstein JN, Guinney J, Saltz J, Holt RA, Rabkin CS, Cancer Genome Atlas Research N, Lazar AJ, Serody JS, Demicco EG, Disis ML, Vincent BG, Shmulevich I (2018) The Immune Landscape of Cancer. Immunity 48 (4):812-830 e814. doi:10.1016/j.immuni.2018.03.023
Prasetyanti PR, Medema JP (2017) Intra-tumor heterogeneity from a cancer stem cell perspective. Mol Cancer 16 (1):41. doi:10.1186/s12943-017-0600-4
Langfelder P, Horvath S (2008) WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9:559. doi:10.1186/1471-2105-9-559
Song WM, Zhang B (2015) Multiscale Embedded Gene Co-expression Network Analysis. PLoS Comput Biol 11 (11):e1004574. doi:10.1371/journal.pcbi.1004574
Huang da W, Sherman BT, Lempicki RA (2009) Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc 4 (1):44-57. doi:10.1038/nprot.2008.211
Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, Simonovic M, Doncheva NT, Morris JH, Bork P, Jensen LJ, Mering CV (2019) STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res 47 (D1):D607-D613. doi:10.1093/nar/gky1131
Chin CH, Chen SH, Wu HH, Ho CW, Ko MT, Lin CY (2014) cytoHubba: identifying hub objects and sub-networks from complex interactome. BMC Syst Biol 8 Suppl 4:S11. doi:10.1186/1752-0509-8-S4-S11
Tang Z, Kang B, Li C, Chen T, Zhang Z (2019) GEPIA2: an enhanced web server for large-scale expression profiling and interactive analysis. Nucleic Acids Res 47 (W1):W556-W560. doi:10.1093/nar/gkz430
Colwill K, Renewable Protein Binder Working G, Graslund S (2011) A roadmap to generate renewable protein binders to the human proteome. Nat Methods 8 (7):551-558. doi:10.1038/nmeth.1607
Friedman J, Hastie T, Tibshirani R (2010) Regularization Paths for Generalized Linear Models via Coordinate Descent. J Stat Softw 33 (1):1-22
Li T, Fan J, Wang B, Traugh N, Chen Q, Liu JS, Li B, Liu XS (2017) TIMER: A Web Server for Comprehensive Analysis of Tumor-Infiltrating Immune Cells. Cancer Res 77 (21):e108-e110. doi:10.1158/0008-5472.CAN-17-0307
Liu CJ, Hu FF, Xia MX, Han L, Zhang Q, Guo AY (2018) GSCALite: a web server for gene set cancer analysis. Bioinformatics 34 (21):3771-3772. doi:10.1093/bioinformatics/bty411
Gorchakov AA, Kulemzin SV, Kochneva GV, Taranin AV (2020) Challenges and Prospects of Chimeric Antigen Receptor T-cell Therapy for Metastatic Prostate Cancer. Eur Urol 77 (3):299-308. doi:10.1016/j.eururo.2019.08.014
Trock BJ, Han M, Freedland SJ, Humphreys EB, DeWeese TL, Partin AW, Walsh PC (2008) Prostate cancer-specific survival following salvage radiotherapy vs observation in men with biochemical recurrence after radical prostatectomy. JAMA 299 (23):2760-2769. doi:10.1001/jama.299.23.2760
Mateo J, Fizazi K, Gillessen S, Heidenreich A, Perez-Lopez R, Oyen WJG, Shore N, Smith M, Sweeney C, Tombal B, Tomlins SA, de Bono JS (2019) Managing Nonmetastatic Castration-resistant Prostate Cancer. Eur Urol 75 (2):285-293. doi:10.1016/j.eururo.2018.07.035
Chang AJ, Autio KA, Roach M, 3rd, Scher HI (2014) High-risk prostate cancer-classification and therapy. Nat Rev Clin Oncol 11 (6):308-323. doi:10.1038/nrclinonc.2014.68
Wang G, Zhao D, Spring DJ, DePinho RA (2018) Genetics and biology of prostate cancer. Genes Dev 32 (17-18):1105-1140. doi:10.1101/gad.315739.118
Zhao SG, Lehrer J, Chang SL, Das R, Erho N, Liu Y, Sjostrom M, Den RB, Freedland SJ, Klein EA, Karnes RJ, Schaeffer EM, Xu M, Speers C, Nguyen PL, Ross AE, Chan JM, Cooperberg MR, Carroll PR, Davicioni E, Fong L, Spratt DE, Feng FY (2019) The Immune Landscape of Prostate Cancer and Nomination of PD-L2 as a Potential Therapeutic Target. J Natl Cancer Inst 111 (3):301-310. doi:10.1093/jnci/djy141
Steegmaier M, Hoffmann M, Baum A, Lenart P, Petronczki M, Krssak M, Gurtler U, Garin-Chesa P, Lieb S, Quant J, Grauert M, Adolf GR, Kraut N, Peters JM, Rettig WJ (2007) BI 2536, a potent and selective inhibitor of polo-like kinase 1, inhibits tumor growth in vivo. Curr Biol 17 (4):316-322. doi:10.1016/j.cub.2006.12.037
Frank T, Tuppi M, Hugle M, Dotsch V, van Wijk SJL, Fulda S (2019) Cell cycle arrest in mitosis promotes interferon-induced necroptosis. Cell Death Differ 26 (10):2046-2060. doi:10.1038/s41418-019-0298-5
Smits VA, Klompmaker R, Arnaud L, Rijksen G, Nigg EA, Medema RH (2000) Polo-like kinase-1 is a target of the DNA damage checkpoint. Nat Cell Biol 2 (9):672-676. doi:10.1038/35023629
Macurek L, Lindqvist A, Lim D, Lampson MA, Klompmaker R, Freire R, Clouin C, Taylor SS, Yaffe MB, Medema RH (2008) Polo-like kinase-1 is activated by aurora A to promote checkpoint recovery. Nature 455 (7209):119-123. doi:10.1038/nature07185
Nieto-Jimenez C, Galan-Moya EM, Corrales-Sanchez V, Noblejas-Lopez MDM, Burgos M, Domingo B, Montero JC, Gomez-Juarez M, Picazo-Martinez MG, Esparis-Ogando A, Pandiella A, Ocana A (2020) Inhibition of the mitotic kinase PLK1 overcomes therapeutic resistance to BET inhibitors in triple negative breast cancer. Cancer Lett 491:50-59. doi:10.1016/j.canlet.2020.06.020
Zhang Z, Hou X, Shao C, Li J, Cheng JX, Kuang S, Ahmad N, Ratliff T, Liu X (2014) Plk1 inhibition enhances the efficacy of androgen signaling blockade in castration-resistant prostate cancer. Cancer Res 74 (22):6635-6647. doi:10.1158/0008-5472.CAN-14-1916
Miki H, Okada Y, Hirokawa N (2005) Analysis of the kinesin superfamily: insights into structure and function. Trends Cell Biol 15 (9):467-476. doi:10.1016/j.tcb.2005.07.006
Sheng L, Hao SL, Yang WX, Sun Y (2018) The multiple functions of kinesin-4 family motor protein KIF4 and its clinical potential. Gene 678:90-99. doi:10.1016/j.gene.2018.08.005
Cao Q, Song Z, Ruan H, Wang C, Yang X, Bao L, Wang K, Cheng G, Xu T, Xiao W, Xiong Z, Liu D, Yang M, Zhou D, Yang H, Chen K, Zhang X (2020) Targeting the KIF4A/AR Axis to Reverse Endocrine Therapy Resistance in Castration-resistant Prostate Cancer. Clin Cancer Res 26 (6):1516-1528. doi:10.1158/1078-0432.CCR-19-0396
Hu G, Yan Z, Zhang C, Cheng M, Yan Y, Wang Y, Deng L, Lu Q, Luo S (2019) FOXM1 promotes hepatocellular carcinoma progression by regulating KIF4A expression. J Exp Clin Cancer Res 38 (1):188. doi:10.1186/s13046-019-1202-3
Taniwaki M, Takano A, Ishikawa N, Yasui W, Inai K, Nishimura H, Tsuchiya E, Kohno N, Nakamura Y, Daigo Y (2007) Activation of KIF4A as a prognostic biomarker and therapeutic target for lung cancer. Clin Cancer Res 13 (22 Pt 1):6624-6631. doi:10.1158/1078-0432.CCR-07-1328
Wang P, Geng J, Gao J, Zhao H, Li J, Shi Y, Yang B, Xiao C, Linghu Y, Sun X, Chen X, Hong L, Qin F, Li X, Yu JS, You H, Yuan Z, Zhou D, Johnson RL, Chen L (2019) Macrophage achieves self-protection against oxidative stress-induced ageing through the Mst-Nrf2 axis. Nat Commun 10 (1):755. doi:10.1038/s41467-019-08680-6
Lee S, Park YH, Chung JS, Yoo YD (2015) Romo1 and the NF-kappaB pathway are involved in oxidative stress-induced tumor cell invasion. Int J Oncol 46 (5):2021-2028. doi:10.3892/ijo.2015.2889
Srinivas US, Tan BWQ, Vellayappan BA, Jeyasekharan AD (2019) ROS and the DNA damage response in cancer. Redox Biol 25:101084. doi:10.1016/j.redox.2018.101084
Idelchik M, Begley U, Begley TJ, Melendez JA (2017) Mitochondrial ROS control of cancer. Semin Cancer Biol 47:57-66. doi:10.1016/j.semcancer.2017.04.005
Jin M, Wang J, Ji X, Cao H, Zhu J, Chen Y, Yang J, Zhao Z, Ren T, Xing J (2019) MCUR1 facilitates epithelial-mesenchymal transition and metastasis via the mitochondrial calcium dependent ROS/Nrf2/Notch pathway in hepatocellular carcinoma. J Exp Clin Cancer Res 38 (1):136. doi:10.1186/s13046-019-1135-x
Sosa V, Moline T, Somoza R, Paciucci R, Kondoh H, ME LL (2013) Oxidative stress and cancer: an overview. Ageing Res Rev 12 (1):376-390. doi:10.1016/j.arr.2012.10.004
Lee SB, Kim JJ, Chung JS, Lee MS, Lee KH, Kim BS, Do Yoo Y (2011) Romo1 is a negative-feedback regulator of Myc. J Cell Sci 124 (Pt 11):1911-1924. doi:10.1242/jcs.079996
Chung JS, Park S, Park SH, Park ER, Cha PH, Kim BY, Chung YM, Woo SR, Han CJ, Kim SB, Suh KS, Jang JJ, Lee K, Choi DW, Lee S, Lee GY, Hahm KB, Shin JA, Kim BS, Noh KH, Kim TW, Lee KH, Yoo YD (2012) Overexpression of Romo1 promotes production of reactive oxygen species and invasiveness of hepatic tumor cells. Gastroenterology 143 (4):1084-1094 e1087. doi:10.1053/j.gastro.2012.06.038
Kim HJ, Jo MJ, Kim BR, Kim JL, Jeong YA, Na YJ, Park SH, Lee SY, Lee DH, Lee HS, Kim BH, Lee SI, Min BW, Yoo YD, Oh SC (2017) Reactive oxygen species modulator-1 (Romo1) predicts unfavorable prognosis in colorectal cancer patients. PLoS One 12 (5):e0176834. doi:10.1371/journal.pone.0176834
Yu MO, Song NH, Park KJ, Park DH, Kim SH, Chae YS, Chung YG, Chi SG, Kang SH (2015) Romo1 is associated with ROS production and cellular growth in human gliomas. J Neurooncol 121 (1):73-81. doi:10.1007/s11060-014-1608-x
Rodriguez FJ, Folpe AL, Giannini C, Perry A (2012) Pathology of peripheral nerve sheath tumors: diagnostic overview and update on selected diagnostic problems. Acta Neuropathol 123 (3):295-319. doi:10.1007/s00401-012-0954-z
Boedtkjer E, Pedersen SF (2020) The Acidic Tumor Microenvironment as a Driver of Cancer. Annu Rev Physiol 82:103-126. doi:10.1146/annurev-physiol-021119-034627
Martinez-Outschoorn UE, Peiris-Pages M, Pestell RG, Sotgia F, Lisanti MP (2017) Cancer metabolism: a therapeutic perspective. Nat Rev Clin Oncol 14 (1):11-31. doi:10.1038/nrclinonc.2016.60

Download PDF

Journal Publication

published 07 Jan, 2022

Read the published version in Scientific Reports →

Version 1

posted

You are reading this latest preprint version

Network Models of Prostate Cancer Immune Microenvironments Identify ROMO1 as Heterogeneity and Prognostic Marker

Status:

Journal Publication

Version 1

Abstract

Figures

Introduction

Materials And Methods

Data acquisition

Determination of cluster count and membership

Identification of DEGs based on consensus clustering and PCA

Gene co-expression network analysis to obtain key gene modules associated with tumorigenesis and malignant progression

Validation and survival analysis of pivotal genes

Development and validation of the hub gene prognostic model

Relationship between key genes and tumor-infiltrating immune cells

Exploring the drug sensitivity of the hub gene

Result

Combined consistency clustering and principal component analysis to obtain sample cohorts

Identification of differentially expressed genes and GO functional annotation in PCa

Functional enrichment analysis of four groups of differentially expressed genes

Co-expression network analysis reveals key modules of PCa network

Validation of key candidate genes

Construction and validation of prognostic gene markers in TCGA

Alterations in the immune microenvironment may result from aberrant regulation of key genes in PCa

Drug sensitivity and resistance of key genes

Discussion

Conclusions

Abbreviations

Declarations

References

Supplementary Files

Status:

Journal Publication

Version 1