SWATH-MS Based Proteomic Profiling of Prostate Cancer Cells Reveals Adaptive Molecular Mechanisms in Response to Anti-Androgen Therapy

Simple Summary Androgen targeted therapy has been foundational in the management of advanced prostate cancer. Nevertheless, responses to the therapy are found to be seldom sustained, in which patients develop resistance and progress to a lethal and incurable castration-resistance stage. Therefore, comprehensive understanding of the molecular basis of treatment induced cellular responses is required to circumvent molecular mechanisms driving castration-resistance. Using an advanced, robust quantitative proteomic analysis, we profiled the prostate cancer cell proteome induced by the two AR antagonists/anti-androgens, Bicalutamide and Enzalutamide. We highlighted key molecular signatures and cellular pathways that provide insights into the anti-androgen induced adaptive cellular programming in prostate cancer cells. Targeting these molecules and associated pathways might be useful in developing novel therapeutic approaches and/or as biomarkers of predicting prostate cancer treatment response. Abstract Prostate cancer (PCa) is the second most common cancer affecting men worldwide. PCa shows a broad-spectrum heterogeneity in its biological and clinical behavior. Although androgen targeted therapy (ATT) has been the mainstay therapy for advanced PCa, it inevitably leads to treatment resistance and progression to castration resistant PCa (CRPC). Thus, greater understanding of the molecular basis of treatment resistance and CRPC progression is needed to improve treatments for this lethal phenotype. The current study interrogated both proteomics and transcriptomic alterations stimulated in AR antagonist/anti-androgen (Bicalutamide and Enzalutamide) treated androgen-dependent cell model (LNCaP) in comparison with androgen-independent/castration-resistant cell model (C4-2B). The analysis highlighted the activation of MYC and PSF/SFPQ oncogenic upstream regulators in response to the anti-androgen treatment. Moreover, the study revealed anti-androgen induced genes/proteins related to transcription/translation regulation, energy metabolism, cell communication and signaling cascades promoting tumor growth and proliferation. In addition, these molecules were found dysregulated in PCa clinical proteomic and transcriptomic datasets, suggesting their potential involvement in PCa progression. In conclusion, our study provides key molecular signatures and associated pathways that might contribute to CRPC progression despite treatment with anti-androgens. Such molecular signatures could be potential therapeutic targets to improve the efficacy of existing therapies and/or predictive/prognostic value in CRPC for treatment response.


Introduction
Despite the advances in prostate cancer (PCa) diagnosis and treatment, temporal trends of PCa incidence and mortality varied significantly globally during the past few years [1]. Although, localized disease is easily managed with curative surgery or radiotherapy, some patients recur systemically with advanced metastatic PCa. Due to the androgen-dependency and the crucial role of the androgen receptor (AR) in PCa progression, androgen deprivation therapy (ADT) has been denoted as the mainstay therapy for relapsed or metastatic patients [2]. ADT has been progressively supplemented with inhibitors of androgen synthesis and AR antagonists/anti-androgens, therefore termed as androgen targeted therapy (ATT) [3]. Although ADT/ATT is initially effective, many patients invariably fail to respond after~13 months and eventually progress to a lethal stage: androgen-independent/castration-resistant PCa (CRPC) [4].
PCa is an extremely heterogeneous disease, hence understanding the cellular basis of cell heterogeneity is important to comprehend the ADT/ATT resistance during CRPC progression. Complex AR-dependent adaptive mechanisms have been identified that cause reactivation of AR signaling enabling the continued growth and progression of CRPC. This includes AR variants with increased sensitivity to residual androgens after castration and overexpression of AR or downstream effectors via ligand-independent modifications [5]. Despite the current findings of AR signaling in CRPC, AR-independent molecular mechanisms involved in treatment resistance and the paradigm shift of hormone-naive PCa leading to CRPC remain largely unknown. Therefore, a broader understanding of the cellular functioning during ADT/ATT is essential as it might lead to novel or co-treatment strategies to improve efficacy of ATT. Although previous studies have investigated the changes in androgen-responsive genes (ARGs) towards PCa progression, only a few attempts have been taken to explain off-target effects of AR antagonists/anti-androgens that bypass AR signaling and other AR-independent mechanisms [6,7].
In the current study, we aimed to delineate AR antagonists/anti-androgens induced molecular alternations in PCa cells. We used androgen-dependent (AD) cell model (LNCaP) treated with 1st generation non-steroidal anti-androgen: Bicalutamide (BIC) and 2nd generation non-steroidal anti-androgen: Enzalutamide (ENZ). In order to profile the PCa proteomic changes, we have used the sequential window acquisition of all theoretical mass spectra (SWATH-MS) label-free quantification method, and RNA sequencing (RNAseq) analysis was used to observe the changes in PCa transcriptome. Comparative analysis between LNCaP cells and the metastatic, androgen-independent (AI) cell model (C4-2B) represented molecular changes acquired during the paradigm shift of androgen-dependent PCa to androgen-independence. Ultimately, the study highlighted key molecular alternations and pathways that emerged during anti-androgen therapy, supporting treatment escape and PCa progression to CRPC.

Summary of the Proteomic and Transcriptomic Profiling of PCa Cell Lines
To mimic the anti-androgen therapy given for patients in the androgen-dependent stage and find alternations in the PCa cell proteome, LNCaP cell line was obtained and supplemented with androgen/dihydrotestosterone (DHT), BIC and ENZ treatments ( Figure S1, see Section 4). Via SWATH-MS proteomic profiling of treated LNCaP cells, a total of 2628 proteins were identified at 1% False discover rate (FDR). According to the differential expression analysis we identified 318, 215 and 182 significantly differentially expressed proteins (DEPs) for DHT, DHT + BIC and DHT + ENZ treatments, respectively ( Figure 1A). A total of 74 DEPs were found to be shared between all three treatments, while 75 and 51 DEPs were found specific to DHT + BIC and DHT + ENZ treatment, respectively. In addition, 319 and 258 DEPs were detected with BIC and ENZ only treatments, providing treatment specific protein profiles stimulated under direct inhibition of AR signaling ( Figure 1B). To mimic the PCa progression from androgen-dependent to castration-resistant stage, we obtained the C4-2B cell line: a metastatic derivative of the LNCaP cell model. SWATH-MS proteomic profiling of the C4-2B cells has identified 2703 proteins at 1% FDR. Differential expression analysis has reported 159 DEPs compared to LNCaP cells cultured in growth factor supplement conditions (LNCaP-FBS) and 208 DEPs compared to LNCaP cells cultured under hormone/androgen depleted conditions (LNCaP-CSS) ( Figure 1C). Finally, the common feature of DEPs was the predominance of up-regulated protein expression in response to the treatments, particularly in BIC (69.4%) and ENZ (75.0%) only treatments ( Figure 1D). Moreover, a majority of the proteins were found up-regulated in C4-2B cells compared to LNCaP-FBS or LNCaP-CSS cells. (Figure 1D).  To mimic the PCa progression from androgen-dependent to castration-resistant stage, we obtained the C4-2B cell line: a metastatic derivative of the LNCaP cell  model. SWATH-MS proteomic profiling of the C4-2B cells has identified 2703 proteins at  1% FDR. Differential expression analysis has reported 159 DEPs compared to LNCaP cells cultured in growth factor supplement conditions (LNCaP-FBS) and 208 DEPs compared to LNCaP cells cultured under hormone/androgen depleted conditions (LNCaP-CSS) ( Figure 1C). Finally, the common feature of DEPs was the predominance of up-regulated protein expression in response to the treatments, particularly in BIC (69.4%) and ENZ (75.0%) only treatments ( Figure 1D). Moreover, a majority of the proteins were found upregulated in C4-2B cells compared to LNCaP-FBS or LNCaP-CSS cells. (Figure 1D). Both principal component analysis (PCA) and hierarchical clustering analysis demonstrated the distinct clustering between DHT and DHT + BIC/ENZ treatment proteomic profiles (Figure 2A,B). DHT + BIC/ENZ treatments were clustered separately from the BIC/ENZ only treatments, suggesting the differential effect of anti-androgens on PCa cell behaviour, depending on androgen availability (Figure 2A,B). Moreover, a clear separation between LNCaP and C4-2B proteomic profiles was observed, signifying distinct molecular programs characterized between AD and AI PCa phenotypes ( Figure 2C,D). Both principal component analysis (PCA) and hierarchical clustering analysis demonstrated the distinct clustering between DHT and DHT + BIC/ENZ treatment proteomic profiles (Figure 2A,B). DHT + BIC/ENZ treatments were clustered separately from the BIC/ENZ only treatments, suggesting the differential effect of anti-androgens on PCa cell behaviour, depending on androgen availability (Figure 2A,B). Moreover, a clear separation between LNCaP and C4-2B proteomic profiles was observed, signifying distinct molecular programs characterized between AD and AI PCa phenotypes ( Figure 2C,D).
RNAseq analysis has identified a total of 19,549 genes in LNCaP cells. A total of 2556 common proteomic and transcriptomic elements were reported allowing a better integration (~97%) between the two datasets ( Figure 3A). Although, linear regression analysis showed a positive correlation between the two datasets, the high variability (low R 2 ) indicated the limited concordance between transcriptomic and proteomic expressions in PCa cells ( Figure 3B). Differential expression analysis has identified a total of 2377, 433 and 503 significantly differentially expressed genes (DEGs) in DHT, DHT + BIC and DHT + ENZ treatments, respectively ( Figure S2A). A total of 81 DEGs were found to be shared between all three treatments, whereas 137 and 201 DEGs were specific to DHT + BIC and DHT + ENZ treatments, respectively. In contrast to the proteomic data, transcriptomic data showed that majority of the genes found downregulated (70.7%) with the DHT treatment ( Figure S2B). In line with the observations of proteomic profiling, hierarchical clustering has given a clear separation between DHT and DHT + BIC/ENZ transcriptomic profiles ( Figure S2C). RNAseq analysis has identified a total of 19,549 genes in LNCaP cells. A total of 2556 common proteomic and transcriptomic elements were reported allowing a better integration (~97%) between the two datasets ( Figure 3A). Although, linear regression analysis showed a positive correlation between the two datasets, the high variability (low R 2 ) indicated the limited concordance between transcriptomic and proteomic expressions in PCa cells ( Figure 3B). Differential expression analysis has identified a total of 2377, 433 and 503 significantly differentially expressed genes (DEGs) in DHT, DHT + BIC and DHT + ENZ treatments, respectively ( Figure S2A). A total of 81 DEGs were found to be shared between all three treatments, whereas 137 and 201 DEGs were specific to DHT + BIC and DHT + ENZ treatments, respectively. In contrast to the proteomic data, transcriptomic data showed that majority of the genes found downregulated (70.7%) with the DHT treatment ( Figure S2B). In line with the observations of proteomic profiling, hierarchical clustering has given a clear separation between DHT and DHT + BIC/ENZ transcriptomic profiles ( Figure S2C).

Effects on AR Signaling during Anti-Androgen Treatment
Since growth and proliferation of PCa cells largely depends on AR signaling, we characterized the AR signaling activities of LNCaP cells in response to anti-androgen treatments. In line with previous studies, our study reported AR expression was unaltered

Effects on AR Signaling during Anti-Androgen Treatment
Since growth and proliferation of PCa cells largely depends on AR signaling, we characterized the AR signaling activities of LNCaP cells in response to anti-androgen treatments. In line with previous studies, our study reported AR expression was unaltered in response to both anti-androgens, at mRNA and protein levels ( Figure 4A,B) [8,9]. However, AR mediated transcriptional activation of ARGs (MSigDB Hallmark "Androgen Response" gene set) were found to be suppressed in DHT + BIC/ENZ treatments compared to the DHT treatment ( Figure 4C). A similar observation was made at the protein level in which protein expression of four ARGs: KLK3, FKBP5, NDRG1, and UAP1 were found to be supressed with DHT + BIC/ENZ and BIC/ENZ only treatments compared to DHT treatment ( Figure 4D). The greater reduction of the ARG expression with ENZ treatment (with or without DHT stimulation) may indicate stronger AR-antagonist activity of ENZ compared to BIC treatment. Although not statistically significant, we observed downregulation of proteins related homologous recombination (HR): ATM, RAD51, FANCA, CHEK2 and upregulated PARP1 protein expression in response to BIC and ENZ, demonstrating suppressed HR activity and activated PARP DNA repair system during anti-androgen treatment ( Figure S3).  All measurements are normalized to control treatment. Error bars represent the means ± SD, * p < 0.05.

Androgen and Anti-Androgen Mediated Molecular Networks and Biological Processes in PCa Cells
Androgen and anti-androgen mediated disrupted gene/protein networks were identified by the IPA analysis, hinting at adaptive molecular pathways supporting PCa cell survival and proliferation. We focused on the top scoring cancer related protein/gene networks induced by androgen/anti-androgen treatments showing a functional relationship between DEPs/DEGs (Table S1). Under DHT treatment, the analysis showed a significant enrichment in the protein network of RNA post-transcriptional modification (Score = 67, DEPs = 30) and the gene network of cellular movement (Score = 34, DEGs = 32) ( Figure  S4A,B). In addition, DHT + BIC treatment showed an overrepresentation in the protein network of RNA-post transcriptional modifications (Score = 60, DEPs = 30) and the gene network of molecular transport (Score = 27, DEGs = 18) ( Figure S4C,D). Moreover, DHT +

Androgen and Anti-Androgen Mediated Molecular Networks and Biological Processes in PCa Cells
Androgen and anti-androgen mediated disrupted gene/protein networks were identified by the IPA analysis, hinting at adaptive molecular pathways supporting PCa cell survival and proliferation. We focused on the top scoring cancer related protein/gene networks induced by androgen/anti-androgen treatments showing a functional relationship between DEPs/DEGs (Table S1). Under DHT treatment, the analysis showed a significant enrichment in the protein network of RNA post-transcriptional modification (Score = 67, DEPs = 30) and the gene network of cellular movement (Score = 34, DEGs = 32) ( Figure S4A,B). In addition, DHT + BIC treatment showed an overrepresentation in the protein network of RNA-post transcriptional modifications (Score = 60, DEPs = 30) and the gene network of molecular transport (Score = 27, DEGs = 18) ( Figure S4C,D). Moreover, DHT + ENZ treatment demonstrated a significant enrichment in the protein network of protein synthesis and RNA damage and repair (Score = 63, DEPs = 30) and the gene network of post-translational modification (Score = 34, DEGs = 22) ( Figure S4E,F).

Upstream Regulators Dysregulated in Response to the Anti-Androgen Treatment
Next, we performed upstream regulatory analysis to identify likely affected transcription regulators during androgen and anti-androgen treatments (Table S3). Intriguingly, the MYC upstream regulator was found activated in both DHT + BIC (activation Z-score: 3.188, p = 1.14 × 10 −09 ) and DHT + ENZ (activation Z-score: 2.849, p = 8.11 × 10 −06 ) treatments, while it was found inhibited in response to DHT treatment (activation Z-score: −2.678, p = 1.35 × 10 −16 ) (Table S3, Figure 6A,B). Both DHT + BIC/ENZ treatments have induced the upregulation of several MYC targets, mainly ribonucleoproteins (HNRNPs, SNRNP), pre-mRNA splicing factors (PRPs) and DEAD-BOX (DDX) RNA helicases that regulate RNA metabolism. Moreover, MYC oncogenic targets such as DEK and YBX1 proteins were found upregulated in response to both DHT + BIC/ENZ treatments. In addition to MY C, activation of proline and glutamine-rich (PSF/SFPQ) RNA binding/splicing factor upstream regulator was observed in both DHT + BIC (activation Z-score: 2.429, p = 3.65 × 10 −08 ) and DHT + ENZ (activation Z-score: 2.216, p = 9.5 × 10 −06 ) while no change in PSF/SFPQ activation/inhibition was observed in DHT treatment (Table S3, Figure 6C,D). Expression of several PSF/SFPQ target splicing factors, such as HNRNPU, SF3B2 and U2AF1/2 were found upregulated in response to both BIC and ENZ treatments. Overall, this analysis shed light on the differential expression of numerous mRNA-splicing factors regulated by MYC and PSF/SFPQ upstream regulators during anti-androgen treatment. regulate RNA metabolism. Moreover, MYC oncogenic targets such as DEK and YBX1 proteins were found upregulated in response to both DHT + BIC/ENZ treatments. In addition to MY C, activation of proline and glutamine-rich (PSF/SFPQ) RNA binding/splicing factor upstream regulator was observed in both DHT + BIC (activation Z-score: 2.429, p = 3.65 × 10 −08 ) and DHT + ENZ (activation Z-score: 2.216, p = 9.5 × 10 −06 ) while no change in PSF/SFPQ activation/inhibition was observed in DHT treatment (Table S3, Figure 6C,D). Expression of several PSF/SFPQ target splicing factors, such as HNRNPU, SF3B2 and U2AF1/2 were found upregulated in response to both BIC and ENZ treatments. Overall, this analysis shed light on the differential expression of numerous mRNA-splicing factors regulated by MYC and PSF/SFPQ upstream regulators during anti-androgen treatment.

Anti-Androgen Induced Molecular Pathways that Denote the Metastatic Androgen Independent PCa Phenotype
Next, we sought to determine enrichment of treatment induced DEPs in molecular pathways that might provide PCa cells with a growth advantage. DEPs induced by DHT Figure 6. Upstream regulatory analysis for the differentially expressed genes/proteins identified by anti-androgen treatments. Regulatory effect of MYC upstream regulator on target genes in response to (A) DHT + BIC treatment, (B) DHT + ENZ treatment compared to control treatment. Regulatory effect of SFPQ upstream regulator on target genes in response to (C) DHT + BIC treatment, (D) DHT + ENZ treatment compared to control treatment. Colour in the target genes represent their upregulation or downregulation while colour in the upstream regulator represent their predicted activation or inactivation. Colour intensity represents the relative magnitude of change in expression. Edges are coloured to indicate the expected direction of regulation.

Anti-Androgen Induced Molecular Pathways That Denote the Metastatic Androgen Independent PCa Phenotype
Next, we sought to determine enrichment of treatment induced DEPs in molecular pathways that might provide PCa cells with a growth advantage. DEPs induced by DHT + BIC/ENZ treatments and DEPs identified in C4-2B vs. LNCaP comparison were mapped in to KEGG pathways. Interestingly, the analysis revealed that pathways dysregulated in the androgen-independent C4-2B cell model were comparable with those dysregulated in DHT + BIC/ENZ treatments (Table S4).
We observed a significant enrichment of DEPs in genetic information processing pathways, such as (i) spliceosome pathway: HNRNPU/K; YBX1; SRSF2/3/7; DDX5, SNRNP70 ( Figure S6A), (ii) ribosome biogenesis: DKC1 ( Figure S6B), (iii) ribosome: RPL17/19/29; RPS10/17/19 ( Figure S6C), (iv) RNA transport pathway: SRRM1; translation initiation factors (eIFs) ( Figure S6D) and (v) Protein processing factors in ER: SSR3; DAD1 ( Figure S6E). Moreover, proteins related to energy metabolism: oxidative phosphorylation (OXPHOS) were found dysregulated, mainly the proteins involved in NADH dehydrogenase activity: NDUFA2/10; NDUFS2 ( Figure S6F). In addition, proteins involved in cell growth, death and cell communication pathways were found dysregulated, such as (i) MCM proteins related to cell cycle maintenance and (ii) PAK4 associated with cellular focal adhesion ( Figure S6G). Furthermore, proteins that belong to the PI3K-AKT and MAPK signaling pathway, such as MAPK1 were found dysregulated ( Figure S6H,I). A predominant involvement of cancer related molecular alterations was observed, for instance, (i) overrepresentation of AR transcriptional activity reported by the upregulation of AR targets: KLK3 and heat shock protein (HSPs): HSP90AB1 ( Figure S6J) and (ii) enrichment of cancer associated proteoglycans, such as Cortactin (CTTN) ( Figure S6K). Overall, the current analysis highlighted a set of proteins related to 11 cancer related signaling pathways that might aggravate treatment resistance and lead to the insurgence of androgen-independence (Table S4, Figure 7). We observed a significant enrichment of DEPs in genetic information processing pathways, such as (i) spliceosome pathway: HNRNPU/K; YBX1; SRSF2/3/7; DDX5, SNRNP70 ( Figure S6A), (ii) ribosome biogenesis: DKC1 ( Figure S6B), (iii) ribosome: RPL17/19/29; RPS10/17/19 ( Figure S6C), (iv) RNA transport pathway: SRRM1; translation initiation factors (eIFs) ( Figure S6D) and (v) Protein processing factors in ER: SSR3; DAD1 ( Figure S6E). Moreover, proteins related to energy metabolism: oxidative phosphorylation (OXPHOS) were found dysregulated, mainly the proteins involved in NADH dehydrogenase activity: NDUFA2/10; NDUFS2 ( Figure S6F). In addition, proteins involved in cell growth, death and cell communication pathways were found dysregulated, such as (i) MCM proteins related to cell cycle maintenance and (ii) PAK4 associated with cellular focal adhesion ( Figure S6G). Furthermore, proteins that belong to the PI3K-AKT and MAPK signaling pathway, such as MAPK1 were found dysregulated ( Figure S6H,I). A predominant involvement of cancer related molecular alterations was observed, for instance, (i) overrepresentation of AR transcriptional activity reported by the upregulation of AR targets: KLK3 and heat shock protein (HSPs): HSP90AB1 ( Figure S6J) and (ii) enrichment of cancer associated proteoglycans, such as Cortactin (CTTN) ( Figure S6K). Overall, the current analysis highlighted a set of proteins related to 11 cancer related signaling pathways that might aggravate treatment resistance and lead to the insurgence of androgen-independence (Table S4, Figure 7).  Anti-androgen induced molecular pathways that denote the insurgence of metastatic, androgen independent PCa phenotype. The chart represents the key dysregulated proteins common to anti-androgen administered androgen-sensitive LNCaP cells and C4-2B androgen-independent cells identified by the comparative bioinformatics analysis of proteomic and transcriptomic data. Molecules were distributed in to 11 KEGG pathways possibly affected in cancer.

Validation of Gene/Protein Expression Using PCa Clinical Data
We have utilized publicly available PCa transcriptomic and immunohistochemical (IHC) datasets to verify the expression of the two upstream regulators and key oncogenic genes/proteins in clinical samples. Based on PCa patient transcriptomic data retrieved from Oncomine TM database, we observed the marked upregulation of the SFPQ upstream regulator in CRPC tissues compared to localized PCa and normal prostate tissues. Although, the MYC upstream regulator was upregulated in localized PCa compared to normal prostate tissues, no significant change was observed in CRPC tissues. In addition, HNRNPU; DKC1; CTTN; MCM2; MAPK1; PAK4 genes were found to be significantly upregulated in CRPC tissues compared to normal prostate tissues (Figure 8). Based on the IHC data retrieved from the HPA database, both MYC and SFPQ upstream regulators showed a moderately higher protein expression in high-grade PCa tissues compared to the lowgrade PCa tissues. Likewise, HNRNPU, MCM2, MAPK1, CTTN protein targets also demonstrated relatively higher expression in high-grade PCa tissues compared to the low-grade PCa tissues ( Figure 9). Intriguingly, DKC1 showed a high expression in both high-/low-grade tissues, whereas PAK4 was expressed at moderate levels in both high-/low-grade PCa tissues. In summary, the above analysis highlighted on the possible involvement of anti-androgen mediated upstream regulators and oncogenic proteins in PCa pathogenesis and their effect on fueling PCa progression to CRPC.

Validation of Gene/Protein Expression Using PCa Clinical Data
We have utilized publicly available PCa transcriptomic and immunohistochemical (IHC) datasets to verify the expression of the two upstream regulators and key oncogenic genes/proteins in clinical samples. Based on PCa patient transcriptomic data retrieved from Oncomine TM database, we observed the marked upregulation of the SFPQ upstream regulator in CRPC tissues compared to localized PCa and normal prostate tissues. Although, the MYC upstream regulator was upregulated in localized PCa compared to normal prostate tissues, no significant change was observed in CRPC tissues. In addition, HNRNPU; DKC1; CTTN; MCM2; MAPK1; PAK4 genes were found to be significantly upregulated in CRPC tissues compared to normal prostate tissues (Figure 8). Based on the IHC data retrieved from the HPA database, both MYC and SFPQ upstream regulators showed a moderately higher protein expression in high-grade PCa tissues compared to the low-grade PCa tissues. Likewise, HNRNPU, MCM2, MAPK1, CTTN protein targets also demonstrated relatively higher expression in high-grade PCa tissues compared to the low-grade PCa tissues ( Figure 9). Intriguingly, DKC1 showed a high expression in both high-/low-grade tissues, whereas PAK4 was expressed at moderate levels in both high-/low-grade PCa tissues. In summary, the above analysis highlighted on the possible involvement of anti-androgen mediated upstream regulators and oncogenic proteins in PCa pathogenesis and their effect on fueling PCa progression to CRPC.

Discussion
Identifying molecular alterations dysregulated during ATT which likely leads to treatment resistance could offer new opportunities for patient selection and management, preventing disease progression to CRPC. The current study investigated both transcriptomic and proteomic changes induced by BIC and ENZ anti-androgen treatments that might be responsible in developing PCa treatment resistance and androgen-independence.
We increased the depth of proteomic profiling in our cell models via SWATH-MS strategy: an accurate, highly reproducible label-free approach in whole proteome quantifi-cation [11,12]. Two recent studies that used SWATH-MS strategy on LNCaP cells validated the total number of proteins identified in our study [13][14][15]. In addition, we observed a limited correlation between PCa transcriptomic and proteomic datasets, in agreement with previous PCa omics-based research [16,17]. Varied post-transcriptional modifications (PTMs) associated with translation regulation and kinetic differences between protein synthesis and turnover may explain the poor correlation between mRNA and protein abundances in complex biological samples [18].
Although, CRPC is frequently associated with the loss of AR expression, recent studies describe that anti-androgen resistance and PCa progression can be acquired via an AR indifferent state, where cancer cells do not depend on AR signaling regardless of persistent AR expression [8,9]. We observed that anti-androgen treatment blocks AR mediated transcriptional activities without affecting AR expression, which provides evidence for the potential of LNCaP cells to acquire AR-independent adaptive mechanisms via an AR-indifferent state. The strong AR-antagonistic activity shown by ENZ can be explained by its higher affinity for AR in comparison to BIC [19]. The clinical benefit of ENZ over BIC in patients with metastatic CRPC and its superior activity in suppressing AR signaling has been demonstrated in clinical trials [20,21]. Nevertheless, AR signaling can be continued in a fraction of cells due to the potential agonistic functioning of BIC and ENZ, which allow LNCaP cells to maintain their dependency on AR signaling [8,22,23]. This observation is clinically manifested by the anti-androgen withdrawal syndrome and clinical improvement observed upon anti-androgen discontinuation [22]. Moreover, our data demonstrated moderate downregulation of HR gene expression and moderate upregulation of PARP activity in response to both BIC and ENZ treatments. Although expression levels are not statistically significant direction of HR and PARP protein regulation is in line with previous observations that show PARP backup DNA repair system is activated in absence of HR which support the survival of PCa cells during ADT [23,24].
Cancer cell survival depends on the coordinated activation of translation initiation factors and major signaling cascades regulating translation that are implicated in a broad range of cellular processes [25,26]. Our analysis highlighted the overrepresentation of molecular networks associated with transcriptional and translational regulation, which might allow them to adapt to anti-androgen induced cellular stress. This finding was further confirmed by the MYC and PSF/SFPQ mediated transcriptional activation of splicing factors that regulate the expression of genes promoting PCa survival. A previous study confirmed that androgen-treatment in hormone-depleted PCa cells reduced MYC levels [27]. MYC plays an important role in regulating AR and its variant expression through transcriptional regulation coupled with modulation of protein stability [28]. Furthermore, it has been demonstrated that BIC has the potential to reverse androgen-mediated reduction of MYC [29]. Therefore, the transcriptional regulation conferred by MYC activation is considered to be a mechanistic insight into the ability of PCa tumors to overcome stressful conditions induced by anti-androgen therapy [30][31][32]. On the other hand, PSF/SFPQ mRNA splicing factor has been previously reported to be involved with the worse prognosis of PCa as shown by the upregulation of PSF/SFPQ downstream targets in metastatic PCa patients [33]. These splicing factors are responsible for AR overexpression and alternative-splicing (AR-v7) which are key drivers towards androgen-independence [33]. Hence, upstream MYC and PSF/SFPQ regulators and downstream transcription factors have an obvious pharmacological importance as their inhibition could be a viable approach to alleviate anti-androgen resistance [28].
We observed the coordinated activation of pathways in LNCaP cells treated with anti-androgens and their overlap with those dysregulated in C4-2B cell populations. Key oncogenic proteins dysregulated in these pathways represent the molecular repertoire that support treatment escape and progression to CRPC. The highly perturbed environment shaped by the androgen-depleted medium supplemented with anti-androgens has triggered molecular pathways related to gene transcription and translation regulation. This includes the dysregulation of spliceosome components, such as HNRNPs, YBXs, SRSFs, DDXs and SNRNPs, that act as proto-oncogenic splicing regulators involved in the adaptive transcriptional programming of hormone-related cancers, such as breast and PCa [34][35][36][37][38]. HNRNPU and HNRNPK are two splicing factors involved in regulating AR and AR-v7 expression by post-transcriptional mechanisms and it has been previously found that BIC treatment alters AR/HNRNPK interaction and AR activation [35,39,40]. Proteins involved in ribosome biogenesis: DKC1 and ribosomal components: RPL19/29, RPS19/20 have been previously reported for their dysregulation and malignant potential in advanced PCa phenotypes [41,42]. Splicing factor: SRRM1 that is involved in the RNA transport pathway has been recently described for its association with PCa proliferation/migration-rate via modulating AR-v7 levels [39]. ER signal sequence receptor subunit: SSR3 that regulates protein entry into ER plays a major role in protein processing in ER and its dysregulation has been found associated with tumor growth in PCa [40,43,44]. DAD1 protein is involved in protein glycosylation and translocation across the ER that exhibit anti-apoptotic functions [45]. Previous data shows the upregulated expression of DAD1 in AI-PCa cell lines and grade-associated expression in PCa clinical tissues, thus, providing functional insights of DAD1 in therapy resistance and PCa progression [46,47].
Our data showed the dysregulation of OXPHOS, particularly characterized by upregulated NADH dehydrogenase activity which is the rate-limiting step in overall respiration. Overrepresentation of OXPHOS provides PCa cells with a growth advantage by overproduction of ROS to cause more DNA damage and mutations [48,49]. This observation is in line with previous data that propose OXPHOS involvement in androgenindependence [50,51]. Overexpression of minichromosome maintenance proteins (MCM) involved in cell cycle regulation are currently used as proliferation markers to determine tumor growth propensities [52][53][54]. Dysregulation of MCM proteins observed in our findings suggest their role in facilitating tumor cell proliferation during anti-androgen treatment. Upregulation of proteins involved in focal adhesion, such as PAK1 indicate the coordinated activation of cancer cell invasion in response to anti-androgen treatment [55][56][57]. Reiterating previous literature, our data showed the dysregulation of proteins belonging to the PI3K-AKT and MAPK pathway, such as MAPK1. Activation of PI3K-AKT/MAPK signaling promotes PCa cell proliferation and survival, hence combining PI3K-AKT/MAPK inhibitors with anti-androgen therapy were found to be promising in clinical trials [58,59]. Upregulation of KLK3/PSA in C4-2B cells confirms the activation of AR-signaling during PCa progression to CRPC via various AR bypassing mechanisms [60][61][62]. Moreover, overexpression of PSA can induce a feedback regulatory loop by which AR expression is maintained during CRPC progression [63]. On the other hand, HSPs, which play a pivotal role in AR stabilization during AR signaling, were found upregulated in response to antiandrogen treatment and in CRPC cells. Previous literature suggests that BIC agonistic AR activation and its downstream interaction with coactivators and transcriptional factors may rely on the stabilization provided by HSP90s [64]. Hence, using HSP inhibitors (e.g., geldanamycin) in combination with anti-androgens would increase the efficiency of PCa therapy. Proteoglycans are key constituents of extracellular matrix capable of modulating cancer cell invasion [65]. Therefore, it signifies the contribution of upregulated proteoglycans, such as CTTN towards PCa cell survival and invasion during anti-androgen therapy.
The limited correlation between the proteomic and transcriptomic data might limit the identification of extensive functional changes in PCa cells in response to stress conditions. Thus, understanding the limitations between different measurement techniques is critical before interpreting the results integrating datasets from multiple sources [66]. SWATH-MS DIA data acquisition is the easiest compared to DDA acquisition methods, yet difficult in terms of data analysis as it requires large informatic resources and sophisticated software tools [67]. Moreover, protein identification in SWATH-MS strategy is largely limited by the spectral library composition; therefore, generating a larger-scale spectral library is essential using both cell line and tissue DDA data. Protein targets identified in the current study are highly likely to be biomarkers for CRPC as well as targets for novel therapies. Nevertheless, our findings require further validation by applying different AD/AI cell lines and molecule/pathway inhibitors, before testing them using in vivo mouse models. Based on the resistance observed for BIC and ENZ, novel anti-androgens (e.g., Apalutamide, Darolutamide) have already made their way through preclinical and clinical studies to treat both metastatic and non-metastatic CRPC [68]. Therefore, the next wave of scientific focus should be identifying possible adaptive mechanisms induced by these novel anti-androgens to alter the course of PCa treatment resistance.

SWATH-MS/MS Data Analysis
SWATH-MS/MS raw data was processed against the above created spectral library. Extracting peak areas and scoring was performed using the PeakView ® SWATH micro app (Version 2.1, AB SCIEX). Shared peptides were excluded while importing the spectral library. Minimum number of 3 transitions (fragment ions) and 3 precursor ions with >90% confidence were used for SWATH-MS/MS data processing. Peak areas were extracted for peptides with <1% FDR. Peak area of fragment ions was used for peptide quantification, whereas the mean value of the peptide quantities was used to quantify proteins.

RNA Isolation
Total RNA was extracted from PCa cells using RNAeasy Mini Kit according to the standard protocol (Qiagen, Hilden, Germany). RNA concentration and purity were measured using NanoDrop TM 1000 (Thermo Scientific).

Real-Time Quantitative PCR (qRT-PCR) Analysis
RNA (1 µg) was reverse transcribed to cDNA using SuperScript III Reverse Transcriptase (Invitrogen, Carlsbad, CA, USA). Quantitative RT-PCR was performed in MicroAmp ® Optical 384-Well Reaction Plate with barcode (Applied Biosystems, Foster City, CA, USA) using the ViiA7 Real-Time PCR system (Applied Biosystems). Each reaction contained 1X final concentration of SYBR Green PCR Master Mix 2X (Life Technologies), 50 nM forward and reverse primer, 2.0 µL of diluted cDNA (1:5) and nuclease-free water to a final volume of 8 µL. The cycling parameters were 95 • C for 10 min, 40 cycles of 95 • C for 15 s and 60 • C for 1 min followed by a dissociation step. Relative fold expression was performed by the comparative CT (∆∆CT) method.

RNAseq Analysis
Total RNA extracted from the androgen and anti-androgen treated LNCaP cell lines used for the RNAseq analysis performed through the Australian Genome Research Facil-ity (AGRF). Three independent RNAs for each treatment were pooled and used for the ribosomal depletion. Paired-end sequencing (100 bp resolution) has performed on the Illumina HiSeq platform using the Illumina TruSeq strand-specific protocol (Life Technologies). RNAseq was carried out at 100 bp resolution, paired end, and strand specific. Raw data was analysed by mapping RNAseq reads using Tophat2 (hg19 assembly) [70], and logFC and p values for the gene expression analysis were determined using the edgeR program [71].

Statistical Analysis
Exported SWATH-MS/MS peak area values were statistically analysed using the MetaboAnalyst 4.0 online tool [72]. First, peak area values of each protein, obtained from the three replicates were averaged, normalized (based on the sum) and log 2 transformed. PCA analysis was performed for the first three principal components using log 2 transformed normalized areas. For the Heatmap visualization, normalized SWATH areas were clustered based on Ward algorithm and row/column orders were set by applying a hierarchical clustering based on Euclidean distance. Fold change (FC) expression between each treatment group and the log 2 transformation of FC (Log 2 FC) was calculated. A threshold of |log 2 FC| ≥ 1 was used to identify DEPs, whereas |logFC| ≥ 1.5 threshold was used identify DEGs. p < 0.05 was considered significant.

Network Analysis
SWATH-MS/MS Proteomic datasets and RNAseq datasets were analysed using Ingenuity Pathway Analysis (IPA) (Qiagen, http://www.ingenuity.com) and list of molecular networks were generated based on comparing imported datasets and Ingenuity ® Knowledge Base. A score (p-score = −log10 (p-value)) according to the fit of the set of supplied protein/genes and molecular networks stored in the Ingenuity Knowledge Base were generated. Both direct and indirect relationships were considered and only DEPs/DEGs with a p < 0.05 were considered for the analysis.

Upstream Regulatory Analysis
IPA upstream regulator analysis was used to predict the activation/inhibition of upstream transcriptional regulators from the datasets based on the literature and compiled in the Ingenuity ® Knowledge Base. A p-value was computed based on significant overlap between known targets regulated by the transcriptional regulator and the molecules in the dataset. The activation Z-score algorithm was used, to make predictions on the upstream regulator. The expression of the upstream regulator itself and only DEPs/DEGs with a p < 0.05 were considered for the analysis.

Gene Set Enrichment Analysis
GSEA was used to determine over-represented GO-BPs enriched by DEPs/DEGs identified in androgen/anti-androgen treated LNCaP cell models [10]. GSEA converted the submitted proteins/genes into Entrez genes on the MSigDB v7.1 [73][74][75]. The analysis evaluated the overlap of provided DEPs/DEGs with human MSigDB GO-BP gene sets and estimated the statistical significance at FDR/q < 0.05.

KEGG Mapper Analysis
KEGG mapper tool was used to map DEPs/DEGs in to KEGG pathways that are dysregulated in AD LNCaP cell model in response to anti-androgen treatment and AI C4-2B cell model. The tool converted the submitted proteins/genes into KEGG GENES ID and matched with the human specific KEGG pathway maps enabling interpretation of cellular functions and other high-level features [76,77]. Pathway lists and categories were downloaded from the KEGG database and different colour codes were used to demonstrate treatment specific and common DEPs/DEGs identified [78].

Validating the Expression of Differentially Regulated Genes/Protein
Differentially expressed protein/genes identified in PCa cell models were revaluated by measuring their expression in Grasso PCa clinical datasets [79] using the online Oncomine TM database (https://www.oncomine.org) [80]. Protein expression in PCa clinical/tissue specimens were analyzed by the PCa immunohistochemistry (IHC) datasets using the Human Protein Atlas (HPA) database (http://www.proteinatlas.org) [81,82].

Conclusions
The current study tested the hypothesis that molecules dysregulated in response to ATT enable the overrepresentation of adaptive molecular pathways supporting PCa treatment escape. Using SWATH-MS strategy, we highlighted anti-androgen induced key proteomic signatures that are involved in transcription/translation regulation, energy metabolism, cell growth and death, cell communication and signaling cascades that might potentially evolve to treatment resistance and androgen-independence in PCa cells. Following further validation, these proteomic signatures and associated pathways can be co-targeted to improve the efficacy of existing treatments and develop novel treatment strategies to attenuate disease progression to CRPC.

Supplementary Materials:
The following are available online at https://www.mdpi.com/2072-669 4/13/4/715/s1. Figure S1: Workflow of the proteomic and transcriptomic analysis of androgen and anti-androgen treated of LNCaP cells and C4-2B cells, Figure S2: Summary of the RNAseq analysis of androgen and anti-androgen treatment of LNCaP cells, Figure S3: Homologous recombination (HR) related protein and PARP protein expression in androgen and anti-androgen treated LNCaP cells, Figure S4: Top scoring gene networks characterized by androgen and anti-androgen treatments, Figure S5: GO-biological processes overrepresented by the genes up-/down-regulated in response to androgen and anti-androgen treatment, Figure S6: KEGG pathways differentially regulated during anti-androgen treatment, Table S1: Dysregulated molecular networks in response to androgen and anti-androgen treatment, Table S2: GO-biological processes affected by androgen and anti-androgen treatment, Table S3: Dysregulated upstream regulators in response to androgen and anti-androgen treatment, Table S4: Top differentially regulated KEGG pathways and associated proteins/genes in LNCaP cells treated with anti-androgens and C4-2B cells.

Conflicts of Interest:
The authors declare no conflict of interest.