The complexity and multigenic nature of cancer has necessitated various genomewide studies to achieve a systems-level understanding of the key genetic mediators involved in prostate cancer ⁽¹⁾ Most diseases are due to the collapse of cellular processes together with interaction networks ⁽²⁾). Therefore, exploring the biological network for complex diseases provides an understanding of the functional alterations in chronic diseases ⁽³⁾.

Network-based approaches contain many clinical applications to explore human diseases systematically. A better understanding of the effects of cellular networks on disease progression may lead to the identification of disease genes which, in turn, may offer better targets for drug development ⁽⁴⁾. One focal point in cancer analysis is the reconstruction of Gene Regulatory Networks (GRN) ⁽⁵⁾. However, cancer progression is a dynamic process with multiple stages; so, reconstruction of one static GRN may not be informative enough for the inference. Instead, reconstructing stage-specific GRNs during cancer progression and then comparing these GRNs would be beneficial to characterize the main genes and interactions involved in cancer progression.

The availability of genomewide gene expression data has helped develop various state-of-art GRN reconstruction methods ^(5–7). These methods seek to identify putative gene regulatory interactions by assuming that alterations in the expression level of a regulator (such as a transcription factor) have a direct effect on the cognate regulated genes.

Empirical evidence of extensive GRN rewiring during cancer progression, along with the availability of GRN reverse engineering approaches, have inspired us to conduct a systematic investigation to characterize the topological changes that occur in a prostate cell's GRN during cancer progression. Therefore, in this study, an attempt was made to reveal candidate disease-associated genes and biomarkers for prostate cancer progression by integrative gene expression profiling and network analysis at a systematic level. In this way, four stage-specific GRNs were reconstructed based on a comprehensive prostate cancer gene expression dataset containing 171 different samples monitoring gene expression at different disease phases. Topological comparison of these four GRNs based on hub type variations recapitulates the previous findings about extensive GRN rewiring ⁽⁸⁾. Through sub-network analysis, it is possible to identify significant genes which were supposed to change their hub type with highly relevant to specific phases of prostate cancer.

Enormous efforts have been made to identify biomarkers for various cancers by the analysis of different transcriptome data ^(9–11). Moreover, there were similar studies for analysis of sub-networks or hub genes which had been helpful for the understanding of the metastasis of cancer at the molecular level ⁽¹²⁾. Nonetheless, there are still few studies on identification of prostate cancer biomarkers for disease progression ⁽¹³⁾. Therefore, in this study, a new integrative network-based approach was developed to detect party hubs and date hubs based on Degree and BN algorithms during cancer progression.

Our analysis led to identification of 38 important genes putatively involved in prostate cancer. Through extensive literature search, experimental evidences revealed the role of 76.3% of candidate genes in prostate cancer (Table 1). This level of experimental confirmation reflects the high accuracy of the proposed approach.

Our study hereby demonstrates a useful approach for analysis of prostate cancer at the systematic level. For the genomewide investigations, this will be a fundamental attempt for future development of the translational medical informatics, which lead to better patient diagnostics with high-throughput data through systems biology ⁽¹⁴⁾.

To reconstruct cell stage specific GRNs, an attempt was made to focus on the available comprehensive transcriptome dataset, originally published ⁽²⁰⁾. This dataset was generated by sampling from four different types of prostate tissues including normal cells (Nor-mal), normal cells adjacent to cancer cells (Adjacent), primary tumor cells (Tumor), and metastatic cells (Metastasis). In our approach, genes were analyzed and prioritized based on the transcriptome data. Hence, it was possible to make reliable predictions only for genes with altered expression level across various conditions. To focus on these genes, only up/down-regulated genes were considered (fold change ≥2 and p < 0.05) in at least one cell stage (978 genes). Also, enrichment of known cancer genes was checked among this set by using a previously curated list of 555 high confidence cancer genes, originally published ⁽³¹⁾. 100 genes were collected and identified as mediators in metastatic prostate cancer from ⁽³²⁾, and 276 genes were added and annotated as either a cancer pathway or prostate cancer gene in the KEGG database. It was found that the cancerrelated genes were about 1.55-fold (hypergeometric two-tailed test, p = 2.52E-6) and the prostate cancerrelated genes were about 2.22-fold enriched (hypergeometric two-tailed test, p = 6.07E-6) in our selected gene set which were fluctuated during prostate cancer.

To identify master regulators and their associated circuits governing cell-specific behavior in each state, the GRNs of prostate cells were compared in different stages with each other. Because the CLR algorithm merely relies on the similarity of expression patterns to infer interactions, constructed networks in this step contain both regulatory interactions (interactions between regulated genes and their putative regulators) as well as interactions between co-regulated genes (non-regulatory interactions). Hereafter, these networks are called coexpression networks. To extract gene regulatory interactions from these networks, only interactions involved at least one human transcription factor were considered and a list of them were extracted ^{(33, 34)}. These networks are referred to as GRNs.

To predict the key genes in the prostate cancer, an attempt was made to find the stage-specific coexpression networks of prostate cancer for high connectivity (hub) or bottleneck genes ^{(22, 23)}. Hub and bottleneck properties are considered important centrality indices because they are major intersections between clusters in the network and if they are removed from a network, the network will be disrupted ⁽²⁴⁾. Han et al suggested the existence of two types of protein hubs in the protein-protein interaction networks, namely party hubs and date hubs ⁽³⁰⁾. Although both interact with many proteins, the difference is that party hubs are proteins that interact with many other proteins simultaneously, whereas date hubs interact with their partners asynchronously ⁽³⁰⁾. By definition, the bottleneck proteins are responsible for the interconnection of clusters in the network, and thus bottlenecks with high degrees are most likely to be date hubs which contain groups of genes that assist in presenting common functions ^{(24, 35)}. The obtained results recapitulate previous findings in which some active sub-networks contained regulatory interactions were supplanted by new interactions which changed their degrees during different conditions ⁽³⁶⁾.

The result also reflected the high level of rewiring of gene regulatory circuits during cancer progression, as suggested elsewhere ⁽⁸⁾. As shown in table 3, for example, more than 2-fold decrease in the number of interactions for KLF6 was observed which controlled cell cycle progression and apoptosis. Indeed, experimental data suggest that KLF6 is inactivated in many cancers such as prostate, ovary and colon ^{(37, 38)}. On the other hand, consistent with more than 2-fold increase in the number of interactions (from 14 in normal stage to 33 in the metastasis stage) for FOXC1 (Table 3), it was indicated that this gene is linked to androgen dependent growth of prostate cancer ⁽³⁹⁾.

Our result led to identification of 38 transcription factors which were bottleneck and changed their interaction during cancer progression. Although the functional role of some famous transcription factors such as AR, SMAD3 and VDR are well known as genes linked to prostate cancer ^(40–42), the 9 transcription factors (CAMTA1, ISL1, MNX1, NHLH2, NKX2-2, STAT2, ZNF146, ZNF205, and ZNF529) are new candidates that may have critical roles in prostate cancer based on topological significance and regulatory changes during cancer progression. Among the remaining 9 transcription factors, 5 of them were associated with other cancer types. ZNF146, CAMTA1, NKX2-2, MNX1, and ISL1 are most prominent in colorectal cancer, neuroblastoma, Ewing's sarcoma, leukemia and breast cancer, and bladder cancers, respectively ^(43–48). No evidence could be found to show the relationship between the 4 remaining transcription factors and any type of cancer.

In this paper, an accurate network-based framework for the analysis of transcriptome data was presented. The analysis of prostate state specific GRNs revealed 38 transcription factors which are critically important for prostate cancer progression. Also, 14 transcription factors were identified to be linked putatively to prostate cancer metastasis stage, so they would be used as key factors for future research in the field of cancer studies. Additionally, experimental evidences revealed the role of 29 of candidate transcription factors in prostate cancer.

The low number of predictions and high degree of overlap with previously known events in the prostate cancer demonstrate the high efficiency of our approach. In addition, the low number of predicted gene sets makes it easy to design follow up experiments to validate the results. In this study, it is believed that the results may provide critical information to gain better understanding of networks dynamics in the cell through complex diseases such as cancer.