Computational Insights on the Competing Effects of Nitric Oxide in Regulating Apoptosis

Despite the establishment of the important role of nitric oxide (NO) on apoptosis, a molecular- level understanding of the origin of its dichotomous pro- and anti-apoptotic effects has been elusive. We propose a new mathematical model for simulating the effects of nitric oxide (NO) on apoptosis. The new model integrates mitochondria-dependent apoptotic pathways with NO-related reactions, to gain insights into the regulatory effect of the reactive NO species N2O3, non-heme iron nitrosyl species (FeLnNO), and peroxynitrite (ONOO−). The biochemical pathways of apoptosis coupled with NO-related reactions are described by ordinary differential equations using mass-action kinetics. In the absence of NO, the model predicts either cell survival or apoptosis (a bistable behavior) with shifts in the onset time of apoptotic response depending on the strength of extracellular stimuli. Computations demonstrate that the relative concentrations of anti- and pro-apoptotic reactive NO species, and their interplay with glutathione, determine the net anti- or pro-apoptotic effects at long time points. Interestingly, transient effects on apoptosis are also observed in these simulations, the duration of which may reach up to hours, despite the eventual convergence to an anti-apoptotic state. Our computations point to the importance of precise timing of NO production and external stimulation in determining the eventual pro- or anti-apoptotic role of NO.


Introduction
The survival of an organism depends on homeostatic mechanisms that establish a balance between cell proliferation and cell death. Apoptosis, a form of programmed cell death, assists in regulating cell proliferation. This process stands in contrast to necrosis, which is thought to be uncontrolled. Dysregulation of apoptosis has been implicated in various disease processes in which the cells apoptose to a higher or lower extent compared to those in healthy tissues [1]. When cells undergo apoptosis, a series of morphological and biochemical changes occur, the mechanisms of which are current topics of broad interest [2].
Apoptosis may be induced by various events, such as binding of extracellular (EC) death signaling ligands to host cell receptors, the lack of pro-survival signals, and genetic damage. These events are usually followed by the activation of caspases, cysteine-dependent aspartate-specific proteases, which initiate and execute apoptosis. Caspases are activated through two major pathways: (a) liganddependent or receptor-induced activation (extrinsic pathway), involving death receptors such as Fas or the members of tumor necrosis factor (TNF) receptor superfamily, and (b) mitochondriadependent activation (intrinsic pathway) via cytochrome c (cyt c) release from mitochondria, triggered by stress, irradiation or inflammation [3,4].
Binding of death ligands such as Fas ligand (FasL), TNF, or tumor necrosis-related apoptosis-inducing ligand (TRAIL) usually induces the oligomerization of associated receptors, followed by binding of adaptor proteins, e.g., Fas-Associated Death Domain proteins (FADD), to the cytoplasmic domains of the receptors [5]. The resulting Death Inducing Signaling Complex (DISC) recruits multiple procaspase-8 molecules that mutually cleave and activate one another into caspases-8 (casp8). In Type I cells, large quantities of casp8 activate other caspases including the executioner caspase-3 (casp3) molecules that ultimately lead to apoptosis. In Type II cells, on the other hand, the amount of casp8 activated at the DISC is small, such that the activation of casp8 does not propagate directly to casp3, but instead is amplified via the mitochondria.
Nitric oxide has opposite, competing effects in regulating apoptosis: it exerts an anti-apoptotic effect on hepatocytes [6][7][8], endothelial cells [9][10][11][12][13] and keratinocytes [14], whereas it is proapoptotic in the case of macrophages [15][16][17][18]. The variability and complexity of the effects of NO on ultimate cellular fate may arise from this molecule's ability to react with oxygen, reactive oxygen species, metal ions, small thiol-containing molecules, and proteins. The resulting reactive NO species can either trigger or suppress apoptosis through various mechanisms. Chief among them is the S-nitrosative suppression of caspase activation, subsequent to the generation of FeL n NO or other species capable of carrying out Snitrosation reactions (see below) [7,19]. Differences in the levels of NO and its reaction products may also arise from diverse inflammatory settings in which the expression of nitric oxide synthases (NOS) is affected. For example, quiescent endothelial cells express constitutive NOS (eNOS) that directly produce NO molecules and mediate so-called ''direct'' effects [20]. Some inflammatory stimuli, on the other hand, lead to inducible NOS (iNOS) expression that subsequently generates reactive NO species, which in turn mediate ''indirect'' effects of NO. The simultaneous presence of oxygen radicals can generate other reactive NO species that mediate further indirect effects of NO [20]. As another example, hepatocytes and macrophages have different amounts of non-heme iron complexes, which affect the levels of iron-nitrosyl species when NO is produced [21]. Finally, different intracellular levels of glutathione (GSH) can also modulate the time evolution of NO-related compounds [22].
Computational approaches have been used previously to help unravel the complex biology of NO. Biotransport of NO was first modeled by Lancaster [23,24] followed by other groups, among them Zhang and Edwards [25] (reviewed by Buerk [26]). Recently, Hu and coworkers focused on detailed reaction mechanism of NO [22]. These models have shed light into the biotransport of NO and the types of chemical reactions that involve NO and related reactive species. Additionally, a number of mathematical models have been proposed for understanding the mechanisms of apoptosis [27][28][29][30][31][32][33][34][35], including in particular the work of Eissing et al., which demonstrated the importance of IAP inhibition for imparting bistability in type I cells [30], and that of Rehm et al. [33] and Legewie et al. [32] that showed the same effect in type II cells. These studies have improved our understanding of the robustness of switch mechanisms for regulating apoptosis, but none of them has addressed the dichotomous effects of NO [27][28][29][30][31][32][33][34][35].
Herein, we propose a mathematical model that may shed light on the pro-and anti-apoptotic effects of NO in specific contexts. The model we propose couples the apoptotic cascade [28] to an extended model of NO reaction pathways initially proposed by Hu et al. [22]. First, we illustrate how identical cells can undergo apoptosis at different time points after being exposed to apoptotic stimuli, in accord with experimental data collected on single cells [36,37]. Then, we examine the apoptotic behavior in response to changes in N 2 O 3 , FeL n NO, ONOO 2 , and GSH levels in the presence of NO production by iNOS. Our simulations provide insights into the origin of the dichotomous effects of NO on apoptosis observed in experiments.

Results
First, we illustrate how different strengths of EC pro-apoptotic signals may result in opposite qualitative responses or different quantitative (time-dependent) responses in the same type of cells [37], using our recently introduced bistable model [28] (illustrated in Figure 1A). Then, we examine the differences in the bistable response of diverse NO producing cells, e.g. cells with different concentrations of GSH and FeL n -and in different settings, i.e., with or without production of superoxide.

Delay in apoptosis induction (Model I)
Tyas et al. [37] showed that cells of the same type simultaneously subjected to EC stimuli initiate their apoptotic The components of the model are procaspase-8 (pro8), procaspase-3 (pro3), procaspase-9 (pro9), caspase-8 (casp8), caspase-9 (casp9), caspase-3 (casp3), IAP (inhibitor of apoptosis), cytochrome c (cyt c), Apaf-1, the heptameric apoptosome complex (apop), the mitochondrial permeability transition pore complex (PTPC), p53, Bcl-2, Bax, Bid, truncated Bid (tBid). The reader is referred to our previous work [28] for more details. Three compounds (N 2 O 3 , FeL n NO and ONOO 2 ) not included in the original Model I [28] are highlighted. These compounds establish the connection with the nitric oxide pathways delineated in panel B. (B) Nitric oxide (NO)-related reactions in Model II. The following compounds are included: ONOO 2 (peroxynitrite), GPX (glutathione peroxidase), O 2 2 (superoxide), GSH (glutathione), GSNO (nitrosoglutathione), GSSG (glutathione disulfide), CcOX (cytochrome c oxidase), SOD (superoxide dismutase), FeL n (nonheme iron compounds), FeL n NO (non-heme iron nitrosyl compounds), NADPH (reduced form of nicotinamide adenine dinucleotide phosphate), NADP+ (oxidized form of nicotinamide adenine dinucleotide phosphate). FeL n NO, ONOO 2 and N 2 O 3 , highlighted in both panels A and B, bridge between Models I to II. Model III integrates both sets of reactions/pathways through these compounds. GSH modulates their concentrations by reacting with them. GSH is converted by these reactions to GSNO, which is then converted to GSSG and finally back to GSH. Those compounds and interactions are shown in blue. See Table 1 for the complete list of reactions and rate constants. doi:10.1371/journal.pone.0002249.g001 response at different times. Figure 2  This analysis shows that cells of the same type may undergo apoptosis at different times due to their different EC microenvironments. Hence, the difference in the onset times among cells of the same type in a given cell culture may be explained without recourse to alterations in the underlying network of biochemical reactions [28].
Nitric oxide-associated network (Model II) ( Figure 1B) The results from our calculations using Model II are shown in Figure 3. Here, we focused on the time evolution of four compounds, GSH, N 2 O 3 , FeL n NO and ONOO 2 , displayed in respective panels A-D. The NO species N 2 O 3 , FeL n NO and ONOO 2 have been proposed to carry out various indirect effects of NO on cellular pathways, including apoptosis, during inflammation [20].
GSH is an anti-oxidant reduced to GSSG by reacting with nitrosative N 2 O 3 and FeL n NO, and with oxidative ONOO 2 ( Table 1). GSH is depleted to low levels in a switch-like manner due to those reactions (panel A). The depletion of GSH is accompanied by increases in N 2 O 3 and FeL n NO concentrations (panels B-C). On the other hand, this switch-like behavior is not that pronounced in [ONOO 2 ] time dependence (panel D). Simulations performed with different initial GSH concentrations (three different curves in each panel) change the steady-state concentrations of all three NO-related compounds that interfere with apoptotic pathways (panels B-D). The switch-like increase in [N 2 O 3 ] and non-switch-like increase in [ONOO 2 ] is in agreement with the results of Hu et al. [22].

Anti-apoptotic and pro-apoptotic effects of NO (Model III)
We analyze here the dynamics of the reduced mitochondriadependent apoptosis model coupled to anti-and pro-apoptotic pathways associated with NO; see Materials and Methods for the list of reactions/interactions/steps that come into play in this model (III). As mentioned above, NO-related pathways are coupled to apoptotic pathways through N 2 O 3 , FeL n NO, and ONOO 2 that are produced by the reaction of NO with O 2 , FeL n and O 2 2 , respectively. For simplicity, those effects of NO mediated by cGMP [38,39] are not included in this initial mathematical model.
Modulating roles of N 2 O 3 and GSH in apoptosis. We initially excluded non-heme iron compounds in order to assess the effect of N 2 O 3 exclusively. The production rate of superoxide was likewise assumed to be zero. N 2 O 3 is produced by reactions (xii) and (xiii) in Table 1. NO production and EC stimulation were initiated simultaneously. These results suggest that N 2 O 3 does not affect the bistable character of the response to EC stimuli, except for modifying the threshold for onset of apoptosis, which is shifted to higher [casp8] 0 (i.e. rendered more difficult) with decreasing [GSH] 0 . However, high initial concentrations of GSH restore the threshold back to 8.35610 25 mM. Therefore, N 2 O 3 can serve as an effective modulator of apoptosis provided that the level of GSH in the system is sufficiently low.
Effect of N 2 O 3 on the threshold degradation rates of Bax for transition from bistable to monostable behavior. In our previous computational study of apoptotic pathways, we observed a bistable behavior (selecting between cell death and survival) for degradation rates of Bax (m Bax ) lower than a threshold value (0.11 s 21 ), while monostable cell survival was predicted when m Bax .0.11 s 21 ( Figure 4A in Ref. [28]). This critical value of m Bax for the transition from bistability to monostability is called a limit point. We explored how the inclusion of NO reactions affects these findings. The limit point value of the Bax degradation rate for monostable cell survival is found to remain unchanged (at 0.11 s 21 ) for the range 10 3 #[GSH] 0 #10 4 mM. However, it decreases to 0.098 s 21 for [GSH] 0 = 10 2 mM and 0.096 s 21 for [GSH] 0 = 0 mM in the present model. The model again predicts that N 2 O 3 is not influential when the GSH level is sufficiently high in the cell.
Roles of non-heme iron complexes and GSH in apoptotic response. One of the important anti-apoptotic effects of NO is presumed to occur via its ability to react with non-heme iron complexes (FeL n ) to form FeL n NO. These species inhibit caspases by S-nitrosating the catalytic cysteine in the active site of these enzymes [19,40,41].
The results are presented in Figure 5, panels A-F, organized similarly to Figure 4 (i.e. using different [casp8] 0 in each row, and different [GSH] 0 in the two columns). Our calculations suggest that when the FeL n concentration is higher than 0.03 mM, there are no longer two stable steady-states at long times: caspase-3 levels always decrease to zero, even though their time evolutions depend on [casp8] 0 and [GSH] 0 . Yet, depending on the level of GSH, both apoptosis and cell survival may be possible. Panels A-C correspond to relatively high [GSH] 0. In panel A, [casp3] decreases to 10 28 mM that is less than 1 molecule per cell, hence zero, from 10 25 mM within the first two hours. However, in panels B and C, [casp3] increases to nanomolar values and remains at those levels for more than three hours. Caspase-3 may cause enough damage to kill the cell before it is depleted at longer times. We note that lower [GSH] 0 (e.g. [GSH] 0 = 10 3 mM, panels D-F and [GSH] 0 ,10 3 mM, data not shown) do not permit the casp3 concentration to reach such pro-apoptotic levels and monostable cell survival is observed irrespective of [casp8] 0 .
Various cell types subject to different intracellular microenvironments, or even the same cells under different settings (e.g. healthy state vs. inflammation or oxidative stress), may produce or experience different reactive NO intermediates [7,20,42]. For example, more FeL n NO may be produced in hepatocytes than in RAW264.7 macrophage-like cells due to the high level of nonheme iron complexes in hepatocytes [21]. In our previous study, RAW264.7 cells underwent apoptosis in the presence of NO; conversely, no casp3 activation was observed in either hepatocytes or iron loaded RAW264.7 cells [21]. The results (Figure 4 and data not shown) suggest that in cells with iron concentrations lower than 0.03 mM (e.g. RAW264.7 cells), both cell survival and apoptosis are possible depending on the strength of apoptotic stimuli (in agreement with our experimental results) [21]. However, a change in the intracellular environment of the same cell can change the response. Figure 5D-F shows that casp3 is not activated in the presence of non-heme iron ([FeL n ] 0 = 0.05 mM) when [GSH] 0 = 10 3 mM and [GSH] 0 ,10 3 mM (data not shown). We also checked if casp3 is activated when [casp8] 0 is as high as 0.1 mM when [GSH] = 10 3 mM. In this case, caspase-3 concentration increased to 0.0007 mM for approximately 5 minutes, an apoptotic stimulus that is likely insufficient for apoptosis. This prediction is in good agreement with our observation that caspase-3 is not activated in non-heme iron-loaded RAW264.7 cells whose [GSH] 0 does not reach 10 4 mM [21].
Roles of ONOO 2 and GSH in apoptotic response. The mechanism by which NO or its reactive species exert pro- apoptotic effects is not well established [43]. In the present study, we assume that the pro-apoptotic effect of NO occurs via formation of ONOO 2 , as has been suggested from a large number of experimental studies both in vitro and in vivo [44,45]. Experimental studies suggest that ONOO 2 may induce the opening of mitochondrial permeability transition pores (MPTPs) and subsequent cyt c release from mitochondria [38]. The possible mechanisms of cyt c release from mitochondria are diverse and controversial [46,47]. In our model, we assume that cyt c release is mediated by activation of MPTPs, independent of Bax channel formation on mitochondria. The complex that forms the MPTPs is called mitochondrial permeability transition pore complex (PTPC). The complex consists of peripheral benzodiazepine receptor, cyclophilin D, adenine nucleotide translocator (ANT), voltage-dependent anion channel (VDAC), and other proteins [48]. ANT is proposed to be converted from a specific transporter to a non-specific pore which then releases cyt c into the cytoplasm and subsequently induces apoptosis. It has been suggested that ONOO 2 acts on PTPC, specifically on ANT, to convert it to a non-specific pore (PTPC act ) [49]. We represent this process as: Cytochrome c is then released from the pore formed by PTPC act cyt c mito zPTPC act ?cyt czPTPC act The results are shown in Figure 6. The initial concentration of PTPC is assumed to be high (0.01 mM). At that value, Model I predicts the response to apoptotic stimuli to be monostable apoptosis (Figure 6 in Ref. [28]). We see a similar response in Figure 6A; a low initial value of casp8 (10 25 mM) results in an increase of [casp3] to nanomolar levels. Casp3 activation was observed with even lower values of [casp8] 0. However, casp3 does not reach nanomolar concentrations when [GSH] 0 = 10 3 mM (Figures 6D-F) and [GSH] 0 ,10 3 mM (data not shown). Initial concentrations [casp8] 0 higher than 1.5610 24 mM did not change this prediction.
These results suggest that in cells with large numbers of MPTPs (probably with high numbers of mitochondria), there are two possible outcomes in the presence of NO and O 2 2 production:  Table 4 is included in the model while those involving FeL n NO and ONOO 2 (reactions (xx, xxiii-xxv) are not, assuming FeL n concentration and rate of formation of superoxide to be zero. The solid curves depict the time evolution of [casp3], and dotted curves refer to [GSH]. The three rows of panels are the counterparts of those in Figure  and ONOO 2 ) that can cause apoptosis. The reason for this paradoxical prediction is that GSH has both protective and pro-apoptotic effects in our simulations: it exerts apoptotic effects via its reaction with anti-apoptotic N 2 O 3 and FeL n , and protective effects due to its reaction with pro-apoptotic O 2 2 and ONOO 2 . Simulations ( Figure 6) suggest that the proapoptotic effect of GSH is stronger than its protective effect using the interactions and parameters adopted in current simulations.
To examine the possibility of an alternative response, we repeated the computations depicted in Figure 6 in the absence of O 2 (so that N 2 O 3 is not produced) and FeL n . We also used initial PTPC concentration of 0.0001 mM, at which Model I predicts bistability ( Figure 6 in ref [28]). As seen in Figure 7 . The present analysis thus shows that the protection by GSH against oxidative stress is possible provided that O 2 and FeL n levels are sufficiently low.

Discussion
We present here the results from simulations that incorporate the main chemical interactions of NO with components of the apoptotic interactions network, with the goal of shedding light on the dichotomous effects of NO on apoptosis. Based on previously published studies, we considered N 2 O 3 and FeL n NO to be antiapoptotic and ONOO 2 pro-apoptotic. The results predict that cell survival or apoptosis is determined by a complex interplay among these reactive NO species and GSH. We observed that relative concentrations of anti-apoptotic and pro-apoptotic species determine the ultimate cell fate at late time points. Interestingly, transient apoptotic effects were observed under specific conditions (e.g. Figure 5 panels B-C). These intriguing findings point to the importance of the timing of NO production and apoptotic stimuli in determining the actual anti-or pro-apoptotic effect, even if steady state conditions favor cell survival, in agreement with our previous observations [50][51][52]. Another interesting effect we observed in our simulations was the time shift/delay in the onset of apoptosis in the presence of weak EC stimulus (panel B-D in Figure 2), consistent with the experiments of Tyas et al. [37].
Our simulations suggest that N 2 O 3 and non-heme iron nitrosyl form in a switch-like manner after depletion of GSH. ONOO 2 formation, on the other hand, hardly shows any switch-like behavior. We further found that N 2 O 3 does not eliminate the bistability between cell survival and apoptosis, but rather increases the threshold [casp8] 0 for onset of apoptosis. However, high initial concentrations of GSH restore the threshold back to its original value. Therefore, we would predict, non-intuitively, that N 2 O 3 does not influence cell survival when [GSH] 0 level is high.
On the other hand, our simulations suggest that there are no longer two stable steady states (cell survival and apoptosis) in the presence of non-heme iron at a level higher than a threshold value. Caspase-3 levels always decrease to zero even though its time evolution may depend on [casp8] 0 and [GSH] 0 . Yet, despite the steady state conditions that favor cell survival, executioner caspase concentrations can reach and retain apoptotic levels for several hours before they level off, when [GSH] 0 is high. When [GSH] 0 is low, on the other hand, our simulations predict resistance to apoptosis, in agreement with experimental observation [21].
In Tiedge et al. [53] have shown that pancreatic beta cells have low anti-oxidant levels (notably, GSH) and that the number of mitochondria is a determining factor in survival. They have also shown that transfection of the cells with a peroxide-inactivating enzyme, catalase, can protect against high-glucose induced apoptosis. An interesting experiment would be to correlate the number of mitochondria in the transfected cells with their survival. Oyadomari et al. [54] have shown that the endoplasmic reticulum (ER) plays a crucial role in the fate of NO-sensitive beta cells via calcium signaling. A natural next step in the present model would be to include these effects via a model which incorporates the effects of NO on the ER.
Our results are subject to several limitations. While we have adopted values for kinetic parameters and concentrations in accord with experimental data whenever available (Tables 1 and  2), many of the true intracellular rate constants for the reactions in our simulations are unknown. Given that the observed apoptotic responses are so sensitive to model parameters, detailed knowledge of reaction mechanisms and accurate values of rate constants are needed in modeling reaction networks as complicated as the ones presented here. Due to an extensive literature basis, we have posited that the pro-apoptotic NO species is ONOO 2 ; however, other species may in fact exert this effect. Additionally, the hypotheses raised by our simulations remain to be tested by further experiments. Some of the predictions could be tested by iron chelation and/or treatment with superoxide donors in a cell-free system or in single-cell studies, though each of these manipulations may have additional, artifactual effects. The hypothesis of bistability with regards to the apoptotic response can be tested as suggested by Legewie et al. [32], either in cell free-systems by adding caspase-3 or in single living cells by microinjecting caspase-3. The time evolution of caspase-3 can be monitored by fluorescent caspase-3 substrates. The time needed for caspase-3 activation will increase abruptly as caspase-3 concentration added will approach threshold value in a bistable system ( Figure 2D). Such combined experimental and computational studies may potentially help us understand and design therapeutics for diseases associated with apoptosis dysregulation.

Models
Three models are considered in this study. Model I, proposed in our earlier work [28], focuses on the pathways involved in mitochondria-dependent apoptosis ( Figure 1A). Model II is an extension of the kinetic model of NO-associated reactions recently proposed by Hu et al. [22] ( Figure 1B). Finally, Model III is the integration of Models I and II, proposed in the present study, to examine the pro-apoptotic and anti-apoptotic effects of NO.  [GSH] 0 = 10 4 mM (or otherwise specified) [22] doi:10.1371/journal.pone.0002249.t002 within a short time interval (,20 minutes) after initiation of the simulations for [GSH] 0 #10 3 mM and within four and half hours for [GSH] 0 = 10 4 mM, whereas five compounds (superoxide dismutase (SOD), glutathione peroxidase (GPX), CO 2 , O 2 , and cyt c) retain their equilibrium concentrations. Table 2 lists the initial and equilibrium concentrations different from zero, adopted in Model II, and the corresponding references.

Model III-Effects of NO-related reactions on apoptotic pathways
Model III combines Models I and II upon inclusion of the additional reactions presented in Table 4. See the highlighted compounds in Figure 1, which point to the species that couple the apoptotic and NO pathways. We note that ONOO 2 has a proapoptotic effect, while N 2 O 3 and FeL n NO (reactions labeled (xxii)-(xxv)) deactivate the caspases, thus inducing anti-apoptotic effects. The associated rate constants and references are given in Table 4. Table 5 provides the rate expressions (rows 30-34) and differential rate equations (rows [35][36][37][38][39][40][41][42][43]