Stemness-related key genes are highly expressed in high Gleason grade prostate cancer samples
The results showed that the stemness of PCa cells of Gleason grade 3-5 was significantly higher than that of Gleason grade 1-2 (Fig. 1a), based on the OCLR's results and transcriptome data of PCa clinical specimen in TCGA. Besides, by using WGCNA (Weighted Gene Co-expression Network Analysis) to analyze the PCa transcriptome data for cancer cell stemness-related gene transcripts, we got 30 transcript modules related to PCa cell stemness (Fig. 1b), and many samples were clustered into different modules (Supplementary Fig. 1a). About 88.2% of the 2158 known cell stemness-related genes  were found in the results of our WGCNA analysis (Supplementary Materials S1; Supplementary Fig. 1c). Among these modules, we found that the MEmagenta module was most related to the stemness of prostate cancer cells (Fig. 1b), and the genes in the MEmagenta module were also positively correlated with stemness (Supplementary Fig. 1b). Therefore, the genes in the MEmagenta module (Supplementary Table S1) might play a vital role in the stemness of prostate cancer.
From the analysis results, we found that the expression of genes in the MEmagenta module in PCa cancer specimen was generally higher than that in the normal prostate specimen. Expression of these stemness genes in PCa increased with Gleason grade, and the highest expression was found at Gleason grade 5 (Fig. 1c). Further, we analyzed the data from 33 samples of the GEO database and found that most of the stemness genes expression in small cell prostate cancer (SCPC) were higher than those in prostate adenocarcinoma (PRAD) (Supplementary Fig. 1d). Subsequently, the transcriptome data of different prostate (normal or PRAD) cell lines and SCPC cell lines were analyzed, and it's found that the expression of stemness gene in PCa cells (PRAD cell lines) was higher than that in normal cells. Simultaneously, the expression of stemness genes in SCPC cell lines was slightly higher than in PRAD cell lines (Supplementary Fig. 1e).
In the genes of the MEmagenta module (Supplementary Table S1), we found that the known genes closely related to cancer stemness in non-prostate cancer, such as BRCA1 , EZH2 , FOXM1 , CDC20 , and CDCA8  were all clustered in this module, indicating these genes were also closely related to the stemness of prostate cancer cells. From supplementary table S1, we also found that EZH2 was the most significant gene in the correlation with the stemness of PCa cells.
The CNV of stemness-related genes elevated in prostate cancer
The segments containing all MEmagenta module genes were obtained from TCGA and then analyzed using the GISTIC (Genome Identification of Significant Targets in Cancer) method. Analysis results showed that most of the genes in the MEmagenta module were found to have high variations, which were distinctly higher than those in normal prostate tissue, in prostate cancer, indicating that the stemness of PCa cells was influenced by the CNVs (Copy Number Variations) of these genes (Fig. 2a and Supplementary Fig. 2a). In two types of CNV (amplification and deletion, Supplementary materials S2-S3), the change of deletion CNVs was more than amplification CNVs of the genes in the MEmagenta module of PCa samples. Both amplification and deletion CNVs in PCa samples were much higher than those in prostate normal samples (Fig. 2b and 2c, and Supplementary Fig. 2b and Supplementary Fig. 2c). Besides, by combining analysis with Gleason grade of PCa, GISTIC results showed that CNVs of genes in the MEmagenta module also increased with the Gleason grade (Fig. 2d), suggesting that CNVs of genes in the MEmagenta module were related to the malignancy of PCa.
In details of CNVs of some genes in the MEmagenta module, it was found that the CNVs (amplification) and expression of SKA3 and RUVBL1, which promoted tumor metastasis and played a role in the development of stem cells [37-40], were all increased with the Gleason grade (Supplementary Fig. 3a and 3b); and the CNVs (deletion) and expression of MCM6 and CENPH, which enhanced cancer cell proliferation, stemness, and metastasis and promoted cancer development [41-45], were also increased with the increase of Gleason grade (Supplementary Fig. 3c and Supplementary Fig. 3d).
Transcription factors of stemness-related genes enhanced the stemness of PCa cells
ATAC-seq data of the MEmagenta module genes in PCa samples obtained from TCGA and the sequence data of the 2kb range of the transcription start site (TSS) were analyzed and displayed. Results showed that both normal and tumor prostate samples had the binding signals of transcription factors. The binding signals of all PCa samples (including the different Gleason grades of PCa and small cell PCa) were weaker than that of normal prostate samples (Fig. 3a), and the binding signals of transcription regulators in small cell PCa was the weakest ones (Fig. 3a). These results indicated that the transcriptional regulation manner and the types of transcription regulators of the MEmagenta module genes in PCa might be different from those in normal prostate tissues.
According to the ATAC-seq data analysis results, it was found that most transcription factors of stemness genes in prostate samples were the same (Supplementary Fig. 4b). After filtering the 238 common transcription factors according to the genome information and transcriptome data, the expression profile of 163 transcription factors was finally obtained. After clustering the transcriptome data of the common 163 transcription factors, 3 clusters were obtained according to the expression of 163 common transcription factors (Supplementary Fig. 4a). The first cluster has 86 transcription factors (such as FOXA1, HOXB13, ERG1, and so on), the second cluster has 44 transcription factors (such as SOX15, FOXM1, SOX2, and so on), and the third cluster has 33 transcription factors (such as NANOG, SOX3, FOXA2 and so on) (Supplementary Fig. 4a). In the top ten results of principal component analysis (PCA) of 163 common transcription factors, it was found that stemness-related PUM1 [46, 47], CLOCK , SP1, and TCF12 played a major positive regulation role on PCa cell stemness (Fig. 3b and Supplementary Fig. 4c). In comparison, IRF3  negatively correlated with other 9 transcription factors and played a negative regulatory role on PCa cell stemness (Fig. 3b). As we know, IRF3 (interferon regulatory factor 3) signaling played an essential role in TLR3-mediated apoptosis in LNCaP cells through the activation of the intrinsic and extrinsic apoptotic pathways , suggesting that the immune system might play a role in suppressing the stemness of PCa cells. Furthermore, our analysis results showed that AR, FOXA1, NFYA, CTCF, and FOXO3 might enhance the stemness of PCa cells, where FOXF1 might be negatively correlated with these transcription factors (Fig. 3c).
In normal prostate samples, we found that the major transcription regulators were HOXB4, NFYB, and TFE3 (Supplementary Fig. 5a and Supplementary Fig. 5c), which differed from those in prostate cancer. In normal samples, the role of FOXA1, NFYA, and FOXP1 in regulating stemness genes was changed by comparing their results of PCA analysis in prostate cancer (Supplementary Fig. 5b and Fig. 3c). These results indicated that the same transcription factors might play different roles in regulating the cell stemness in prostate normal and cancer samples.
Our transcriptional regulator analysis results showed that the upstream of the EZH2, the most relevant gene to the stemness of PCa cells, could be significantly bound by NFY (Fig. 3d). As a transcription factor, NFY regulates the self-renewal of hematopoietic stem cells  and promotes the self-renewal and expansion of prostate cancer cells and their stemness . Hence EZH2 might play a stemness role in prostate cancer.
Immunological microenvironment negatively related to the stemness of PCa cells
By scoring the stemness and immunity of the immune microenvironment of PCa in different Gleason grades of clinical samples, there was a negative correlation between the PCa cell stemness and the immunity of the microenvironment of PCa clinical samples in all Gleason grade. The correlation coefficient of Gleason grade 1-2 was almost the same as that of Gleason grade 3-5 (Fig. 4a). By analyzing and scoring the stromal and immune cells of the PCa microenvironment, we found that the scores of stromal and immune cells in the PCa microenvironment were all inversely related to the stemness of PCa cells (Supplementary Fig. 6a). Although the immune score increased with the increase of Gleason grade, in each type of Gleason grade, the stemness of the low-immunity group was significantly increased than the stemness of the high-immunity group. (Supplementary Fig. 6b). After clustering the MEmangenta module genes based on cell stemness and immune scores of PCa, we found that PCa cell stemness was negatively related to the immunity of the PCa microenvironment in the clinical samples with high expression of genes of the MEmagenta module (Supplementary Fig. 6d).
Furthermore, by analyzing and scoring the immune infiltration of PCa clinical samples from TCGA, the number of 22 types of immune cells in prostate cancer with different Gleason grades was obtained (Fig. 4b). The number of most types of immune cells in the microenvironment of PCa was increased with the Gleason grade increase; and CD8+ T cells  and macrophage M1  were the most significantly increased in all types of immune cells, while the plasma cells  were reduced considerably with the Gleason grade increase (Fig. 4b), indicating that plasma cells (B cells) in the microenvironment of PCa played a vital role in anti-PCa immunity.
According to the results of immune cell correlation, immune cells were divided into 4 clusters by clustering. Fifteen immune cells were clustered in cluster-A, three immune cells were clustered in cluster-B, two immune cells were clustered in cluster-C, and two immune cells were clustered in cluster-D. And there were twelve edges in cluster-A, and the edges in other clusters were all connected to cluster-A. Therefore, cluster-A was more complicated in this immune cell network. For examples, the activated NK cell was negatively correlated with resting NK cell, which was most significant in all correlation among immune cells (Supplementary Fig. 6c), and the number of activated NK cell increased with the Gleason grade increase in PCa (Fig. 4b), indicating that activated NK cell might play a significant inhibitory function on prostate cancer stemness; the most significant positive correlation was between the activated dendritic cell and the memory B cell, and the number of both activated dendritic cells and memory B cells was all increased with the Gleason grade increase (Fig. 4c), which indicated that these two types of cells might play essential roles in inhibiting the stemness of PCa cells.
The protein interaction network of stemness-related genes and its relationship with prostate cancer immune cells
By screening and analyzing human protein interaction data containing MEmagenta module genes from the STRING database, we found the important protein-protein interaction network in the proteins of MEmagenta module genes and found that EZH2 interacted directly with 17 proteins. In the EZH2-related 17 protein-protein interactions, we also found that EZH2 could regulate the entire protein interaction network of MEmagenta module stemness genes by mainly interacting with CENPA, BUB1B, and PARP1 (Fig. 5a) [52-55]. Furthermore, function enrichment analysis of genes of the MEmagenta module revealed that most function of these genes was concentrated in cell mitosis, and the most significant functional pathway was related cell cycle (Supplementary Fig. 7a and Supplementary Fig. 7b), suggesting that these stemness genes might involve in the regulation of PCa stem cell mitosis.
From the analysis results of the relationship between the expression of the MEmagenta module gene and immune infiltration of PCa, we found that different immune cells had different effects on the expression of stemness genes in PCa. Different types of immune cells could affect the same stemness gene expression, and one stemness gene expression could also reversely affect different types of immune cells (Supplementary materials S5). In correlation of expression of stemness genes and immune cells, we found that expression of most stemness genes was positively correlated to memory B cells and naive B cells and negatively correlated to plasma cells (Fig. 5b). Further, the number of these B cells increased, while the number of plasma cells decreased, increasing the PCa Gleason grade (Fig. 4b). These results indicated that B cells might play the opposite effects on PCa cell stemness in different conditions (it's consistent with reference ). Besides, the number of activated NK cells and memory CD4+ T cells were all increased with the Gleason grade; and the expression of stemness genes of the MEmagenta module was positively correlated to resting NK cells and memory CD4+ T cells, while negatively correlated to activated NK cells (Fig. 4b and 5b).
Expression of EZH2, the most relevant gene to PCa cell stemness, was most positively correlated to the activated memory CD4+ T cells and negatively correlated to the resting Mast cells (Fig. 5c and 5d). Most types of cells positively correlated to the expression of the EZH2 were T cells and B cells (Fig. 5c), suggesting that T cells and B cells were the critical immune cells in regulating PCa cell stemness by controlling the expression of the EZH2. Further, the number of resting mast cells decreased, and the number of activated mast cells increased with the increased PCa Gleason grade (Fig. 4b); and the resting mast cells were the most negatively correlated to the expression of the EZH2 (Fig. 5c and 5d). These results indicated that the immunity and immune cells of the microenvironment of PCa played an essential role in the tumorigenesis and development of prostate cancer.