 Research
 Open Access
 Published:
Towards targeted combinatorial therapy design for the treatment of castrationresistant prostate cancer
BMC Bioinformatics volume 18, Article number: 134 (2017)
Abstract
Background
Prostate cancer is one of the most prevalent cancers in males in the United States and amongst the leading causes of cancer related deaths. A particularly virulent form of this disease is castrationresistant prostate cancer (CRPC), where patients no longer respond to medical or surgical castration. CRPC is a complex, multifaceted and heterogeneous malady with limited standard treatment options.
Results
The growth and progression of prostate cancer is a complicated process that involves multiple pathways. The signaling network comprising the integral constituents of the signature pathways involved in the development and progression of prostate cancer is modeled as a combinatorial circuit. The failures in the gene regulatory network that lead to cancer are abstracted as faults in the equivalent circuit and the Boolean circuit model is then used to design therapies tailored to counteract the effect of each molecular abnormality and to propose potentially efficacious combinatorial therapy regimens. Furthermore, stochastic computational modeling is utilized to identify potentially vulnerable components in the network that may serve as viable candidates for drug development.
Conclusion
The results presented herein can aid in the design of scientifically wellgrounded targeted therapies that can be employed for the treatment of prostate cancer patients.
Background
Prostate cancer is the most common noncutaneous male malignancy and one of the leading causes of cancer mortality in the western world [1]. The growth and progression of prostate cancer is stimulated by androgens [2]. Androgens are male sex steroid hormones that are responsible for the development of male characteristics. Testosterone is the most important androgen in men. The effects of androgens are mediated through the androgen receptor (AR) [3]. The androgen receptor is a nuclear receptor, which is activated in response to the binding of androgens. Upon activation, it mediates transcription of target genes that modulate growth and differentiation of prostate epithelial cells. In malignant prostate cells, androgen signaling is deregulated and the homeostatic balance between the rate of cell proliferation and programmed cell death is lost. As prostate cancer relies on androgens for growth, the main line of treatment focuses on abrogating the action of androgens. Androgen deprivation therapy (ADT) in the form of surgical or medical castration is the cornerstone of treatment for prostate cancer [4]. Initially, androgen ablation induces significant regression of the tumor. However, the response to ADT is temporary and prostate cancer invariably stops responding to this treatment regimen, leading to a clinical condition that is known as hormonerefractory prostate cancer, androgenindependent prostate cancer or castrationresistant prostate cancer (CRPC). CRPC is a more aggressive and typically lethal phenotype where the tumor continues to grow in spite of the very low levels (<50 ng/ml) of circulating serum testosterone. Standard treatment options are limited and palliative docetaxelbased chemotherapy is generally used for patients who have become refractory to hormone treatment. However, median survival time for patients following firstline chemotherapeutic treatment is just eighteen to twentyfour months [5]. There is therefore a clear rationale for advances in alternative therapeutics in order to evolve and expand the landscape of treatment options for malignant forms of prostate cancer that recur after abatement.
Over recent years, there has been a significant effort towards furthering our understanding of the molecular mechanisms underpinning tumor development, growth and progression. It is now appreciated that in spite of castrate levels of androgens, the cancer cells are able to maintain persistent androgen receptor signaling through a variety of contributory mechanisms including AR gene amplification that results in overexpression of AR, gainoffunction mutations in AR which enable promiscuous activation of the receptor through other steroids or even in the absence of ligand binding, changes in AR coactivators and the expression of AR splice variants [6]. This compensatory response allows cancer cells to survive in a low testosterone environment and the reactivated AR signaling axis continues to play a role after neoplastic transformation. Additionally, certain androgenindependent cellular signaling pathways that promote proliferation and inhibit apoptosis, have been critically implicated as drivers of continued progression of prostate cancer. Hence, accumulating evidence indicates that the growth and progression of prostate cancer is a complicated process that involves interaction between multiple pathways. Advances in our knowledge of the biology of prostate cancer has led to the development of a number of novel therapies designed to target signaling pathways involved in disease progression. With the exception of certain androgen synthesis and AR signaling antagonists that have received regulatory approval, these advanced agents are under various stages of clinical trials [7].
Castrationresistant prostate cancer is a complex malady. Given the inherent complexity of the CRPC signaling cascade, there is no one dominant molecular driver across all tumors and hence no single drug can act as a “magic bullet” by being uniformly effective for treating the malignancy [8, 9]. At best, limited benefit will be derived from targeting a single molecule. Rational combinations of signalmodulating therapeutic agents have higher likelihood of yielding better outcomes. While there are several drugs being tested on cell lines, most of these studies focus on a single pharmaceutical agent and very few of those experiments involve trying out drug combinations. Furthermore, prostate cancer is a markedly heterogeneous disease, with different tumors varying in their composition and makeup. In other words, different tumors will harbor different malfunctions in the signaling pathways. Thus, tailored targeted therapies based on individual tumor characteristics are required to maximize the potential benefits from treatment.
Mathematical and computational modeling plays a pivotal role in systems biology in elucidating biological insights from largescale biomolecular signaling networks that are not amenable to straightforward intuitive interpretation. A diverse array of formalisms have been proposed in this domain as suitable representations for complex multicomponent networks such as cellular signaling pathways [10]. Amongst these frameworks, Boolean network models [11, 12] have emerged as an extremely useful parameterfree approach to capture the qualitative behavior of extensive genetic networks wherein knowledge of kinetic parameters is scarce. Boolean logic models have been successfully applied to study biological signaling networks and cellular processes [13, 14], for instance the cell cycle [15], apoptosis [16], the T cell survival network [17], hypoxia stress response pathways [18] and the gene regulatory network regulating cortical area development [19]. In this paper, we use Boolean logic modeling of the key signaling pathways implicated in the development and progression of prostate cancer to simultaneously test various combinations of agents for their efficacy in attenuating cancer growth and design targeted therapies for the management of the disease. In addition, we attempt to delineate components in the signaling network that can be pharmacologically manipulated to therapeutic advantage.
Methods
Prostate cancer signal transduction network
Cellular processes such as growth and division are regulated by an interconnected network of molecules referred to as signaling pathways. Key cellular signal transduction pathways known to play a major role in cell survival, growth, differentiation and the development of castrationresistance in prostate cancer are the Androgen Receptor (AR), PI3K/AKT/mTOR and MitogenActivated Protein Kinase (MAPK) pathways. The aberrant behavior of prostate cancer cells is characterized by dysfunction in these selective oncogenic signaling pathways promoting malignant characteristics. These pathways play a role in a diverse range of essential physiological cellular processes such as differentiation, survival, proliferation, protein synthesis and metabolism. Malfunctions in these pathways are common in prostate cancer malignancies. For example, approximately 70% of advanced prostate cancers have genomic alterations in the PI3K/AKT/mTOR pathway [20]. These three pathways are the most frequently overactivated pathways increasing survival of cancer cells and promoting cancer progression [21]. A schematic representation of these pathways is shown in Fig. 1 [22–24]. The pharmacologic agents depicted in red boxes in the figure are highly specific pathway inhibitors. These reagents modulate growthfactor receptors and the downstream pathways abnormally activated in CRPC by targeting with great specificity certain signaling nodes in the network.
Boolean modeling of prostate cancer signaling
In the context of methodologies that are applied to model cellular signal transduction networks, Boolean networks are probably the simplest where the state of each node in the network is either active (on) or inactive (off). In a Boolean network, the nodes are the genes and the edges represent the interaction amongst the genes. Since the molecules in a generegulatorynetwork (GRN) exhibit switchlike behavior, genes may be regarded as binary devices where a gene can be considered to be active if it is being transcribed and inactive if it is not. Moreover, the relationships amongst the genes may be represented by means of logical functions. Thus, a GRN is amenable to such a representation. The Boolean formalism is analogous to a digital circuit where logic gates can be used to represent the regulatory relationships amongst the nodes and the activation level of the nodes is indicated by binary logic. The biological interactions amongst the various nodes (genes) represented in the gene regulatory network of Fig. 1 can therefore be translated to an equivalent Boolean circuit [25]. Let’s say either gene X or Y can activate a third gene Z, then we can model this component of the signaling network with an OR gate with two inputs, namely X and Y and with output Z. Thus, the signaling network of Fig. 1 can be mapped to the combinational circuit shown in Fig. 2. This digital logic circuit represents our multiinput multioutput (MIMO) systems model of the prostate cancer signaling transduction network.
Cancer is a disease of abnormal cell signaling caused by a breakdown in the normal signaling pathways leading to the loss of cell cycle control and uncontrolled cell proliferation. These abnormalities in the signaling network can be represented as stuckat faults [26]. A stuckat fault is said to occur when a line in the network is permanently set to a fixed value of one (stuckatone fault) or zero (stuckatzero fault) with the result that the state of the line is stuck at the faulty value and no longer depends on the state of the signaling network upstream that drives that line i.e. the faulty line has a constant (1/0) value independent of other signal values in the circuit. A stuckatfault can occur either at the input or output of a gate. An example of a stuckatfault is given in Fig. 3. Suppose the input vector is <abcd>= 1100. In this case, the output is 0. However, if there is a stuckatone fault at the output of the NAND gate with the same input vector as before, the output of the faulty circuit is one instead of zero. This notion of stuckatfaults has immediate biological relevance: on account of mutations or other structural abnormalities, a gene might become dysfunctional and hence stuck at a particular state irrespective of the signals that it is receiving from surrounding genes [27]. These biological defects can be abstracted as stuckat faults. For instance, as discussed earlier, a diverse array of mechanisms engender persistent AR signaling in CRPC even with castrate serum levels of androgen. This constitutive (permanent) activation of the androgen receptor where the receptor remains active i.e. continues to signal downstream even in the absence of androgens can be represented as a stuckat1 fault. By the same token, the inactivation in cancer of a tumor suppressor, which acts as a molecular brake on cell growth in a normal cell, can be represented as a stuckat0 fault. From our Boolean circuit model, we can explicitly enumerate the different locations where a fault can occur. These fault locations are numbered in Fig. 2 with the stuckat0 and stuckat1 faults in red and black numerals respectively. There is a total number of 24 possible fault locations.
The objective is to counteract the effect of these faults by targeted drug intervention, so we incorporate the drugs in our model. The drug intervention points are illustrated in Fig. 2 which are the locations of the molecules that these prostate cancer drugs are known to target. Since the drugs inhibit the activity of their target i.e. the main mechanism of action of the anticancer drugs is to cut off downstream signaling, their action is incorporated in our model as an inverted input to an AND gate with the result that whenever the drug is applied, the gene that it targets is turned off.
Simulation for fault mitigation with drug intervention
We can now use our Boolean model to test different combination therapies in terms of their efficacy in mitigating the effects of the faults. For each fault, we would like to intervene with the best possible drug combination i.e. we want to determine which set of drugs would be most effective in attempting to nullify the effect of that fault, thereby providing us with a targeted therapy based on the tumor signature. Define, the input vector as follows:
The first four components of this vector are growth factors, which are external signals that stimulate a cell to grow and replicate. The next two input components, namely PTEN and NKX3.1 are tumor suppressors which act as molecular brakes on cell division. The last input vector component consists of the external hormones that stimulate the AR pathway in a normal prostate cell. The input vector is set to be [0000110]. This corresponds to all the external signals that stimulate cell growth being absent and the molecular brakes being active i.e. this input vector corresponds to a nonproliferative input which produces a nonproliferative output in the faultfree case. The output vector is defined to be:
The output vector consists of key markers of cell growth and proliferation in prostate cancer. In the faultfree scenario, a nonproliferative input to the regulatory network should produce a nonproliferative output characterized by the allzero vector. However, faults in the network will produce a nonzero (proliferative) output even when the input is nonproliferative. The objective is to drive the faulty network’s output as close as possible to that of the faultfree circuit i.e. towards the allzero vector through targeted drug intervention. Define, the drug vector as:
Each component of the drug vector is one if the corresponding drug is applied and is zero otherwise i.e. the i^{th} bit of the drug vector is one if the drug is selected and zero if it is not. Thus, for example, the drug vector [0010010] represents the combination of AZD6244 and Temsirolimus. Since, the total number of drugs is seven, the number of possible drug combinations is 128. The objective is to determine the best possible therapy for each fault. Each fault represents a different molecular abnormality and hence a tumor with a different profile.
For each of the faults, the problem is to find the drug selection that can rectify the fault i.e. change the faulty output to the correct output. If that is not possible, the best drug vector will drive the output as close as possible to the faultfree output. A simple metric that can be used as a distance measure to determine how far the output vector is from the faultfree vector is Hamming distance. Faults that produce an output vector with a greater Hamming distance from the correct output have more of the proliferative genes active and presumably a greater proliferative effect. Since the correct output is the allzero vector, the Hamming distance of the output vector from the correct output is simply the Hamming weight of the output vector (for binary vectors Hamming weight is equivalent to the L _{1}norm). For each fault, we determine the output under every possible drug vector. The best therapy for that fault is the drug vector that produces the output with the smallest Hamming weight. In addition, since the drugs have deleterious sideeffects, we would like to choose a drug combination with the fewest number of drugs. Thus, the best targeted therapy for each of the cancerinducing faults is the one that under the presence of the fault, produces the best output with the smallest Hamming weight with the minimal number of drugs.
To determine the best combination therapy across all faults, for each drug combination we determine the sum of the Hamming weights of the output vector across all possible combinations of faults and choose the drug combination that yields the smallest total. In order to keep the computation tractable, we restrict the number of possible faults in any fault combination to be no more than three i.e. up to three genes can be faulty simultaneously. We constrain the cardinality of the drug vector to be less than or equal to three, in essence limiting the number of drugs in the combination to three since on account of the harmful sideeffects of the drugs, administering four or more cancer drugs simultaneously might not be prudent.
Let us formalize the qualitative description above of the selection of best therapy for each fault and that of the overall optimal drug vector. For the Boolean network (BN) of Fig. 2, let N,M and P be the total number of primary inputs, primary outputs and fault locations respectively, then N=7, M=6 and P=24. Let \(\mathbf {x} \in \mathcal {X}\) and \(\mathbf {z} \in \mathcal {Z}\) be the input and output vectors respectively where \(\mathcal {X}\) and \(\mathcal {Z}\) represent the space of all binary vectors of dimensions N and M respectively. Let x ^{∗}=[ 0,0,0,0,1,1,0] be the input vector corresponding to the nonproliferative input.
Let D represent the total number of drug combinations (vectors) with no more than three drugs in any combination, then \(D=\sum \limits _{k=0}^{3}\binom {7}{k}\). Denote each drug vector in the drug space as d _{ i }with i=0,…,D−1 (d _{0} is the allzero drug vector meaning no drug is applied). Let \(\mathcal {D}\) be this space of drug vectors.
Let C be the total number of fault combinations with no more than three faults in any combination, then \(C=\sum \limits _{k=0}^{3}\binom {P}{k}\). Assign each fault combination in the fault space a label f _{ j } with j=0,…,C−1 (f _{0} represents the faultfree case). Let \(\mathcal {F}\) be this set of faults.
Let ψ denote the mapping from a given input vector, drug combination and fault combination to an output vector: \(\mathbf {x} \in \mathcal {X}, \mathbf {d} \in \mathcal {D}, f \in \mathcal {F} \xrightarrow {\boldsymbol {\psi }} \mathbf {z} \in \mathcal {Z}\) i.e. ψ represents the output of the BN for a given input x when a drug combination d is applied under fault scenario f. Let ψ _{ i } be the i^{th} component of this Mdimensional vector ψ.
The best drug vector d _{ i }, i∈{0,1,…,D−1} for each single fault f _{ j }, j∈{1,2,…,P} is the vector of smallest Hamming weight that minimizes ∥ψ(x ^{∗},d _{ i },f _{ j })∥_{1}.
The optimal drug combination across all faults is:
\(\mathbf {d_{i}^{*}}\) is determined by exhaustive enumeration by explicitly searching for the drug combination that for a nonproliferative input, minimizes the sum of Hamming weights (L _{1}norms) of the output vector across all possible combinations of faults.
Node vulnerability assessment
In electronic circuits, reliability refers to the probability of a circuit functioning as intended i.e. producing the correct output. Reliability assessment is used to determine the vulnerability of a circuit to faults. A number of different techniques have been proposed for reliability analysis in digital circuits [28]. Recently, in [29] a scalable, efficient and accurate simulationbased framework based on stochastic computations was introduced for logic circuit reliability evaluation. In biological systems, dysfunctions in nodes in the signaling network cause deviation from normative behavior. Reliability assessment methodologies can be leveraged on Boolean network models of pathways to determine the vulnerability of the network to the dysfunction of each node [30, 31]. In this section we conduct a stochastic logic based vulnerability analysis of the prostate cancer signal transduction network in order to discover the most vulnerable nodes thereby allowing us to prioritize such segments in the network whose perturbation has the greatest potential to yield the most clinical benefit.
In stochastic logic, signal probabilities are encoded in random binary bit streams (the signal probability of a node corresponds to the likelihood of that node having logic value one). For example, the binary sequence 0110010100 of length ten encodes the probability 0.4 since the proportion of ones in this sequence is \(\textstyle \frac {4}{10}\). In practice, the length of the stochastic sequences typically used is much larger. Since the biological literature is devoid of precise ligand binding probabilities, each primary input is assumed equally likely to be 0 or 1 i.e. all primary input signal probabilities are taken to be 0.5.
Stochastic logic often makes use of Bernoulli sequences for the random binary streams where each bit in the stream is generated independently from a Bernoulli random variable with a probability of one equal to p. The use of probabilistic sequences inevitably introduces stochastic fluctuations which implies that the result produced is nondeterministic. These fluctuations can be significantly reduced by representing the initial input probabilities by nonBernoulli sequences [32] defined as random permutations of sequences containing a fixed number of ones and zeros. For a given probability p and sequence length L, a nonBernoulli sequence contains a fixed number pL of ones, with the positions of the ones determined by a random permutation. Thus, for example, to represent the probability 0.5 by a nonBernoulli stream of length 10, we could randomly permute the sequence 1111100000 which has five ones (instead of generating each bit from a Bernoulli random variable with p=0.5 as would have been done to represent the same probability by a Bernoulli sequence). We use nonBernoulli sequences of random permutations of fixed number of ones and zeros in order to encode the initial input probabilities.
A logic circuit operating on stochastic bit streams (see Fig. 4 for an example), accepts as input random sequences representing the probability of each input being one and produces ones and zeros like any digital circuit [33] i.e. a stochastic logic circuit uses Boolean gates to operate on sequences of random bits. Each bitstream represents a stochastic number interpreted as the probability of seeing a one in an arbitrary position. Thus, the computations performed by such a circuit are probabilistic in nature. The output bit stream produced can be decoded as the probability of the output being one by counting the number of ones in the stream and dividing by its length.
The vulnerability of a node is defined as the probability that the system produces incorrect output if that particular node is dysfunctional (faulty) i.e. it is the probability that the output of the network is different when that node is dysfunctional and is the complement of reliability. The procedure to determine the node vulnerabilities is illustrated in Fig. 5 is as follows. We generate nonBernoulli sequences of length L=1,000,000 in which exactly half of the bits are set to one at each of the seven initial inputs. The input stochastic sequences are propagated through both the original errorfree circuit and the circuit in which the node of interest is dysfunctional. As discussed in the previous section, the dysfunction of a node is represented by a corresponding stuckat fault of the requisite type at the particular location. This produces two sets of stochastic bit streams, one at each of the primary outputs of the faultfree circuit and the other at the primary outputs of the unreliable circuit. The proportion of ones in the output bit stream encodes the output signal probabilities i.e. the probability of the output being one. Since the reliability of the circuit under the fault is the probability that the circuit output is same as that of the faultfree circuit, the sequence encoding the output reliability can be obtained from the output sequence of the faulty circuit by comparing it to the output sequence of the faultfree circuit and setting each bit to one whenever the corresponding bits in the sequences are the same and zero if they are different. The proportion of ones in this resulting sequence will then correspond to the reliability of that output. Thus, we can obtain the stochastic sequence representing the reliability of each output by taking the XOR of each output bit stream of the faulty circuit with the complement of the corresponding output bitstreams of the faultfree circuit. For a circuit with multiple primary outputs as is the case here, the stochastic sequence encoding the joint output reliability can be obtained by taking the stochastic AND of the outputs of the XOR gates as the stochastic AND operation on the output of XOR gates produces a one only if all the corresponding bits at each XOR gate are one i.e. if all the corresponding bits in the respective outputs of the faultfree and faulty circuit are same. We then take the complement of the bit stream at the output of this AND gate to obtain the stream that encodes vulnerability. This bit stream can then be decoded to determine the node vulnerability with the proportion of ones in this stream equivalent to the vulnerability of the node.
The procedure for computing the vulnerability of a node described above and depicted in Fig. 5 is summarized as follows:

Generate nonBernoulli streams encoding input probabilities at each of the primary inputs.

Propagate the input binary streams through the faultfree circuit and obtain a random bit sequence for each output.

Propagate the same input binary streams through the circuit with a stuckat fault at the location of the node whose vulnerability we want to determine and again obtain a random bit sequence for each output.

XOR each primary output sequence from the faulty circuit obtained in step 3 with the complement of the corresponding primary output sequence from the faultfree circuit.

AND all the sequences obtained from each XOR gate. Take the complement of the stream so obtained. The vulnerability of the node is the fraction of ones in the resulting bit stream.
Thus, in a nutshell, the node vulnerabilities are obtained by propagating the initial input stochastic bit streams encoding the input probabilities through both the faulty and faultfree circuit, comparing the respective outputs obtained from each and decoding probabilities from the resulting streams.
Let x _{1},x _{2},…,x _{ N } represent input nonBernoulli sequences of length L with each sequence represented as a vector of length L whose i^{th} component is equal to the i^{th} bit in the sequence. Define the L×N matrix \(X=(\mathbf {x}_{1}^{\top }\;\mathbf {x}_{2}^{\top }\;\dotsc \;\mathbf {x}_{N}^{\top })\). Thus, each row of this matrix contains the corresponding bits of each of the primary input streams. The vulnerability v _{ j } of node j∈{1,2,…,P} is given by:
where ^{′} is the bitcomplement operator and ⊕ is the binary XOR operator.
Results and discussion
Simulation results for drug intervention
We use the Boolean network model to determine an apposite therapy for each fault. As described in the methods section, the best targeted therapy for each of the cancerinducing faults is the one that under the presence of the fault, produces the output with the smallest Hamming weight with the minimal number of drugs. The best therapy for each of the faults is shown in table 1 with the drug vector defined as before. Note that for certain faults, no drug vector can improve the output. Such faults are said to be untestable since no test (drug vector in this case) can rectify the fault. This is because there are no drugs on the fanout of these genes. However, all these faults with the exception of fault 18 are minimally proliferative as they produce a faulty output with the least possible Hamming weight of one.
Thus, there are many locations in the gene regulatory network of prostate cancer where malfunctions can occur resulting in a cancer that is different, requiring a specific targeted therapy. The table facilitates arriving at such a therapy as it maps each malfunction to an appropriate set of drugs. The lookup table can be used to devise therapies that have a higher likelihood of success since they are tailored specifically to the molecular abnormalities in critical pathways and thereby facilitates an individualized approach to therapy design.
In order to find the best combination therapy across all possible faults, as discussed in the methods section, for each drug combination we determine the sum of the Hamming weights of the output vector across all possible combinations of faults and choose the drug combination that yields the smallest total. This gives us the drug cocktail of AZD6244, AZD5363 and Enzalutamide as a combination therapy for advanced prostate cancer. In a recent study, the drug combination of AZD5363 and Enzalutamide has demonstrated an impressive response in prostate cancer models [34]. Moreover, AZD6244 in partnership with an AKT pathway inhibitor (analogous to AZD5363), has been proposed as a strategy for the treatment of CRPC [35]. Thus, we propose that the aforementioned drug triad which represents a horizontal blockade approach, wherein combination therapy is used for the concerted pharmacologic inhibition of multiple compensatory pathways, as a therapeutic modality that may attenuate prostate cancer survival and growth.
Node vulnerabilities
Using the framework delineated in the methods section, we quantify the vulnerability of different nodes. The vulnerability values obtained are given in Table 2. Vulnerability assessment can be used to identify candidates for targeted drug development. Nodes whose vulnerabilities are higher should be presumably better targets for drugs since potentially therapeutic benefit is more likely for nodes which are more vulnerable. We observe that the ARmediated signaling axis remains a valid target. Furthermore, we see that dysfunction in the AKT nexus and the loss of tumorsuppressors have higher vulnerability values so drugs that attempt to alleviate these aberrations should be efficacious in attenuating tumor growth. The design of anticancer therapeutics directed at the loss of tumor suppressors has been difficult [36]. Additionally, AKTselective drug development is challenging due to its homology with other kinases [37]. These complications notwithstanding, accelerated development of novel agents that target these aberrations is warranted. In contrast, the vulnerabilities for certain nodes such as those in the mTOR axis are low indicating that they might not be attractive targets for drug development. Indeed, marginal clinical activity has been observed for mTOR inhibition with agents such as everolimus and temsirolimus failing to impact tumor proliferation in men with prostate cancer [4, 38]. Finally, in terms of the key pathways implicated in the disease we see that castrationresistant prostate cancer shows most vulnerability on aggregate to dysfunction in the AKT pathway. In a study it was demonstrated that the AKT pathway dominates AR signaling in CRPC [39].
Conclusion
Castrationresistant prostate cancer is a hormone refractory phenotype of significant morbidity and mortality in the prostate cancer disease continuum where patients no longer respond to androgen ablation therapy. The biomolecular network representing the signaling pathways involved in the pathogenesis of this lethal malignancy is translated to a digital circuit. The locations of possible malfunctions in the digital circuit are identified and computer simulation of the equivalent model is used to predict effective therapies that mitigate the effect of different faults. A prospectively attractive combinatorial therapeutic strategy for the constellation of abnormalities is to leverage an AR axis targeted agent in conjunction with reciprocal inhibitors of other dysregulated pathways that are fundamental in coordinately driving oncogenesis. Proof of principle of clinical use for the proposed regimen remains to be demonstrated. A reliability (vulnerability) analysis methodology of digital circuits premised on stochastic logic modeling is utilized to quantify the vulnerability of the network to the dysfunction in discrete components in the signaling cascade thereby identifying key variables as targets for intervention that conceivably might be exploited by a new generation of novel therapeutics. These findings can contribute to the development of new rational approaches for the possible treatment of androgenrefractory prostate cancer. There is however a paucity of companion predictive biomarkers that can be used for the stratification of patients based on molecular aberrations in order to prescribe the apposite treatment. Furthermore, the histological and clinical heterogeneity of CRPC and the inherent redundancy along with the presence of feedback loops in pathways whose molecular underpinnings in the context of the disease induction and development are not yet fully understood, tender any potential translation into objective clinical efficacy of therapeutic implications derived from computations fraught with challenges.
Abbreviations
 ADT:

Androgen deprivation therapy
 AR:

Androgen receptor
 BCL2:

Bcell lymphoma 2
 BN:

Boolean network
 CDK2:

Cyclin dependent kinase 2
 CRPC:

Castration resistant prostate cancer
 EGF:

Epidermal growth factor
 GRN:

Gene regulatory network
 HBEGF:

Heparin binding EGFlike growth factor
 IGF:

Insulinlike growth factor
 MAPK:

Mitogen activated protein kinase
 MIMO:

Multiinput multioutput
 NRG1:

Neuregulin 1
 PSA:

Prostate specific antigen
 PTEN:

Phosphatase and tensin homolog
 SP1:

Specificity protein 1
 TMPRSS2:

Transmembrane protease serine 2
References
 1
Siegel RL, Miller KD, Jemal A. Cancer statistics 2015. CA: a cancer journal for clinicians. 2015; 65:5–29.
 2
Feng J, Zheng SL, Liu W, Isaacs WB, Xu J. Androgen receptor signaling in prostate cancer: new twists for an old pathway. J Steroids Hormon Sci. 2011.
 3
Boyd LK, Mao X, Lu YJ. The complexity of prostate cancer: genomic alterations and heterogeneity. Nat Rev Urol. 2012; 9(11):652–64.
 4
Derleth CL, Evan YY. Targeted therapy in the treatment of castrationresistant prostate cancer. Oncology. 2013; 27(7):620–30.
 5
Leo S, Accettura C, Lorusso V. Castrationresistant prostate cancer: targeted therapies. Chemotherapy. 2010; 57(2):115–27.
 6
Chen Y, Clegg NJ, Scher HI. Antiandrogens and androgendepleting therapies in prostate cancer: new agents for an established target. Lancet Oncol. 2009; 10:981–91.
 7
Patel JC, Maughan BL, Agarwal AM. Batten AM. Zhang TY: Agarwal N. Emerging molecularly targeted therapies in castration refractory prostate cancer. Prostate Cancer; 2013.
 8
Agarwal N, Sonpavde G, Sternberg CN. Novel molecular targets for the therapy of castrationresistant prostate cancer. Eur Urol. 2012; 61(5):950–60.
 9
Aggarwal R, Ryan CJ. Castrationresistant prostate cancer: targeted therapies and individualized treatment. Oncologist. 2011; 16(3):264–75.
 10
Nedumparambathmarath V, Chakrabarti SK, Sreekumar J. Modeling of gene regulatory networks: a review. J Biomed Sci Eng. 2013; 6:223–31.
 11
Kaufmann SA. The origins of order: selforganization and selection in evolution. New York:Oxford University Press; 1993.
 12
Helikar T, Konvalina J, Heidel J, Rogers JA. Emergent decisionmaking in biological signal transduction networks. Proc Natl Acad Sci. 2008; 105(6):1913–8.
 13
Albert R, Thakar J. Boolean modeling: a logicbased dynamic approach for understanding signaling and regulatory networks and for making useful predictions. Wiley Interdiscip Rev Syst Biol Med. 2014; 6:353–69.
 14
Wang R, Saadatpour A, Albert R. Boolean modeling in systems biology: an overview of methodology and applications. Phys Biol. 2012;9(5).
 15
Fauré A, Naldi A, Chaouiya C, Thieffry D. Dynamical analysis of a generic Boolean model for the control of the mammalian cell cycle. Bioinformatics. 2006; 22(14):e124–e31.
 16
Schlatter R, Schmich K, Vizcarra IA, Scheurich P, Sauter T, Borner C, Ederer M, Merfort I, Sawodny O. ON/OFF and beyonda Boolean model of apoptosis. PLoS Comput Biol. 2009; 5(12):e1000595.
 17
Zhang R, Shah MV, Yang J, Nyland SB, Liu X, Yun JK, Albert R, Loughran TP. Network model of survival signaling in large granular lymphocyte leukemia. Proc Natl Acad Sci. 2008; 105(42):16308–13.
 18
Sridharan S, Varghese R, Venkatraj V, Datta A. Hypoxia stress response pathways: modeling and targeted therapy. IEEE J Biomed Health Inform. in press. http://ieeexplore.ieee.org/document/7460895/.
 19
Giacomantonio CE, Goodhill GJ. A Boolean model of the gene regulatory network underlying mammalian cortical area development. PLoS Comput Biol. 2010; 6(9):e1000936.
 20
Carver BS, Chapinski C, Wongvipat J, Hieronymus H, Chen Y, Chandarlapaty S, Arora VK, Le C, Koutcher J, Scher H, Scardino PT, Rosen N. Reciprocal feedback regulation of PI3K and androgen receptor signaling in PTENdeficient prostate cancer. Cancer Cell. 2011; 19(5):575–86.
 21
Yuen HF, Abramcyzk O, Montgomery G, Chan KK, Huang YH, Sasazuki T, Shirasawa S, Gopesh S, Chan KW, Fennell D, Janne P, ElTanani M, Murray JT. Impact of oncogenic driver mutations on feedback between the PI3K and MEK pathways in cancer cells. Biosci Rep. 2012; 32(4):413–22.
 22
Kanehisa M, Goto S, Furumichi M, Tanabe M, Hirakawa M. KEGG for representation and analysis of molecular networks involving diseases and drugs. Nucleic Acids Res. 2010; 38:D355–D60.
 23
Kanehisa M, Goto S, Hattori M, AokiKinoshita KF, Itoh M, Kawashima S, Katayama T, Araki M, Hirakawa M. From genomics to chemical genomics: new developments in KEGG. Nucleic Acids Res. 2006; 34:D354–D7.
 24
Kanehisa M, Goto S. KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Res. 2000; 28:27–30.
 25
Watterson S, Marshall S, Ghazal P. Logic models of pathway biology. Drug Discov Today. 2008; 13(9):447–56.
 26
Abramovici M, Breuer MA, Friedman AD. Digital systems testing and testable design. New York: WileyIEEE Press; 1994.
 27
Layek R, Datta A, Bittner M, Dougherty E. Cancer therapy design based on pathway logic. Bioinformatics. 2011; 27(4):548–55.
 28
Choudhury MR, Mohanram K. Reliability analysis of logic circuits. Comput Aided Des Integr Circ Syst, IEEE Trans. 2009; 28(3):392–405.
 29
Han J, Chen H, Liang J, Zhu P, Yang Z, Lombardi F. A stochastic computational approach for accurate and efficient reliability evaluation. IEEE Trans Comput. 2014; 63(6):1336–50.
 30
Abdi A, Tahoori MB, Emamian ES. Fault Diagnosis Engineering of Digital Circuits Can Identify Vulnerable Molecules in Complex Cellular Pathways. Sci Signal. 2008; 1(42):ra10.
 31
Zhu P, Aliabadi HM, Uludağ H, Han J. Identification of Potential Drug Targets in Cancer Signaling Pathways using Stochastic Logical Models. Nat Sci Rep. 2016; 6.
 32
Liang J, Han J. Stochastic Boolean networks: an efficient approach to modeling gene regulatory networks. BMC Syst Biol. 2012; 6.
 33
Zhu P, Han J. Asynchronous stochastic Boolean networks as gene network models. J Comput Biol. 2014; 21(10):771–83.
 34
Toren P, Kim S, Cordonnier T, Crafter C, Davies BR, Fazli L, Gleave ME, Zoubeidi A. Combination AZD5363 with enzalutamide significantly delays enzalutamideresistant prostate cancer in preclinical models. Eur Urol. 2015; 67(6):986–90.
 35
Park H, Kim Y, Sul JW, Jeong IG, Yi HJ, Ahn JB, Kang JS, Yun J, Hwang JJ, Kim CS. Synergistic anticancer efficacy of MEK inhibition and dual PI3K/mTOR inhibition in castrationresistant prostate cancer. Prostate. 2015; 75(15):1747–59.
 36
Dillon LM, Miller TW. Therapeutic targeting of cancers with loss of PTEN function. Curr Drug Targets. 2014; 15:65–79.
 37
Bitting RL, Armstrong AJ. Targeting the PI3K/Akt/mTOR pathway in castrationresistant prostate cancer. Endocr relat Cancer. 2013; 20(3):R83–R99.
 38
Sarker D, Reid AHM, Yap TA. de Bono JS. Targeting the PI3K/AKT pathway for the treatment of prostate cancer. Clin Cancer Res. 2009; 15(15):4799–805.
 39
Kaarbø M, Mikkelsen ØL, Malerød L, Qu S, Lobert VH, Akgul G, Halvorsen T, Mælandsmo GM, Saatcioglu F. PI3KAKTmTOR pathway is dominant over androgen receptor signaling in prostate cancer cells. Anal Cell Pathol. 2010; 32:11–27.
Acknowledgements
Not applicable.
Funding
Publication costs for this article have been funded by the National Science Foundation (NSF) under grant ECCS1404314. This work was supported in part by NSF under grant ECCS1404314 and the Texas Engineering Experiment Station (TEES)  Agrilife Center for Bioinformatics and Genomics Systems Engineering.
Availability of data and materials
Data sharing is not applicable to this article as no datasets were generated or analyzed during the current study.
Authors’ contributions
OAA developed the computational modeling and simulations and wrote the manuscript. AD conceived the main idea. Both authors read and approved the final manuscript.
Competing interests
The authors declare that they have no competing interests.
Consent for publication
Not applicable.
Ethics approval and consent to participate
Not applicable.
About this supplement
This article has been published as part of BMC Bioinformatics Volume 18 Supplement 4, 2017: Selected original research articles from the Third International Workshop on Computational Network Biology: Modeling, Analysis, and Control (CNBMAC 2016): bioinformatics. The full contents of the supplement are available online at https://bmcbioinformatics.biomedcentral.com/articles/supplements/volume18supplement4.
Author information
Affiliations
Corresponding author
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Arshad, O.A., Datta, A. Towards targeted combinatorial therapy design for the treatment of castrationresistant prostate cancer. BMC Bioinformatics 18, 134 (2017). https://doi.org/10.1186/s1285901715222
Published:
Keywords
 Prostate cancer
 Gene regulatory networks
 Boolean modeling
 Combination therapy
 Stochastic logic
 Vulnerability assessment