Shan-Qi Guo,Hao Deng,Ying-Jie Jia,Xiao-Jiang Li＊,Xing-Kang Jiang

1Department of Oncology,First Teaching Hospital of Tianjin University of Traditional Chinese Medicine,Tianjin,300381,China;

2Department of Urology,The Second Hospital of Tianjin Medical University,Tianjin,300211,China.

Abstract Background:The Genotype-Tissue Expression was used to expanded normal tissue of the Cancer Genome Atlas database.This study aimed to investigate genes associated with the pathogenesis and prognosis of prostate cancer.Methods:We conducted prognostic related genes for prostate cancer by using transcriptome data from the Genotype-Tissue Expression Project and the Cancer Genome Atlas data sources,which were analyzed using an integrated bioinformatics strategy.Clinically significant modules were distinguished,and GO and KEGG analysis were used to Database for Annotation,Visualization and Integrated Discovery.Further annotation was performed through Gene set enrichment analysis.Logistic regression was carried out to analyze the associations between clinicopathologic characteristics and the hub genes.Logistic regression model and survival analysis were performed. Results:By using data available from the Cancer Genome Atlas and the Genotype-Tissue Expression databases,we here show that 53 differential expression genes were identified.Through GO and KEGG analysis a prognostic related gene signature consisted of GOLM1,EIF4A1,ABCC4,RPL7P16,NPIPB12 and PCA3 was constructed with a good performance in predicting overall survivals.The majority of the six hub genes were associated with clinical characteristics of prostate cancer. Conclusion:These genes might be considered as new targets for further investigating the diagnostic and prognostic biomarkers to facilitate the molecular targeting therapy since they showed differently expressed in prostate cancer and correlate with overall survival prognosis.

Key words:Prostate cancer,prognosis,bioinformatic analysis,the Cancer Genome Atlas;Genotype-Tissue Expression

Background

Prostate cancer (PCa) is the second most common male malignancy tumor in the world,and the occurrence rate and the mortality of PCa keeps increasing [1,2].It’s worth noting that there still quite a few PCa patients suffering local recurrence or distant metastases and die of their disease.At present,Gleason grade,clinical stage and serum PSA level are mostly used items for managing the risk degrees of PCa in clinical [3].However,clinical and scientific researches indicated that potential mechanisms and molecular pathways contributing to tumorigenesis and tumor progression is essential to enhance diagnostic,prognostic,therapeutic approaches of PCa.The exploration of accurate biomarker to cognit invasive tumors,and to predict the clinical outcome of PCa is of great significance[4,5].

Nowadays,global shared databases really do us a favor to have a deeper understanding on the prognostic factors of PCa.Noticeably,most bioinformation studies on PCa prognosis are based on the Cancer Genome Atlas (TCGA) database or Gene Expression Omnibus.Based on TCGA data,Chen Xuan et al [6]concluded CCNB1,TTK,CNN1,and ACTG2 were correlated with prognosis of PCa,and the protein levels of CCNB1,TTK,and ACTG2 had a degree of differences between tumor tissues and normal tissues.As we known,updated TCGA contains 499 PCa samples while only 52 normal samples of prostate.Based on inadequate control sample,Mu Haiqi and colleagues [7] chose GSE16560 and TCGA data showed a prognostic prompt function of the nine-gene signature,including GART,RRM2,CACNA1D,PIK3R1,MYH11,MYLK,TNC,MYL9;and in their study,KCNMA1 was an independent predictor of prognosis in PCa;while there were some limitations in the study,for instant,the microarray dataset GSE55945 only including 13 malignant and 8 normal fresh-frozen prostate tissue specimens,which was used to identify the differential expression genes (DEGs) between the two groups.He Z et al [8] verified that ABCC4 and SLPI might be utilized as prognostic indicators of PCa filtered out from the GSE103512 dataset,which included 60 PCa samples,and 7 matched normal samples.As you can see from the studies above,lack of normal tissue samples may be an obstacle to the accuracy of prognostic genes.

As we know,TCGA includes non-tumor samples from the same cancer participants,and normal samples in TCGA is often limited or lacking.Considering a limitation of some of the TCGA datasets is paucity of corresponding normal samples,we integrated into TCGA data to use normal samples from the Genotype-Tissue Expression (GTEx) project,which is another large-scale repository cataloging gene expression from healthy individuals,to expanded the information of the normal samples,in order to obtaining more accurate different genes related to prognosis.The GTEx project has generated a systematic,multi-tissue reference for identifying genetic variants associated with changes gene expression[9,10].What attracted to us was that GTEx had normal samples of related organs,which could complement studies and expanding molecular data from ‘normal’ individuals.In light of the above situation,we combined the prostate tissues either with carcinoma or normal samples from TCGA and GTEx,to identify the prognostic signature of hub genes.

Materials and methods

Data collection

TCGA PCa dataset contained normal prostate samples(n=52) and PCa samples (n=499).The transcriptional data and clinical data were downloaded from TCGA)(http://tcga-data.nci.nih.gov/tcga/database.Moreover,100 normal prostate samples were obtained from the GTEx Project which was downloaded from the UCSC Genome Browser(https://genome.ucsc.edu).

Identification of DEGs

We applied the expression profiles to the edgeR package in R language and calculated for differential expression between tumor and normal group samples after normalization and filtration.The adjusted P values (adj.P) were applied to correct the false-positive results.Then the significant DEGs (adj.P<0.05 and fold-change value larger than 3) were selected out for the next step analysis.

Survival analysis

Overall survival (OS) Kaplan–Meier estimate was performed in this study.We performed Kaplan–Meier analysis (R package “survival”) to assess survival and relapse difference across cases with DEGs.For each gene,patients were divided into high-expression group and low-expression group according to the median expression level.Based on the survival curves of each group,these upregulated and downregulated DEGs with a poorer survival rate were within our consideration.Pvalue <0.05 was set as the cutoff point.

Functional and pathway enrichment analysis of DEGs

The Database for Annotation,Visualization and Integrated Discovery (DAVID,https://david.ncifcrf.gov/) is a biological database regularly used to facilitate functional annotation and pathway analysis.KEGG pathway analysis aims to identify and visualize significantly enriched pathways of molecular interactions,reactions,and relations.GOanalysis uses hypergeometric tests to perform enrichment analysis on gene sets.We uploaded selected DEGs associated to DAVID to perform KEGG pathway and GO enrichment analysis.

Gene set enrichment analysis(GSEA)

GSEA,a computational method,determines whether a predefined set of genes shows significant differences between two biological phenotypes under study.Gene set permutations were performed 1,000 times for each analysis.The nominalPvalue and normalized enrichment score were used to sort the pathways enriched in each phenotype.

Results

Data processing and screening of DEGs

The PCa microarray expression dataset of TCGA contained 499 PC samples and 52 normal prostate samples.GTEx Project contains 100 normal prostate tissues.A total of 55187 gene symbols were discerned and the gene expression matrix was constructed.The DEGs were filtered using the limma R package(criteria:adjustedP<0.05 and|log 2 fold change|>3).A total of 53 DEGs were identified between PCa and NP samples,including 15 upregulated genes and 38 downregulated genes.The cluster heatmap of the 53 DEGs were collected in Figure 1.

Figure 1 Cluster heatmap of the differentially expressed genes.

Functional and pathway enrichment analyses of DEGs

GO function annotation and KEGG pathway analyses of the identified DEGs was obtained using the DAVID database and its online analysis tool (Table 1).GO function annotation included the following three portions:biological process(BP),cell component(CC),and molecular function (MF).The results were statistically significant ifP<0.05.As shown in Figure 2A,the upregulated genes were mainly enriched in translation factor activity,RNA binding (ontology:MF).The relationship between GO and relative genes were showed in Figure 2B.The ensuing KEGG pathway analyses results were collected and shown in Figure 2C.The upregulated DEGs enriched in Ribosome biogenesis in eukaryotes pathway and downregulated DEGs enriched in Ribosome and RNA transport pathway.The heatmap of the relationship between DEGs and pathways was also displayed.The Z-score of enriched pathways less than zero indicated that most of the cancer pathways were more likely to be decreased.

Figure 2 GO and KEGG analysis for the differential expression genes.

Table 1 Enriched gene ontology terms and pathways of differentially expressed genes

Verification core prognostic genes

To reveal the relationship between DEGs and patient survival,we draw Kaplan-Meier curves by conducting the “survival” package in R language according to the median expression level of DEGs among 651 patients.These upregulated and downregulated genes with a significantly poorer survival (P<0.05),respectively,were determined.Kaplan-Meier curves of the partial genes between the high-expression group and low-expression group were shown in Figure 3.Golgi membrane protein 1 (GOLM1),eukaryotic translation initiation factor 4A1 (EIF4A1),ATP binding cassette subfamily C member 4 (ABCC4),ribosomal protein L7 pseudogene 16 (RPL7P16),nuclear pore complex interacting protein family member B12(NPIPB12)and PCa associated 3 (PCA3) were considered to be prognosis related genes with thePvalue <0.05 above 53 candidate genes.The results showed that GOLM1 had the best prognostic value and the most significantPvalue (P=0.003).Based on the above analysis,the six genes were predicted to be the core factors differentially expressed between PCa and normal prostate tissues.Meanwhile,the co-relationship of six core prognostic genes in PCa and normal prostate samples were concluded in Figure 4.

Figure 3 Kaplan-Meier survival curves of the prostate cancer patients with high and low expression of the six core prognostic genes.

Figure 4 The gene expression levels of the six genes were detected in prostate cancer samples and normal group.

Clinical relevance and GSEA

We could clearly identify the core factors as cancer risk factors,the expression level increased as Gleason score elevated and higher T,N stage.Logistic regression was carried out to analyze the associations between clinicopathologic characteristics and the core prognostic genes (Table 2).Higher expression of GOLM1,ABCC4 and PCA3 predicts higher T stage.Five(GOLM1,EIF4A1,ABCC4,RPL7P16 and PCA3)of six core gene hinted with N1 stage compared with N0 stage.With regard to Gleason score,the expression of GOLM1,EIF4A1,ABCC4,RPL7P16 and PCA3 were associated with higher scores of more than 7.Moreover,univariate and multivariate analysis were executed to determine significant prognostic genes that related to poor outcome,we found GOLM1 and NPIPB12 acted as independent prognostic factors in TCGA PCa (Figure 5).GSEA analysis suggested that amino acid metabolism pathway and Perocxisome signaling pathway were enriched in GOLM1 and RPL7P16 high expression group jointly,other signaling pathway enriched with different genes were also showed in Figure 6.

Figure 5 Hazard ratio of relative clinicopathologic characteristics and the expression of hub genes.

Figure 6 Gene set enrichment analysis pathway enrichment analysis.Pathway with FDR <0.05 were considered to be significance.

Table 2 Clinicopathologic characteristics of core prognostic genes

Discussion

PCa is one of the most common cancers in men worldwide [1].Patients were considered well treatable if diagnosed at early stage.Appropriate biomarkers especially those related to prognosis were of great essential to make clinical decisions [11].The aim of our study was to delineate the prognostic biomarker of PCa.For this purpose,we not only used TCGA datasets,but also GTEx,which could expand the information of the normal samples.Thus,we finally obtained more accurate different genes related to prognosis of PCa.

Results showed that PCa in TCGA cohort represented 499 tumor and 52 normal tissue samples,while GTEx Project contains 100 normal prostate tissues.Thus,we finally got 499 tumor samples and 152 normal tissue sample.Thus,using GTEx data and TCGA database could perfectly integrate and enhance the diversity and accuracy of data resources.At the same time,GTEx well complemented the lack of normal samples of TCGA database,making our results more credible.Results indicated that a total of 53 DEGs were identified between PCa and NP samples,including 15 upregulated genes and 38 downregulated genes with the combination of the above two databases.

Meaningfully,as for the prognostic role of these genes during the statistical stages of Kaplan-Meier,six genes between the high-expression group and low-expression group including GOLM1,EIF4A1,ABCC4,RPL7P16,NPIPB12 and PCA3 were considered to be prognosis related with thePvalue <0.05 above 53 DEGs.Furthermore,the expression of six genes were detected of significant difference in PCa samples and normal group.As for the prognostic role of these genes,we hereby provide a comprehensive analysis;clinical relevance revealed higher expression of GOLM1,ABCC4 and PCA3 related to higher T stage.GOLM1,EIF4A1,ABCC4,RPL7P16 and PCA3 hinted with N1 stage compared with N0 stage.Moreover,the expression of GOLM1,EIF4A1,ABCC4,RPL7P16 and PCA3 were associated with higher scores of more than 7.Moreover,GOLM1 and NPIPB12 acted as independent prognostic factors in TCGA PCa predicting poor outcome.Chen et al investigated Golgi Membrane Protein 1(GOLM1)in PCa,and they found that GOLM1 was considered as tumor promoter and miR-3935 was identified as the upstream regulator of GOLM1,which were similar to our results that GOLM1 was upregulated in PCa[12].

Zhang and colleagues used Gene Expression Omnibus dataset which included 60 PCa samples and 7 normal prostate samples concluded that GOLM1 and ACLY in PCa may be the potential diagnostic and prognostic biomarkers [13].Yan Guang et al used tissue microarrays and real-time RT-PCR,Western blot,and immunohistochemistry analyses to investigate GOLM1 functions in vitro and in vivo,they substantiated that GOLM1 played oncogenic functions mainly through activating PI3K-AKT-mTOR signaling pathway [14].This remind us that agents block PI3K-AKT-mTOR signaling pathway might be used in PCa patients with GOLM1 up-regulation.Our results showed that EIF4A1 was closely related to the prognosis of PCa.Li Zengxia et al provide evidence that USP9X was a novel regulator of the translation initiation process via deubiquitination of EIF4A1,which offers new insight in understanding the pivotal role of USP9X in human malignancies and neurodevelopmental disorders,thus interpreted the possible mechanism of EIF4A1in PCa [15].ABCC4 is positively correlated with androgen levels in PCa patients.Addition of androgens in PCa cells significantly upregulated ABCC4 and the addition of anti-androgens decreased expression [16].These verified that ABCC4 gene and its association with androgen levels as potential therapeutic targets in PCa[17].As we known,ATP-binding cassette (ABC)transporters are transmembrane proteins responsible for the efflux of a wide variety of substrates.For better understand of the role of ABC transporters in PCa,Demidenko Ritaet al.analyzed the profile of ABC transporter gene expression in PCa and noncancerous prostate tissues.The results showed that eight ABC transporter genes (ABCA8,ABCB1,ABCC6,ABCC9,ABCC10,ABCD2,ABCG2,and ABCG4) were down regulated in PCa compared with normal samples,and only two genes (ABCC4 and ABCG1) were up-regulated.In the current study,we also found that high expression of ABCC4 was correlated with poor outcome [18].Gao et al studied on the side effects of antiandrogen flutamide,which was used to treat PCa,and they found that aryl hydrocarbon receptor may mediate the regulation of the ABCC4 by flutamide,potentially leading to FXR signaling inhibition and resulting in cholestatic liver injury [19].Indicating the potential pathways that ABCC4 involved in PCa treatment.Many researches have declared that PCA3 played a crucial role in PCa [20].Rodon Natalia et al clarified that PCA3 score may serve as an additional marker to reduce the indication for biopsies and avoid overdiagnosis and overtreatment of patients with suspected PCa.Meanwhile,the prognostic significance of s-PCA3 was confirmed and results implicated that PCA3 was closely related to the prognosis of PCa[21].

As for get knowledge of relative pathways of hub genes,GSEA for differentially expressed genes allowed us to detect many biological processes or pathways altered in PCa.We performed a detailed analysis of the six genes showing that amino acid metabolism pathway and Perocxisome signaling pathway were enriched in GOLM1 and RPL7P16 high expression group jointly.Additionally,drug metabolism cytochrome P450 and arachidonic acid metabolism pathway were low expressed in multiple GSEA analysis.As we all know,the above pathways were closely bound up with PCa,whether it refers to the mechanism or the treatment of PCa.This gives us a broad perspective of how much relative signaling pathways we could invest to cover all potentially relevant targets,which have not been in the scope of PCa researchers yet.

In order to identification of appropriate biomarkers,we represent two very useful reference databases to validate the results of different levels of gene expression and provided possible pathways for further study.We believe that our study should also be a starting point for important researches and had particular advantage of prognostic biomarkers with the consideration of normal tissue samples.

Conclusion

Our study provides a comprehensive overview of the prognostic biomarkers in patients with PCa,GOLM1,EIF4A1,ABCC4,RPL7P16,NPIPB12 and PCA3 was constructed with a good performance in predicting overall survivals.And the majority of the six hub genes were associated with clinical characteristics of PCa.These findings could provide a reference point for further biomarker research and validation data for ongoing projects by most PCa researchers.