
==== Front
Comput Struct Biotechnol J
Comput Struct Biotechnol J
Computational and Structural Biotechnology Journal
2001-0370
Research Network of Computational and Structural Biotechnology

S2001-0370(24)00285-X
10.1016/j.csbj.2024.08.026
Research Article
UQCRB and LBH are correlated with Gleason score progression in prostate cancer: Spatial transcriptomics and experimental validation☆
Quan Yongjun qyongjun@yeah.net
a
Zhang Hong b
Wang Mingdong a
Ping Hao pinghaotrh@ccmu.edu.cn
a⁎
a Department of Urology, Beijing Tongren Hospital, Capital Medical University, Beijing 100176, China
b Department of Pathology, Beijing Tongren Hospital, Capital Medical University, Beijing 100176, China
⁎ Corresponding author. pinghaotrh@ccmu.edu.cn
30 8 2024
12 2024
30 8 2024
23 33153326
8 6 2024
9 8 2024
29 8 2024
© 2024 The Authors
2024
https://creativecommons.org/licenses/by-nc-nd/4.0/ This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
Prostate cancer (PCa) is a multifocal disease characterized by genomic and phenotypic heterogeneity within a single gland. In this study, Visium spatial transcriptomics (ST) analysis was applied to PCa tissues with different histological structures to infer the molecular events involved in Gleason score (GS) progression. The spots in tissue sections were classified into various groups using Principal Component Analysis (PCA) and Louvain clustering analysis based on transcriptome data. Anotation of the spots according to GS revealed notable similarities between transcriptomic profiles and histologically identifiable structures. The accuracy of macroscopic GS determination was bioinformatically verified through malignancy-related feature analysis, specifically inferred copy number variation (inferCNV), as well as developmental trajectory analyses, such as diffusion pseudotime (DPT) and partition-based graph abstraction (PAGA). Genes related to GS progression were identified from the differentially expressed genes (DEGs) through pairwise comparisons of groups along a GS gradient. The proteins encoded by the representative oncogenes UQCRB and LBH were found to be highly expressed in advanced-stage PCa tissues. Knockdown of their mRNAs significantly suppressed PCa cell proliferation and invasion. These findings were validated using The Cancer Genome Atlas Prostate Adenocarcinoma (TCGA-PRAD) dataset, as well as through histological and cytological experiments. The results presented here establish a foundation for ST-based evaluation of GS progression and provide valuable insights into the GS progression-related genes UQCRB and LBH.

Graphical Abstract

ga1

Keywords

Prostate cancer (PCa)
Visium spatial transcriptomics (Visium ST)
Gleason score (GS)
Glandular epithelial (GE) cells
UQCRB
LBH
==== Body
pmc1 Introduction

Prostate cancer (PCa) is the most common malignancy and the leading cause of tumor-related mortality among men more than 50 years of age, with an estimated 299,010 newly diagnosed cases and 35,250 associated deaths in the United States in 2024 [1]. Therefore, PCa represents a serious health concern, and additional research is needed to clarify its underlying mechanism and make it possible to design definitive therapies.

PCa is a multifocal disease that is characterized by genomic and phenotypic heterogeneity within a single gland [2], [3]. Recent studies indicate that the substantial intertumoral heterogeneity of PCa arises through the evolution of a monoclonal precursor [4], [5], [6], [7], [8] and that each lesion present in multifocal PCa represents a different developmental stage in the evolutionary process. Compared with that in primary PCa patients, a substantially greater percentage (19.3 % overall) of metastatic castration-resistant PCa (mCRPC) patients exhibit aberrations in homologous recombination repair-related genes such as BRCA1, BRCA2, and ATM [9], and aberrations of this type are linked to therapeutic sensitivity to PARP inhibitors (olaparib) [10], [11]. Therefore, exploration of the genetic differences among multiple PCa lesions can facilitate the identification of the origins of PCa and the triggers of PCa progression and thus may have clinical relevance.

Glandular epithelial (GE) cells, such as luminal cells and basal cells, are reported to be the cells of origin of PCa [12], [13], [14], [15], [16], [17], [18]. PCa cells are diffusely distributed in prostate tissue and form multiple foci, and certain regions of tissue include both tumor tissue and nonglandular [mainly anterior fibromuscular stroma (AFMS)] tissue. Because these foci lack clear visible boundaries with nontumor cells, accurately identifying the gene expression patterns of GE cells is the main challenge in the transcriptomic analysis of PCa.

Single-cell analysis [19], [20], [21], [22], [23], [24] and spatial transcriptomics (ST) [8], [25], [26], [27] have been used to detect drivers of PCa transdifferentiation to malignant stages. Although single-cell analysis provides transcriptome information at the cellular level, this approach does not provide information on spatial patterns and depends to a certain extent on inferences made via bioinformatics [28], [29]. ST combines high-throughput RNA sequencing data with spatial dimensions [8], [25], [26], [27], [30] and makes it possible to quantify total mRNA in the spatial context of histological structures. Visium ST analysis offers high spatial resolution, employing spatial microarrays with 4992 spatially barcoded spots, each with a diameter of 55 μm within a capture area measuring 6.5 × 6.5 mm2 [8], [27]. This high-resolution platform is highly effective for distinguishing between different histological structures by enabling precise spatial mapping of gene expression across tissue sections.

A fundamental prerequisite for analyzing the transcriptome-wide expression patterns that occur in multifocal PCa is defining criteria for the malignant staging of individual lesions using histological images. The Gleason grading system is extensively used in the pathological diagnosis of PCa, and it is the most reliable architectural grading system for evaluating disease severity; the assigned scores range from 1 to 5, representing well-differentiated to poorly differentiated GE cells [31], [32], [33], [34], [35]. The use of the Gleason grading system enables standardized risk assessment of PCa patients with different prognostic outcomes, exemplified by biochemical recurrence and overall survival [31], [32], [34], [35], [36], [37]. The Gleason score (GS) is based on primary and secondary Gleason patterns and ranges from 1 + 1 to 5 + 5; PCa patients with GSs with sums ≥ 7 are considered to have aggressive disease with a high risk of extraprostatic extension and a poor prognosis [38], [39].

In this study, Visium ST analysis was performed on PCa tissues, the transcriptome-wide expression patterns of different histological structures of GE cells were comprehensively analyzed, and the molecular events involved in GS progression were inferred. We identified key genes that are positively associated with PCa malignancy, specifically encode ubiquinol–cytochrome c reductase-binding protein (UQCRB) and limb-bud and heart (LBH). These associations were validated using The Cancer Genome Atlas Prostate Adenocarcinoma (TCGA-PRAD) dataset, as well as through histological and cytological experiments.

2 Materials and methods

2.1 Collection of prostate cancer (PCa) tissues

Forty-five PCa patients who underwent prostatectomy were recruited; their clinicopathological information is shown in Table S1. PCa and prostatic intraepithelial neoplasia (PIN) tissues from these patients were embedded in paraffin for histological analysis. Two PCa tissues with GSs of 4 + 4 and 3 + 4 (obtained from Patient 2 and Patient 1 in our previous study [27]) were subjected to Visium spatial transcriptome sequencing (Table S2). The extraction sites for analysis were determined by prostate magnetic resonance imaging (MRI) and pathological diagnosis of prostate needle core biopsy samples obtained before radical surgery. The fresh tissues were snap-frozen on dry ice, embedded in optimal cutting temperature (OCT) compound, and stored at −80 °C until use.

2.2 Visium spatial transcriptomics (ST) (10x Genomics Visium)

Visium ST analysis was performed as previously described [8], [27]. OCT-embedded tissue sections with RNA integrity numbers (RINs) ≥ 7 were cryosectioned at a thickness of 10 µm and placed on Visium Spatial slide. Spatially barcoded cDNA derived from tissue sections was generated using the Visium Spatial Gene Expression Slide & Reagent Kit (10 × Genomics, CA, USA). The tissue sections were prepared for methanol fixation and hematoxylin—eosin (H&E) staining and subsequently subjected to optimized permeabilization conditions (permeation time of 18 min) using the Visium Spatial Tissue Optimization Slide & Reagent Kit (10 × Genomics). A library was constructed on the basis of the 10x Genomics Visium library preparation protocol and sequenced on the Illumina NovaSeq600 platform (Illumina, CA, USA). Splicing-aware alignment of the reads was performed by matching to a reference genome using the STAR aligner of Space Ranger (10 × Genomics).

2.3 Analysis of gene expression using OmniAnalyzer Pro

The transcriptome profiling data output from Space Ranger were analyzed via OmniAnalyzer Pro [Analytical BioSciences Limited (Abiosciences), Beijing, China] as previously reported [27] to determine spatial gene expression patterns. OmniAnalyzer Pro is based on Scanpy and on the algorithms developed by Abiosciences. Following quality control measures such as spot filtering (with minimum count among all spots: 15.0, maximum count among all spots: 58,208.0, minimum number of detected genes per spot: 14, maximum percentage of mitochondrial gene expression: 20.00 %), gene filtering (with minimum total expression count of all genes: 0.0, maximum total expression count of all genes: 560,533.0, minimum number of spots in which any gene is expressed: 0, maximum number of spots in which any gene is expressed: 4,271), and identification of highly variable genes (HVGs) (with minimum mean value: 0.0125, maximum mean value: 3, minimum dispersion: 0.5, and maximum dispersion: 1,000), the read count data were normalized to ensure uniform total expression across all spots. The data were subsequently transformed using the natural logarithm [ln (1 + normalized value)] to increase the normality of the distribution and mitigate the impact of highly expressed genes. The image processing algorithm was applied to the capture area of the PCa section according to the spatial barcode information in the clean data and the reads for each spot.

2.3.1 Dimensionality reduction

To reduce the dimensionality of the various transcriptome features of PCa sections, principal component analysis (PCA) (https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html) was performed using HVGs with a minimum mean value of 0.0125, a maximum mean value of 3, a minimum dispersion of 0.5, and a maximum dispersion of 1,000. The pipeline uses a python implementation of scipy.sparse.linalg ARPACK for singular value decomposition (SVD). Fifty components were identified from the PCA of the transcriptome profiling data. To visualize the data in a 2-dimensional space, the PCA-reduced data were subsequently passed to uniform manifold approximation and projection (UMAP). The spots were then divided into groups through Louvain clustering analysis. On the basis of this classification, the spots overlaid on the histological image were assigned different colors so that transcriptomic patterns could be simultaneously aligned with histological structures.

2.3.2 Inferred copy number variation (inferCNV) analysis

The acquisition of some genomic aberrations is correlated with PCa progression to advanced stages [4], [5], [6], [7], [8], [9], [40]. CNV analysis facilitates the prediction of clonal hierarchies and tumor malignancy [4], [5], [6], [7], [8], [20], [21]. In this study, somatic large-scale chromosomal CNVs were inferred via Python reimplementation of inferCNV using transcriptomic data from the ST analysis. The steps used for the inferCNV analysis were described previously (https://infercnvpy.readthedocs.io/en/latest/infercnv.html). The intensity of expression of genes across positions in the whole genome was explored in comparison to a set of “normal” reference spots. Because of the absence of definite normal prostate gland cells in the PCa sections, we used cells with a GS of 3 + 3, which have been proposed to represent precancerous lesions [41], [42], as a reference group because of their well-differentiated nature and low capacity for metastasis [43], [44], [45]. The CNV values of local genome fragments across 22 autosomes were visualized in a heatmap, and the median CNV signals in each cluster were illustrated in boxplots to determine the cell of origin [4], [5], [6], [7], [8] and malignancy-related features of each cluster [20], [21].

2.3.3 Diffusion pseudotime (DPT) analysis

DPT enables the reconstruction of the developmental progression of PCa cells in a biological process from snapshot data by illustrating geodesic distance along the graph [46], making it possible to identify branching decisions and differentiation endpoints [46]. It is implemented within Scanpy [47]; the computational algorithm of DPT was described previously [46] (https://scanpy.readthedocs.io/en/stable/generated/scanpy.tl.dpt.html).

2.3.4 Partition-based graph abstraction (PAGA) analysis

PAGA is an extended version of DPT that reconciles clustering with trajectory inference through a topology-preserving map of single cells [48]. In ST analysis, an interpretable abstracted graph-like map (PAGA graph) of continuously connected and discrete disconnected spot-to-spot variations is generated. It consistently predicts developmental trajectories and gene expression changes using a transcriptomic dataset [48]. The relative algorithms used to obtain the PAGA were described in a previous study [48] (https://scanpy.readthedocs.io/en/stable/generated/scanpy.tl.paga.html).

2.3.5 Analysis of differentially expressed genes (DEGs)

DEG analysis in ST facilitates exploration of the heterogeneity of gene expression during the development of histological clusters. The DEGs in the comparisons between the specified clusters were identified through the Wilcoxon test, in which an absolute value of natural logarithm (ln) fold change (FC) (|ln FC|) > 0.1 and a Q value [the false discovery rate (FDR) as calculated by the Benjamini–Hochberg procedure] < 0.001 were used as the filtering thresholds. The distributions of the DEGs in various clusters were summarized and illustrated in volcano plots, Venn diagrams, and dot plots via OmniAnalyzer Pro and R software (4.3.1) with the “VennDiagram”, “grid”, “futile.logger”, and “formatR” packages.

2.4 Enrichment analysis

The biological functions, signaling pathways, and diseases associated with the identified DEGs were identified via Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses using R software (4.3.1) with the “clusterProfiler” and “enrichplot” packages.

2.5 Bioinformatics analysis of The Cancer Genome Atlas Prostate Adenocarcinoma (TCGA-PRAD) dataset

The input genes were verified using the TCGA-PRAD dataset. The clinicopathological and follow-up phenotype files were downloaded from the Genomic Data Commons (GDC) Data Portal (https://portal.gdc.cancer.gov/) and the UCSC Xena website (https://xenabrowser.net/datapages/) [49]. The read count data of genes from TCGA-PRAD patients were transformed to transcripts per kilobase of exon model per million (TPM)-mapped read values, and the clinicopathological characteristics of the patients [GS, pathological T stage (pT), and TP53 mutation status] were then compared via R software (4.3.1) with the “limma”, “reshape2”, and “ggpubr” packages. Kaplan–Meier survival curve analysis and Cox regression analysis were performed on the basis of the following prognostic outcomes: recurrence-free survival (RFS) and progression-free survival (PFS). The RFS and PFS data were extracted from the phenotype files as previously described [49], [50]. Survival analysis was conducted via R software (4.3.1) with the “survival” and “survminer” packages.

2.6 Immunohistochemistry (IHC)

Paraffin-embedded tissues were sectioned at 5 µm, subjected to antigen retrieval and incubated with primary antibodies against UQCRB and LBH and with the appropriate secondary antibodies. Information on these antibodies is presented in Table S3.

For statistical analysis, we selected representative files at 36.4 × magnification and assigned a staining intensity score (1, weak; 2, moderate; and 3, strong) and a staining percentage score (0, ≤ 5 % positive cells; 1, 6–25 % positive cells; 2, 26–50 % positive cells; 3, 51–75 % positive cells; and 4, ≥ 76 % positive cells); we then calculated the staining index (SI, staining intensity score × staining percentage score) as previously described [51], [52].

2.7 Cytological analysis

2.7.1 Cell culture

The C4–2, DU145, and PC3 cell lines [American Type Culture Collection (ATCC), Rockville, MD, USA] were cultured in RPMI-1640 or DMEM (Gibco, Grand Island, NY, USA) supplemented with 10 % fetal bovine serum (FBS) (HyClone, South Logan, UT, USA) and 1 % antibiotic-antimycotic (AA) (Gibco, Grand Island, NY, USA) in a humidified atmosphere at 37 °C and 5 % CO2.

2.7.2 Gene regulation in PCa cell lines

To stably knockdown the expression of UQCRB and LBH, PCa cells were infected with lentiviruses harboring short hairpin RNA (shRNA) sequences (GENECHEM, Shanghai, China) targeting specific genes in combination with 5 μg/ml polybrene and subsequently selected with 1–2 μg/ml puromycin. The sequences of the shRNAs are shown in Table S4.

2.7.3 Reverse transcription (RT) and real-time quantitative polymerase chain reaction (qPCR) analysis

TRIzol™ reagent (Invitrogen, Carlsbad, CA, USA) was used to isolate total RNA. Complementary DNA (cDNA) synthesis and qPCR were conducted using One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen Biotech, Beijing, China) with anchored oligo (dT) primers and Top Green qPCR SuperMix (TransGen Biotech) on an SDS 7500 FAST Real-Time PCR system (Applied Biosystems, Foster City, CA, USA). Endogenous reference genes, such as GAPDH and 18S ribosomal RNA, were used for normalization. The sequences of primers used for the related genes are shown in Supplementary Table S5.

2.7.4 Cell viability assay

Cell Counting Kit-8 (CCK-8) [MedChemExpress (MCE), Monmouth Junction, NJ, USA] assay was used to determine cell viability according to the manufacturer’s instructions. A total of 2 × 103 cells were seeded in 96-well plates, and cell viability was assessed at 48 h, 72 h, and 96 h.

2.7.5 Colony formation assay

A total of 2 × 103 cells per well were evenly seeded and cultured in 6-well plates for 2 weeks. The resulting colonies were stained with a 0.01 % crystal violet solution (Beyotime, Shanghai, China) and fully decolorized.

2.7.6 Transwell invasion assay

Cell invasion was measured using Transwell chambers (Corning, NY, USA). A mixture of 20 % FBS in 500 μl of medium was added to the lower chamber to induce cell invasion. A total of 5 × 104 cells in serum-free medium were seeded in the upper chamber, the upper membrane of which was coated with Matrigel (BD Biosciences, San Jose, CA, USA). The cells were then incubated in the chambers for 2 days; unattached cells on the upper surface of the Transwell membranes were removed, and the remaining cells were stained with a 0.01 % crystal violet solution (Beyotime, Shanghai, China).

2.8 Statistical analysis

Statistical analysis was conducted using OmniAnalyzer Pro, R (4.0.3), SPSS version 23 (IBM, Armonk, New York, USA), GraphPad Prism 9.0.0 (GraphPad Software, La Jolla, CA, USA), and Microsoft Excel 2019 (Microsoft Corp., Redmond, WA, USA). The Wilcoxon test, the Mann—Whitney U test, or Student’s t test with P- or Q-value calculations was used to analyze continuous data for comparisons between two groups. The Kruskal—Wallis test or one-way or two-way analysis of variance (ANOVA) followed by post hoc tests (Tukey's test or Dunn's multiple comparison test) was used to analyze continuous data for comparisons among ≥ 3 groups. Survival analysis was performed using Kaplan—Meier curves with the log-rank test and univariate Cox regression models with hazard ratios (HRs). P- or Q-values < 0.05 indicated statistical significance.

3 Results

3.1 Study workflow

Visium ST analysis was performed on PCa tissues collected from patients (Fig. 1A). On the basis of the transcriptomic data obtained from the ST procedure, the spots were dimensionally reduced and classified into groups through PCA and Louvain clustering analysis. These groups were then annotated to various GSs for GE tissues and AFMS tissues (Fig. 1B-C). Subsequently, inferCNV, DPT, and PAGA analyses were conducted to assess the malignancy and developmental trajectories of the different histological clusters (Fig. 1D-E). Genes related to GS progression were identified through DEG analysis by comparing the groups along a GS gradient (Fig. 1F). GO and KEGG enrichment analyses of the identified genes were performed to evaluate their biological functions during GS progression (Fig. 1G). The genes were further verified and screened using the TCGA-PRAD dataset (Fig. 1H). Finally, the selected genes, UQCRB and LBH, were validated through histological and cytological analyses (Fig. 1I).Fig. 1 Study workflow. (A) Visium spatial transcriptome sequencing was performed using PCa tissues on the Visium Spatial Slide. The spatial microarrays contained 4992 spatially barcoded spots, each with a diameter of 55 µm, in capture areas measuring 6.5 × 6.5 mm². (B) The full transcriptome of each spot was captured via the Visium ST procedure. (C) The spots were divided into subgroups through PCA and Louvain clustering analysis. These subgroups were annotated to various GS or AFMS tissues on the basis of their histological structures. (D-E) The malignancy and developmental trajectories of the different histological clusters were assessed via inferCNV (D), DPT and PAGA (E) analyses. (F-G) Genes related to GS progression were identified through DEG analysis (F) and subsequently subjected to enrichment analyses (G) to evaluate their biological functions during GS progression. (H) The DEGs were verified and further screened using the TCGA-PRAD dataset. (I) The selected genes UQCRB and LBH were investigated through histological and cytological analyses.

Fig. 1

3.2 Transcriptomic heterogeneity of PCa tissue

GS heterogeneity is evident in multifocal PCa [53]. In the H&E-stained image of a selected PCa section (patient 2), we annotated the regions with GSs of 3 + 3, 5 + 4, (5 + 5)_1, and (5 + 5)_2 (Fig. 2A). The tissue sections were then classified into various groups through PCA and Louvain clustering analysis according to the transcriptome information. Overlaying the classified spots on the histological image revealed a notable correlation between the transcriptomic profiles and the histologically identifiable structures (Fig. 2B). Regions with a GS of 3 + 3 have been proposed to represent precancerous lesions [41], [42] because of their good differentiation and low capacity for metastasis [43], [44], [45]. InferCNV analysis was performed using GS (3 + 3) as a reference group, allowing the identification of clonal hierarchies and the degree of tumor malignancy [4], [5], [6], [7], [8], [20], [21]. The median inferCNV scores were estimated and found to be higher in regions with high GSs, i.e., GS (5 + 4), GS (5 + 5)_1, and GS (5 + 5)_2, than for the regions GS (3 + 3) (Fig. 2C). In the unsupervised hierarchical clustering analysis, the GS (5 + 5)_1 and GS (5 + 5)_2 groups were clustered together adjacent to the GS (5 + 4) group (Fig. 2D), indicating their similarity and monoclonal origin.Fig. 2 Transcriptomic heterogeneity of PCa tissue. (A) Four different histologically identifiable lesions within H&E-stained PCa sections were enlarged and annotated with GSs of 3 + 3, 5 + 4, (5 + 5)_1, and (5 + 5)_2. (B) PCa sections were classified into various groups through PCA and Louvain clustering analysis and visualized via UMAP (left). The spots observed on the histological images were annotated as GE or AFMS tissues (right). (C) An inferCNV analysis was conducted, and the inferCNV scores for GE clusters are presented in a boxplot. The median CNV scores for individual GE clusters are marked in red, while the first and third quartiles are indicated in black. (D) The chromosomal landscape of inferCNV values for each cluster is shown in a heatmap. The unsupervised hierarchical clustering is shown on the right vertical axis.

Fig. 2

3.3 Spatial transcriptome-wide expression patterns during GS progression

We then performed DPT and PAGA analyses to reconstruct the developmental trajectories of the classified clusters. The geodesic distances of these clusters along the graph illustrate GS progression (from low to high). We also found that highly aggressive tumors in the GSs of 5 + 4 and 5 + 5 groups were densely distributed and that they displayed a certain separation from those in the GS (3 + 3) group (Fig. 3A-B).Fig. 3 Spatial transcriptome-wide expression patterns during GS progression. (A) DPT analysis was performed, and appropriate coordinates of the GE clusters of different GSs were constructed (top). A diffusion map of dpt_pseudotime (bottom) was used to infer the developmental path. (B) The GS developmental trajectories were inferred through PAGA analysis and are illustrated as a static PAGA plot in which the GE subclusters (left) and each spot (right) are retained. (C) The DEGs in the comparisons among GS (5 + 5)_1 vs. GS (5 + 4), GS (5 + 5)_2 vs. GS (5 + 4), and GS (5 + 4) vs. GS (3 + 3) were visualized via volcano plots. Crucial DEGs (including UQCRB and LBH) are highlighted in these plots. (D) The numbers of genes that were significantly and simultaneously upregulated (left) or downregulated (right) in the above comparisons are shown in Venn diagrams. (E) The average expression patterns of these genes in each GS subcluster are shown in a bubble plot. (F) The intersecting genes were subjected to GO (up) and KEGG (down) enrichment analyses and are shown in dot plots. BP: biological process; CC: cellular component; MF: molecular function.

Fig. 3

To investigate the molecular basis of GS progression, we analyzed DEGs among the comparisons of the following groups: GS (5 + 5)_1 vs. GS (5 + 4), GS (5 + 5)_2 vs. GS (5 + 4), and GS (5 + 4) vs. GS (3 + 3) (Fig. 3C). This gradient of comparisons of GSs (from high to low) led to the identification of 63 genes (including JUN, KLK12, and NPY) that were positively associated with GS progression and 26 genes that were negatively associated with it (these genes had |ln FC| > 0.1 and Q values < 0.001 in the Wilcoxon test) (Fig. 3D-E).

We subsequently performed GO and KEGG enrichment analyses of the 89 genes to evaluate their biological functions during GS progression. The most frequently enriched biological functions associated with these genes were immune-like metabolic processes, especially antigen processing and presentation (Fig. 3F).

3.4 LBH expression was positively associated with GS progression in Patient 2 and Patient 1

The same process used to analyze the tissue obtained from Patient 2 was used to analyze the PCa tissue obtained from Patient 1. According to the pathological diagnosis after radical prostatectomy, the GS for Patient 1 was 3 + 4 (Table S2); the tissue sample was further classified into histologically identifiable structures that were marked GSs of 3 + 3, 3 + 4, and 4 + 3 through PCA and Louvain clustering analysis. The numbers of genes whose expression was simultaneously upregulated or downregulated during GS progression are listed (Fig. 4A). Only the LBH gene was found to be upregulated during GS progression in both Patient 2 and Patient 1, suggesting that there is individual specificity of gene expression during GS progression in different patients (Fig. 4A).Fig. 4 Verification of GS progression-associated genes using the TCGA-PRAD dataset. (A) The numbers of genes that were significantly upregulated (top) or downregulated (bottom) during GS progression in Patient 2 and Patient 1 are shown in Venn diagrams. LBH gene expression was upregulated during GS progression in both patients. (B) The 86 genes screened by ST analysis were verified through comparisons of tumor vs. normal, GS > 7 vs. GS < 7, pT3 vs. pT2, and TP53 mutation vs. wild-type tissues, and survival analysis (RFS and PFS) in the TCGA-PRAD dataset. The numbers of dysregulated genes in these comparisons are shown in a Venn diagram. Representative commonly dysregulated genes (CNKSR3, LSM8, UQCRB, and LBH) are shown in the corresponding regions. (C-D) Distributions of the 4 genes in the comparisons of tumor tissue vs. normal tissue (C) and tissues with various GSs (GS < 7 vs. GS = 7 vs. GS > 7) (top of D) and pT stages (pT2 vs. pT3 vs. pT4) (bottom of D) are shown in boxplots. Each box represents the median ± interquartile range. The P value obtained via the Wilcoxon test is shown. (E) Kaplan−Meier survival analysis of RFS (top) and PFS (bottom) was performed for patients grouped according to the optimum cutoff value for each gene (lowest log-rank P value). The Kaplan−Meier survival curves with log-rank P values and the number of patients at risk are shown in the graph.

Fig. 4

3.5 Verification of GS progression-associated genes using the TCGA-PRAD dataset

We evaluated the 89 intersecting genes and identified 86 genes via genomic alignment of the TCGA-PRAD cohort, and we examined how the expression of these genes correlated with various clinicopathological characteristics of TCGA-PRAD patients. In the comparisons of tumor vs. normal tissues, GS > 7 regions vs. GS < 7 regions, and pathological T3 (pT3) tissues vs. pT2 tissues, 61, 39, and 42 genes, respectively (including CRISP3, LBH, NPY, and UQCRB), were significantly differentially expressed (P < 0.05 in the Wilcoxon test) (Tables S6-S8). In PCa, TP53 mutation has been associated with shorter radiographic PFS, tumor metastasis, and time to CRPC [54]. In the comparison between the TP53 mutation group and the wild-type group of TCGA-PARD patients, 24 genes were significantly differentially expressed (Table S9).

The prognostic impacts of the changes in the expression of the 86 identified genes were evaluated through survival analysis. The TCGA-PRAD patients were divided into two groups according to the optimum threshold segmentation of gene expression, which was determined according to the lowest log-rank P value in the Kaplan−Meier analysis. The results revealed that 17 genes were significantly associated with RFS and that 17 genes were significantly associated with PFS (P value < 0.05 in the Cox regression analysis); 13 and 12 of these genes, respectively (including CNKSR3, LSM8, and PPFIA2), were linked to a poor prognosis (HR > 1 and P value < 0.05 in the Cox regression analysis) (Tables S10 and S11).

We subsequently analyzed the simultaneously dysregulated genes identified in the above comparisons (P < 0.05 according to the Wilcoxon test or Cox regression analysis) (Fig. 4B and Table S12). The results revealed that genes such as CNKSR3, LSM8, UQCRB, and LBH were commonly upregulated as GS increased [from GS (3 + 3) to GS (5 + 5)] in the ST analysis and in advanced stages of the disease (in comparisons of tumor vs. normal tissues, GS > 7 vs. GS < 7, pT3 vs. pT2, TP53 mutation vs. wild-type, and HR > 1 in the RFS and PFS analysis in the TCGA-PRAD dataset) (Fig. 4B-E, Fig. 5A, and Table S12). We also found different expression tendencies of some genes (represented by REPS2, AEBP1, C1QA, and C1QB) during GS progression in the ST analysis and in advanced stages in the TCGA-PRAD dataset. This discrepancy may be attributed to the histological atypia of tumor tissues extracted from patients with multifocal PCa, resulting in discordance between the transcriptome information and the actual pathological diagnosis reported in the TCGA-PRAD dataset.Fig. 5 Expression patterns of UQCRB and LBH during GS progression. (A) The expression patterns of four genes (CNKSR3, LSM8, UQCRB, and LBH) are shown in a bubble plot. (B) The expression levels of UQCRB and LBH are shown in spatial activity maps. (C) The distribution and expression levels of UQCRB and LBH are shown in PAGA plots based on the GE subclusters (left) and spots (right).

Fig. 5

We selected UQCRB and LBH for further verification through histological and cytological experiments, as antibodies to the UQCRB and LBH proteins that can be used in IHC are commercially available. Bubble plots of gene expression, spatial activity maps, and PAGA graphs of the expression of the two genes in samples with different GSs revealed a pattern of gradual upregulation of their expression during GS progression (Fig. 5A-C).

3.6 The UQCRB and LBH proteins are highly expressed in advanced-stage PCa tissues

The expression levels of the proteins encoded by UQCRB and LBH genes in prostatic intraepithelial neoplasia (PIN) and PCa tissues were assessed through IHC. Representative IHC images revealed that the expression of UQCRB and LBH proteins increased gradually during GS progression [from GS (3 + 3) to GS (5 + 4)] (Fig. 6A). The percentages of stained cells and the staining indices for UQCRB and LBH proteins were quantified in PCa tissues with different GSs (GS < 7, GS = 7, and GS > 7) and different pTs [pT2, pT3a and pT3b (pT3a/b), and pT4]. The results indicated that the expression of both proteins was positively associated with GS and pT progression (Fig. 6B).Fig. 6 Expression patterns of UQCRB and LBH proteins in PCa tissues stratified by various clinicopathological features, i.e., GS and pT stage. (A) Representative images showing IHC staining for UQCRB (left) and LBH (right) proteins in PCa tissues with PIN and different GSs [GS (3 + 3), GS (3 + 4), GS (4 + 4), and GS (5 + 4)]. (B) Quantitative analysis of the staining percentage and staining index of regions with various GSs (top) and at various pT stages (bottom) are shown in the histograms. The mean ± standard error of the mean (SEM) is shown in each column. ns, P > 0.05; * , P < 0.05; * *, P < 0.01; * ** , P < 0.001.

Fig. 6

3.7 Knockdown of UQCRB and LBH mRNAs suppresses PCa cell proliferation and invasion

To evaluate the effects of UQCRB and LBH gene expression on cell viability, we conducted cytological analyses using C4–2, DU145, and PC3 cells. UQCRB and LBH were knocked down in these cell lines, and the knockdown efficiency was assessed through qPCR analysis (Figs. S1A-B). The results of the CCK-8, colony formation, and Transwell assays indicated that downregulation of UQCRB and LBH expression significantly suppressed the proliferation and invasion of PCa cell lines (Fig. 7A-C).Fig. 7 Knockdown of UQCRB and LBH mRNAs suppresses PCa cell proliferation and invasion. (A-C) UQCRB and LBH were knocked down in C4–2, DU145, and PC3 cells via lentiviral transduction. CCK-8 (A), colony formation (B), and Transwell (C) assays were conducted to assess PCa cell proliferation and invasion. Ctrl, the parental cells; sh_NC, negative control cells; sh_UQCRB or sh_LBH, UQCRB or LBH stably downregulated cells.

Fig. 7

4 Discussion

PCa cells are diffusely distributed and form multiple foci in prostate tissue. PCa most commonly originates from basal or luminal cells of the GE and is surrounded by nonglandular tissues, mainly AFMS [12], [13], [14], [15], [16], [17], [18]. For this reason, it is difficult to distinguish PCa tissues through visual assessment, making pathological diagnosis and sequencing analysis challenging.

Visium ST technology makes it possible to quantify the entire transcriptome at high spatial resolution; it can be used to precisely determine the gene expression patterns in specific histological structures and is therefore particularly appropriate for analyzing multifocal PCa.

Few studies have focused on the ST analysis of PCa [8], [25], [26], [27]. Berglund et al. published spatial maps of PCa transcriptomes in different histological structures [25]. However, they used conventional low-resolution ST (spatial microarrays with 1007 spots); thus, it was difficult to precisely determine the histological structures and cell features [25]. Erickson et al. observed gene copy number alterations during PCa progression through Visium ST analysis. Their work focused on genome integrity during tumor evolution and did not include analysis of transcriptome-wide expression patterns [8]. Watanabe et al. performed Visium CytAssist ST analysis and revealed the gene expression patterns of neuroendocrine prostate carcinoma (NEPC) [26]. However, their study screened DEGs between androgen receptor pathway-positive prostate carcinoma (ARPC) and NEPC without elucidating the developmental trajectories from ARPC to NEPC [26]. We also previously analyzed the transcriptome-wide expression patterns of PCa lesions with different GSs through Visium ST analysis [27]; PCA was performed on the whole genome, but the screened GS progression-related genes were not validated through histological experiments [27].

To analyze the transcriptome-wide expression patterns that occur during the evolutionary process, it is important to determine the cell origin and the relationships among disparate lesions within the same gland. Acquired specific genetic changes have been reported to accumulate during the progression of PCa to advanced stages [4], [5], [6], [7], [8], [9], [40]. While such additional genomic instabilities trigger the formation of distinct subclonal phenotypes, they are distinct from the dominant conserved genomic alterations that are necessarily shared by all lesions. High-resolution genome-wide analysis has shown that separate foci from the same patient always trace to a single genomically aberrant PCa cell whose macrogenomic CNV values stably replicate with each cell division [4], [5], [8], suggesting that multifocal PCa has a monoclonal origin [4], [5], [6], [7], [8].

On the basis of the above theory, we consider that the multiple lesions present in multifocal PCa represent relative developmental stages that evolve from the same initial PCa precursor. A previous study using ST analysis on PCa tissue sections revealed that PCA based on the 500 most variable genes clearly distinguished the different histological structures in PCa [25], [27]. In the PCA conducted in our study, we used HVGs rather than all genes, as in our previous study [27]. HVGs exhibit highly dynamic patterns of expression across several biological conditions, and these patterns enable the precise identification of different histological structures [55]. Consistent with our expectations, the classified subgroups overlaid on H&E-stained histological images closely mirrored structures with GSs of 3 + 3, 5 + 4, and 5 + 5.

We subsequently assessed the reliability of the PCA-determined histological characteristics through bioinformatic analysis. The average inferCNV value was used as the criterion for identifying malignant stages of PCa [20], [21]. We found that the median inferCNV value was significantly greater for aggressive lesions with GSs of 5 + 4 and 5 + 5 than for lesions with GS (3 + 3), a finding that confirms the accuracy of macroscopically identified GS when verified through malignant criteria. The nearest clonal hierarchies of lesions whose GSs indicated aggressiveness in the inferCNV clustering analysis provided additional evidence of their histological similarities and monoclonal origin. We then evaluated the biological processes associated with the identified subgroups by predicting their developmental trajectories through DPT and PAGA analyses. The geodesic distances along the graph illustrated that the progression patterns of the classified GE regions matched those of the GS regions. The results of all these analyses confirmed the agreement between the results of visual assessment of histological structures and transcriptome-wide expression pattern analysis.

After preliminarily defining our annotated GS clusters as groups of lesions representing tumors that exhibit different degrees of malignancy and are at different stages of progression, we performed pairwise comparisons of groups along a GS gradient to investigate GS progression-related genes. We discovered 89 DEGs that were dysregulated during GS progression and found that these genes were enriched in immune-like metabolic processes such as antigen processing and presentation. Recent studies have shown that the immune checkpoint molecule B7-H3 (CD276) is correlated with GS progression and poor outcomes in prostatectomy cohorts [56] and is a potential therapeutic target in advanced PCa [57], [58]. However, arguably due to the presence of immunosuppressive mechanisms [24], [59], [60], [61], the efficacy of immune checkpoint inhibitor monotherapy is still controversial [62]. Therefore, novel therapeutic targets need to be identified through further study of GS progression-related genes.

By verifying the associations of the identified 86 DEGs with clinicopathological parameters and prognosis in the TCGA-PRAD dataset, we found that genes such as CNKSR3, LSM8, UQCRB, and LBH were commonly and consistently upregulated during progression to advanced stages in both the ST analysis and the TCGA-PRAD dataset analysis. The changes in the expression of UQCRB and LBH were particularly appropriate for histological and cytological verification because of the feasibility of antibody acquisition and RNA interference. Further experiments confirmed a progressive increase in the expression of the proteins encoded by UQCRB and LBH as GS and pT stage advanced, as well as the suppressive effects of mRNA knockdown of these genes on PCa cell proliferation and invasion.

UQCRB is a functional nuclear gene that encodes a 13.3-kDa subunit of mitochondrial complex III [63]. To our knowledge, no study has elucidated the role of UQCRB in PCa. In breast and colon cancer cells, silencing UQCRB expression reversed the epithelial–mesenchymal transition (EMT) and autophagic flux by blocking the production of mitochondrial reactive oxygen species (ROS) [64], [65]. In glioblastoma stem‑like cells, knockdown of UQCRB significantly inhibited cancer stem cell-like phenotypes, as well as the expression of stemness markers, by blocking the mitochondrial ROS/HIF‑1α/c‑Met pathway [66]. A study revealed that mitochondrial UQCRB positively regulates VEGFR2 signaling in endothelial cells (ECs) and that the UQCRB-targeting agent terpestacin could be applied in new therapeutic approaches for human cancer [67].

LBH is a highly conserved and novel transcription cofactor that plays important roles in heart and limb development in vertebrates [68]. LBH is dysregulated by DNA hypomethylation, and its expression is correlated with WNT activation in multiple cancers [69], [70]; LBH expression also indicates a poor prognosis and promotes gastric cancer cell proliferation and invasion by activating the integrin/FAK/Akt pathway [71]. In triple-negative [ER (-), PR (-), and HER2 (-)] breast cancer (TNBC) cells, blocking LBH expression disrupts genome integrity by inducing replicative stress and activating the ATR/CHK1 DNA damage response pathway [72]. In human glioma, LBH promotes angiogenesis via VEGFA-mediated ERK signaling under hypoxia; this process is transcriptionally regulated by HIF-1 and results in a self-reinforcing cycle [73]. However, LBH is reported to inhibit cell proliferation and invasion in nasopharyngeal carcinoma by inducing G1/S cell cycle arrest and downregulating αB-crystallin expression [74], [75]. Therefore, owing to the various functions of LBH in different cancers, further studies are needed to clarify its role in the development of PCa.

Several limitations of this study need to be addressed in future investigations. First, most of the genes identified in the ST analysis in this study were derived from one PCa patient, which is a potential source of error due to individual heterogeneity. Second, the identification of two distinct GS (5 + 5) clusters with different transcriptome expression patterns suggests that relying solely on visual GS classification may have limited value due to the inherent imperfections of histological classification. The optimal analysis requires precisely categorizing PCa sections based solely on their Gleason grades (e.g., Gleason pattern 3 to Gleason pattern 5), a task that is challenging to accomplish. Finally, further research is needed to clarify the roles of the identified genes (such as UQCRB and LBH) in the molecular mechanisms and signaling pathways that are active during GS progression in PCa.

5 Conclusions

In this study, we analyzed the transcriptome-wide expression patterns of GE cells with different GSs using Visium ST technology. Diverse histological structures corresponding to different GSs were observed across spots subgrouped through PCA and Louvain clustering analysis, alongside bioinformatics analysis linked to their respective malignant developmental stages. We discovered crucial oncogenes, two of which are UQCRB and LBH, whose expression gradually increased during GS progression; the carcinogenic properties of these genes were verified using the TCGA-PRAD dataset, as well as through histological and cytological experiments. In summary, the results of this PCa spatial analysis establish a foundation for ST-based studies of GS progression and provide valuable insights into the GS progression-related genes UQCRB and LBH.

Ethics approval and consent to participate

Human prostate tissue samples were obtained from Beijing Tongren Hospital. Ethical consent for the use of these samples in the study was obtained from the Ethics Committee of Beijing Tongren Hospital, which is affiliated with Capital Medical University.

Funding

This work was supported by the 10.13039/501100001809 National Natural Science Foundation of China (grant numbers: 82072833 and 82272864 to Hao Ping) and by the Capital Health Research and Development of Special Fund, Beijing, China (grant number: 2024-2-2059 ).

CRediT authorship contribution statement

Hong Zhang: Methodology, Investigation. Yongjun Quan: Supervision, Software, Resources, Methodology, Investigation, Formal analysis, Data curation. Hao Ping: Visualization, Validation, Funding acquisition, Data curation. Mingdong Wang: Software, Project administration, Investigation.

Declaration of Competing Interest

The authors declare that they have no competing interests.

Appendix A Supplementary material

Table S1

Supplementary material.

.

Fig S1

Verification of the knockdown efficiency of UQCRB and LBH in C4–2, PC3, and DU145 cells. (A-B) The efficiencies of UQCRB (A) and LBH (B) knockdown in C4–2, DU145, and PC3 cells were assessed through qPCR analysis. Ctrl, the parental cells; sh_NC, negative control cells; sh_UQCRB or sh_LBH, UQCRB or LBH stably downregulated cells.

.

Acknowledgments

We thank the members of the Urology and Pathology Department of Beijing Tongren Hospital, Capital Medical University, for assisting with sample collection and providing guidance on pathological diagnoses.

Author contributions

Hao Ping and Yongjun Quan designed the study. Hao Ping provided the clinical PCa samples. Yongjun Quan and Mingdong Wang collected the samples. Yongjun Quan performed the bioinformatic analysis and histological and cytological experiments. Hong Zhang performed the pathological diagnosis. Yongjun Quan wrote the manuscript, which was approved by all the authors.

Consent for publication

All the authors have read and approved the final manuscript for publication.

☆ A full list of author information is provided in the electronic submission form.

Appendix A Supplementary data associated with this article can be found in the online version at doi:10.1016/j.csbj.2024.08.026.
==== Refs
References

1 Siegel R.L. Giaquinto A.N. Jemal A. Cancer statistics, 2024 CA: a Cancer J Clin 74 2024 12 49 10.3322/caac.21820
2 Andreoiu M. Cheng L. Multifocal prostate cancer: biologic, prognostic, and therapeutic implications Hum Pathol 41 2010 781 793 10.1016/j.humpath.2010.02.011 20466122
3 Haffner M.C. Zwart W. Roudier M.P. True L.D. Nelson W.G. Epstein J.I. Genomic and phenotypic heterogeneity in prostate cancer Nat Rev Urol 18 2021 79 92 10.1038/s41585-020-00400-w 33328650
4 Boyd L.K. Mao X. Xue L. Lin D. Chaplin T. Kudahetti S.C. High-resolution genome-wide copy-number analysis suggests a monoclonal origin of multifocal prostate cancer Genes, Chromosomes Cancer 51 2012 579 589 10.1002/gcc.21944 22334418
5 Gundem G. Van Loo P. Kremeyer B. Alexandrov L.B. Tubio J.M.C. Papaemmanuil E. The evolutionary history of lethal metastatic prostate cancer Nature 520 2015 353 357 10.1038/nature14347 25830880
6 Liu W. Laitinen S. Khan S. Vihinen M. Kowalski J. Yu G. Copy number analysis indicates monoclonal origin of lethal metastatic prostate cancer Nat Med 15 2009 559 565 10.1038/nm.1944 19363497
7 Hong M.K. Macintyre G. Wedge D.C. Van Loo P. Patel K. Lunke S. Tracking the origins and drivers of subclonal metastatic expansion in prostate cancer Nat Commun 6 2015 6605 10.1038/ncomms7605 25827447
8 Erickson A. He M. Berglund E. Marklund M. Mirzazadeh R. Schultz N. Spatially resolved clonal copy number alterations in benign and malignant tissue Nature 608 2022 360 367 10.1038/s41586-022-05023-2 35948708
9 Robinson D. Van Allen E.M. Wu Y.M. Schultz N. Lonigro R.J. Mosquera J.M. Integrative clinical genomics of advanced prostate cancer Cell 161 2015 1215 1228 10.1016/j.cell.2015.05.001 26000489
10 Mateo J. Carreira S. Sandhu S. Miranda S. Mossop H. Perez-Lopez R. DNA-Repair defects and olaparib in metastatic prostate cancer N Engl J Med 373 2015 1697 1708 10.1056/NEJMoa1506859 26510020
11 de Bono J. Mateo J. Fizazi K. Saad F. Shore N. Sandhu S. Olaparib for metastatic castration-resistant prostate cancer N Engl J Med 382 2020 2091 2102 10.1056/NEJMoa1911440 32343890
12 Goldstein A.S. Huang J. Guo C. Garraway I.P. Witte O.N. Identification of a cell of origin for human prostate cancer Science 329 2010 568 571 10.1126/science.1189992 20671189
13 Wang X. Kruithof-de Julio M. Economides K.D. Walker D. Yu H. Halili M.V. A luminal epithelial stem cell that is a cell of origin for prostate cancer Nature 461 2009 495 500 10.1038/nature08361 19741607
14 Wang Z.A. Mitrofanova A. Bergren S.K. Abate-Shen C. Cardiff R.D. Califano A. Lineage analysis of basal epithelial cells reveals their unexpected plasticity and supports a cell-of-origin model for prostate cancer heterogeneity Nat Cell Biol 15 2013 274 283 10.1038/ncb2697 23434823
15 Lawson D.A. Zong Y. Memarzadeh S. Xin L. Huang J. Witte O.N. Basal epithelial stem cells are efficient targets for prostate cancer initiation Proc Natl Acad Sci USA 107 2010 2610 2615 10.1073/pnas.0913873107 20133806
16 Lee S.H. Shen M.M. Cell types of origin for prostate cancer Curr Opin Cell Biol 37 2015 35 41 10.1016/j.ceb.2015.10.002 26506127
17 Xin L. Cells of origin for cancer: an updated view from prostate cancer Oncogene 32 2013 3655 3663 10.1038/onc.2012.541 23178496
18 Taylor R.A. Toivanen R. Frydenberg M. Pedersen J. Harewood L. Collins A.T. Human epithelial basal cells are cells of origin of prostate cancer, independent of CD133 status Stem Cells (Dayt, Ohio) 30 2012 1087 1096 10.1002/stem.1094
19 Chen S. Zhu G. Yang Y. Wang F. Xiao Y.T. Zhang N. Single-cell analysis reveals transcriptomic remodellings in distinct cell types that contribute to human prostate cancer progression Nat Cell Biol 23 2021 87 98 10.1038/s41556-020-00613-6 33420488
20 Dong B. Miao J. Wang Y. Luo W. Ji Z. Lai H. Single-cell analysis supports a luminal-neuroendocrine transdifferentiation in human prostate cancer Commun Biol 3 2020 778 10.1038/s42003-020-01476-1 33328604
21 Ma X. Guo J. Liu K. Chen L. Liu D. Dong S. Identification of a distinct luminal subgroup diagnosing and stratifying early stage prostate cancer by tissue-based single-cell RNA sequencing Mol Cancer 19 2020 147 10.1186/s12943-020-01264-9 33032611
22 Taavitsainen S. Engedal N. Cao S. Handle F. Erickson A. Prekovic S. Single-cell ATAC and RNA sequencing reveal pre-existing and persistent cells associated with prostate cancer relapse Nat Commun 12 2021 5307 10.1038/s41467-021-25624-1 34489465
23 Barros-Silva J.D. Linn D.E. Steiner I. Guo G. Ali A. Pakula H. Single-cell analysis identifies LY6D as a marker linking castration-resistant prostate luminal cells to prostate progenitors and cancer e3506 Cell Rep 25 2018 3504 3518 10.1016/j.celrep.2018.11.069 30566873
24 Wong H.Y. Sheng Q. Hesterberg A.B. Croessmann S. Rios B.L. Giri K. Single cell analysis of cribriform prostate cancer reveals cell intrinsic and tumor microenvironmental pathways of aggressive disease Nat Commun 13 2022 6036 10.1038/s41467-022-33780-1 36229464
25 Berglund E. Maaskola J. Schultz N. Friedrich S. Marklund M. Bergenstråhle J. Spatial maps of prostate cancer transcriptomes reveal an unexplored landscape of heterogeneity Nat Commun 9 2018 2419 10.1038/s41467-018-04724-5 29925878
26 Watanabe R. Miura N. Kurata M. Kitazawa R. Kikugawa T. Saika T. Spatial gene expression analysis reveals characteristic gene expression patterns of De novo neuroendocrine prostate cancer coexisting with androgen receptor pathway prostate cancer Int J Mol Sci 24 2023 10.3390/ijms24108955
27 Quan Y. Zhang H. Wang M. Ping H. Visium spatial transcriptomics reveals intratumor heterogeneity and profiles of Gleason score progression in prostate cancer iScience 26 2023 108429 10.1016/j.isci.2023.108429
28 Satija R. Farrell J.A. Gennert D. Schier A.F. Regev A. Spatial reconstruction of single-cell gene expression data Nat Biotechnol 33 2015 495 502 10.1038/nbt.3192 25867923
29 Achim K. Pettit J.B. Saraiva L.R. Gavriouchkina D. Larsson T. Arendt D. High-throughput spatial mapping of single-cell RNA-seq data to tissue of origin Nat Biotechnol 33 2015 503 509 10.1038/nbt.3209 25867922
30 Ståhl P.L. Salmén F. Vickovic S. Lundmark A. Navarro J.F. Magnusson J. Visualization and analysis of gene expression in tissue sections by spatial transcriptomics Science 353 2016 78 82 10.1126/science.aaf2403 27365449
31 Gleason D.F. Mellinger G.T. Prediction of prognosis for prostatic adenocarcinoma by combined histological grading and clinical staging J Urol 111 1974 58 64 10.1016/s0022-5347(17)59889-4 4813554
32 Vis A.N. Roemeling S. Kranse R. Schröder F.H. van der Kwast T.H. Should we replace the Gleason score with the amount of high-grade prostate cancer? Eur Urol 51 2007 931 939 10.1016/j.eururo.2006.07.051 16935413
33 Gleason D.F. Classification of prostatic carcinomas Cancer Chemother Rep 50 1966 125 128 5948714
34 Andrén O. Fall K. Franzén L. Andersson S.O. Johansson J.E. Rubin M.A. How well does the Gleason score predict prostate cancer death? A 20-year followup of a population based cohort in Sweden J Urol 175 2006 1337 1340 10.1016/s0022-5347(05)00734-2 16515993
35 Egevad L. Granfors T. Karlberg L. Bergh A. Stattin P. Prognostic value of the Gleason score in prostate cancer BJU Int 89 2002 538 542 10.1046/j.1464-410x.2002.02669.x 11942960
36 Dong F. Wang C. Farris A.B. Wu S. Lee H. Olumi A.F. Impact on the clinical outcome of prostate cancer by the 2005 international society of urological pathology modified Gleason grading system Am J Surg Pathol 36 2012 838 843 10.1097/PAS.0b013e3182486faf 22592143
37 Popiolek M. Rider J.R. Andrén O. Andersson S.O. Holmberg L. Adami H.O. Natural history of early, localized prostate cancer: a final report from three decades of follow-up Eur Urol 63 2013 428 435 10.1016/j.eururo.2012.10.002 23084329
38 Eggener S.E. Scardino P.T. Walsh P.C. Han M. Partin A.W. Trock B.J. Predicting 15-year prostate cancer specific mortality after radical prostatectomy J Urol 185 2011 869 875 21239008
39 Kozminski M.A. Tomlins S. Cole A. Singhal U. Lu L. Skolarus T.A. Standardizing the definition of adverse pathology for lower risk men undergoing radical prostatectomy e411-416 Urol Oncol 34 2016 415 10.1016/j.urolonc.2016.03.019
40 Ruiz C. Lenkiewicz E. Evers L. Holley T. Robeson A. Kiefer J. Advancing a clinically relevant perspective of the clonal nature of cancer Proc Natl Acad Sci USA 108 2011 12054 12059 10.1073/pnas.1104009108 21730190
41 Eggener S.E. Berlin A. Vickers A.J. Paner G.P. Wolinsky H. Cooperberg M.R. Low-Grade Prostate Cancer: time to stop calling it cancer J Clin Oncol: J Am Soc Clin Oncol 40 2022 3110 3114 10.1200/jco.22.00123
42 Labbate C.V. Paner G.P. Eggener S.E. Should Grade Group 1 (GG1) be called cancer? World J Urol 40 2022 15 19 10.1007/s00345-020-03583-4 33432506
43 Ross H.M. Kryvenko O.N. Cowan J.E. Simko J.P. Wheeler T.M. Epstein J.I. Do adenocarcinomas of the prostate with Gleason score (GS) ≤6 have the potential to metastasize to lymph nodes? Am J Surg Pathol 36 2012 1346 1352 10.1097/PAS.0b013e3182556dcd 22531173
44 Donin N.M. Laze J. Zhou M. Ren Q. Lepor H. Gleason 6 prostate tumors diagnosed in the PSA era do not demonstrate the capacity for metastatic spread at the time of radical prostatectomy Urology 82 2013 148 152 10.1016/j.urology.2013.03.054 23706588
45 Kweldam C.F. Wildhagen M.F. Bangma C.H. van Leenders G.J. Disease-specific death and metastasis do not occur in patients with Gleason score ≤6 at radical prostatectomy BJU Int 116 2015 230 235 10.1111/bju.12879 25060593
46 Haghverdi L. Büttner M. Wolf F.A. Buettner F. Theis F.J. Diffusion pseudotime robustly reconstructs lineage branching Nat Methods 13 2016 845 848 10.1038/nmeth.3971 27571553
47 Wolf F.A. Angerer P. Theis F.J. SCANPY: large-scale single-cell gene expression data analysis Genome Biol 19 2018 15 10.1186/s13059-017-1382-0 29409532
48 Wolf F.A. Hamey F.K. Plass M. Solana J. Dahlin J.S. Göttgens B. PAGA: graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells Genome Biol 20 2019 59 10.1186/s13059-019-1663-x 30890159
49 Liu J. Lichtenberg T. Hoadley K.A. Poisson L.M. Lazar A.J. Cherniack A.D. An integrated TCGA pan-cancer clinical data resource to drive high-quality survival outcome analytics e411 Cell 173 2018 400 416 10.1016/j.cell.2018.02.052 29625055
50 Quan Y. Zhang X. Ping H. Construction of a risk prediction model using m6A RNA methylation regulators in prostate cancer: comprehensive bioinformatic analysis and histological validation Cancer Cell Int 22 2022 33 10.1186/s12935-021-02438-1 35045837
51 Wang L. Li J. Hou J. Li M. Cui X. Li S. p53 expression but not p16(INK4A) correlates with human papillomavirus-associated esophageal squamous cell carcinoma in Kazakh population. Infect Agent Cancer 11 2016 19 27076841
52 Wang L. Zhang Z. Yu X. Huang X. Liu Z. Chai Y. Unbalanced YAP-SOX9 circuit drives stemness and malignant progression in esophageal squamous cell carcinoma Oncogene 38 2019 2042 2055 10.1038/s41388-018-0476-9 30401982
53 Arora R. Koch M.O. Eble J.N. Ulbright T.M. Li L. Cheng L. Heterogeneity of Gleason grade in multifocal adenocarcinoma of the prostate Cancer 100 2004 2362 2366 10.1002/cncr.20243 15160339
54 Deek M.P. Van der Eecken K. Phillips R. Parikh N.R. Isaacsson Velho P. Lotan T.L. The mutational landscape of metastatic castration-sensitive prostate cancer: the spectrum theory revisited Eur Urol 80 2021 632 640 10.1016/j.eururo.2020.12.040 33419682
55 Arora J.K. Opasawatchai A. Teichmann S.A. Matangkasombut P. Charoensawan V. Computational workflow for investigating highly variable genes in single-cell RNA-seq across multiple time points and cell types STAR Protoc 4 2023 102387 10.1016/j.xpro.2023.102387
56 Benzon B. Zhao S.G. Haffner M.C. Takhar M. Erho N. Yousefi K. Correlation of B7-H3 with androgen receptor, immune pathways and poor outcome in prostate cancer: an expression-based analysis Prostate Cancer Prostatic Dis 20 2017 28 35 10.1038/pcan.2016.49 27801901
57 Guo C. Figueiredo I. Gurel B. Neeb A. Seed G. Crespo M. B7-H3 as a therapeutic target in advanced prostate cancer Eur Urol 83 2023 224 238 10.1016/j.eururo.2022.09.004 36114082
58 Shenderov E. De Marzo A.M. Lotan T.L. Wang H. Chan S. Lim S.J. Neoadjuvant enoblituzumab in localized prostate cancer: a single-arm, phase 2 trial Nat Med 29 2023 888 897 10.1038/s41591-023-02284-w 37012549
59 Kiniwa Y. Miyahara Y. Wang H.Y. Peng W. Peng G. Wheeler T.M. CD8+ Foxp3+ regulatory T cells mediate immunosuppression in prostate cancer Clin Cancer Res: J Am Assoc Cancer Res 13 2007 6947 6958 10.1158/1078-0432.ccr-07-0842
60 Runcie K.D. Dallos M.C. Prostate cancer immunotherapy-finally in from the cold? Curr Oncol Rep 23 2021 88 10.1007/s11912-021-01084-0 34125308
61 Cha H.R. Lee J.H. Ponnazhagan S. Revisiting immunotherapy: a focus on prostate cancer Cancer Res 80 2020 1615 1623 10.1158/0008-5472.can-19-2948 32066566
62 Kwon E.D. Drake C.G. Scher H.I. Fizazi K. Bossi A. van den Eertwegh A.J. Ipilimumab versus placebo after radiotherapy in patients with metastatic castration-resistant prostate cancer that had progressed after docetaxel chemotherapy (CA184-043): a multicentre, randomised, double-blind, phase 3 trial Lancet Oncol 15 2014 700 712 10.1016/s1470-2045(14)70189-5 24831977
63 Malaney S. Heng H.H. Tsui L.C. Shi X.M. Robinson B.H. Localization of the human gene encoding the 13.3-kDa subunit of mitochondrial complex III (UQCRB) to 8q22 by in situ hybridization Cytogenet Cell Genet 73 1996 297 299 10.1159/000134360 8751380
64 Monti E. Mancini A. Marras E. Gariboldi M.B. Targeting mitochondrial ROS production to reverse the epithelial-mesenchymal transition in breast cancer cells Curr Issues Mol Biol 44 2022 5277 5293 10.3390/cimb44110359 36354671
65 Gwak E.J. Kim D. Hwang H.Y. Kwon H.J. Mitochondrial ROS produced in human colon carcinoma associated with cell survival via autophagy Cancers 14 2022 10.3390/cancers14081883
66 Jung N. Kwon H.J. Jung H.J. Downregulation of mitochondrial UQCRB inhibits cancer stem cell-like properties in glioblastoma Int J Oncol 52 2018 241 251 10.3892/ijo.2017.4191 29115404
67 Jung H.J. Kim Y. Chang J. Kang S.W. Kim J.H. Kwon H.J. Mitochondrial UQCRB regulates VEGFR2 signaling in endothelial cells J Mol Med (Berl, Ger) 91 2013 1117 1128 10.1007/s00109-013-1049-6
68 Briegel K.J. Joyner A.L. Identification and characterization of Lbh, a novel conserved nuclear protein expressed during early limb and heart development Dev Biol 233 2001 291 304 10.1006/dbio.2001.0225 11336496
69 Young I.C. Brabletz T. Lindley L.E. Abreu M. Nagathihalli N. Zaika A. Multi-cancer analysis reveals universal association of oncogenic LBH expression with DNA hypomethylation and WNT-Integrin signaling pathways Cancer Gene Ther 30 2023 1234 1248 10.1038/s41417-023-00633-y 37268816
70 Rieger M.E. Sims A.H. Coats E.R. Clarke R.B. Briegel K.J. The embryonic transcription cofactor LBH is a direct target of the Wnt signaling pathway in epithelial development and in aggressive basal subtype breast cancers Mol Cell Biol 30 2010 4267 4279 10.1128/mcb.01418-09 20606007
71 Yu R. Li Z. Zhang C. Song H. Deng M. Sun L. Elevated limb-bud and heart development (LBH) expression indicates poor prognosis and promotes gastric cancer cell proliferation and invasion via upregulating Integrin/FAK/Akt pathway PeerJ 7 2019 e6885 10.7717/peerj.6885
72 Garikapati K. Young I.C. Hong S. Rai P. Jain C. Briegel K.J. Blocking LBH expression causes replication stress and sensitizes triple-negative breast cancer cells to ATR inhibitor treatment Oncogene 2024 10.1038/s41388-024-02951-3
73 Jiang Y. Zhou J. Zou D. Hou D. Zhang H. Zhao J. Overexpression of Limb-Bud and Heart (LBH) promotes angiogenesis in human glioma via VEGFA-mediated ERK signalling under hypoxia EBioMedicine 48 2019 36 48 10.1016/j.ebiom.2019.09.037 31631037
74 Liu Q. Guan X. Lv J. Li X. Wang Y. Li L. Limb-bud and Heart (LBH) functions as a tumor suppressor of nasopharyngeal carcinoma by inducing G1/S cell cycle arrest Sci Rep 5 2015 7626 10.1038/srep07626 25557837
75 Wu A. Zhang L. Luo N. Zhang L. Li L. Liu Q. Limb-bud and heart (LBH) inhibits cellular migration, invasion and epithelial-mesenchymal transition in nasopharyngeal carcinoma via downregulating αB-crystallin expression Cell Signal 85 2021 110045 10.1016/j.cellsig.2021.110045
