Apoptosis is a highly regulated cell death mechanism involved in many physiological processes. A key component of extrinsically activated apoptosis is the death receptor Fas which, on binding to its cognate ligand FasL, oligomerize to form the death-inducing signaling complex. Motivated by recent experimental data, we propose a mathematical model of death ligand-receptor dynamics where FasL acts as a clustering agent for Fas, which form locally stable signaling platforms through proximity-induced receptor interactions. Significantly, the model exhibits hysteresis, providing an upstream mechanism for bistability and robustness. At low receptor concentrations, the bistability is contingent on the trimerism of FasL. Moreover, irreversible bistability, representing a committed cell death decision, emerges at high concentrations which may be achieved through receptor pre-association or localization onto membrane lipid rafts. Thus, our model provides a novel theory for these observed biological phenomena within the unified context of bistability. Importantly, as Fas interactions initiate the extrinsic apoptotic pathway, our model also suggests a mechanism by which cells may function as bistable life/death switches independently of any such dynamics in their downstream components. Our results highlight the role of death receptors in deciding cell fate and add to the signal processing capabilities attributed to receptor clustering.
Apoptosis is a coordinated cell death program employed by multicellular organisms that plays a central role in many physiological processes. Normal function of apoptosis is critical for development, tissue homeostasis, cell termination, and immune response, and its disruption is associated with pathological conditions such as developmental defects, neurodegenerative disorders, autoimmune disorders, and tumorigenesis [1]–[5]. Due to its biological significance, much effort has been devoted to uncovering the pathways governing apoptosis. Indeed, recent progress has enabled the proliferation of mathematical models, both mechanistic and integrative [e. g. , 6]–[14], which together have offered profound insights into the underlying molecular interactions. The current work takes a similarly mathematical approach and hence inherits from this legacy. There are two main pathways of apoptotic activation: the extrinsic (receptor-mediated) pathway and the intrinsic (mitochondrial) pathway, both of which are highly regulated [15], [16]. In this study, we focus on the core machinery of the extrinsic pathway, which is initiated upon detection of an extracellular death signal, e. g. , FasL, a homotrimeric ligand that binds to its cognate transmembrane death receptor, Fas (CD95/Apo-1), in a 1∶3 ratio. This clusters the intracellular receptor death domains and promotes the ligation of FADD, forming the death-inducing signaling complex (DISC) [17]–[19]. The DISC catalyzes the activation of initiator caspases, e. g. , caspase-8, through death effector domain interactions. Initiator caspases then activate effector caspases, e. g. , caspase-3, which ultimately execute cell death by direct cleavage of cellular targets [20]–[23]. Apoptosis is typically viewed as a bistable system, with a sharp all-or-none switch between attracting life and death states. This bistability is important for conferring robustness [24]. Consequently, researchers have used computational models to identify and study potential sources of bistability in apoptosis, including positive caspase feedback [8], inhibition of DISC by cFLIP [7], cooperativity in apoptosome formation [10], double-negative caspase feedback through XIAP [11], and double-negative feedback in Bcl-2 protein interactions [25]. In this work, we propose that bistability may be induced upstream by the death receptors themselves. The current model of death ligand-receptor dynamics assumes that FasL activates Fas by direct crosslinking, producing a DISC concentration that varies smoothly with the ligand input [26]. However, recent structural data [27] suggests a different view. In particular, Fas was found in both closed and open forms, only the latter of which allowed FADD binding and hence transduction of the apoptotic signal. Moreover, open Fas were observed to pair-stabilize through stem helix interactions. This affords a mechanism for bistability, similar to the Ising model in ferromagnetism [28], where open Fas, presumably disfavored relative to their native closed forms [29], are able to sustain their conformations even after removal of the initial stimulus promoting receptor opening, past a certain critical density of open Fas. This induces hysteresis in the concentration of active, signaling receptors and therefore in apoptosis. We studied this proposed mechanism by formulating and analyzing a mathematical model. The essential interpretation is that FasL acts as a clustering platform for Fas, which establish contacts with other Fas through pairwise and higher-order interactions to form units capable of hysteresis (Figure 1). At low receptor concentrations, the model exhibits bistability provided that the number of receptors that each ligand can coordinate is at least three. This hence gives a theory for the trimeric character of FasL. Furthermore, at high concentrations, for example, through receptor pre-association [30]–[32] or localization onto lipid rafts [33], irreversible bistability is achieved, implementing a permanent cell death decision. Thus, our model suggests a primary role for death receptors in deciding cell fate. Moreover, our results offer novel functional interpretations of ligand trimerism and receptor pre-association and localization within the unified context of bistability.
Constructing a mathematical model of Fas dynamics is not entirely straightforward as receptors can form highly oligomeric clusters [27], [33]. A standard dynamical systems description would therefore require an exponentially large number of state variables to account for all combinatorial configurations. To circumvent this, we considered the problem at the level of individual clusters. Each cluster can be represented by a tuple denoting the numbers of its molecular constituents, the cluster association being implicit, so only these molecule numbers need be tracked. In our model, a cluster is indexed by a tuple, where represents FasL and, , and are three posited forms of Fas, denoting closed, open and unstable, and open and stable, i. e. , active and signaling, receptors, respectively. Within a cluster, we assumed a complete interaction graph and defined the reactions (1a) (1b) (1c) (1d) The first reaction describes spontaneous receptor opening and closing; the second, constitutive destabilization of open Fas; the third, ligand-independent receptor cluster-stabilization; and the fourth, ligand-dependent receptor cluster-stabilization (Figure 2). The orders of the cluster-stabilization events are limited by the parameters and, which capture the effects of receptor density and Fas coordination by FasL, respectively. Although only pair-stabilization () has been observed experimentally [27], higher-order analogues, for example, as facilitated by globular interactions, are not unreasonable. Formally, these reactions are to be interpreted as state transitions on the space of cluster tuples. However, the reaction notation is suggestive, highlighting the contribution of each elementary event, which we modeled using constant reaction rates (for simplicity, we set uniform rate constants and for all ligand-independent and -dependent cluster-stabilization reactions of molecularity, respectively). Then on making a continuum approximation, we reinterpreted the molecule numbers as local concentrations and applied the law of mass action to produce a dynamical system for each cluster in the concentrations of. Validity of the model requires that the molecular concentrations are not too low and that the timescale of receptor conformational change is short compared to that of cluster dissociation. To study the long-term behavior of the model, we solved the system at steady state (denoted by the subscript). Introducing the nondimensionalizations (2a) (2b) (2c) (2d) (2e) where is a characteristic concentration and is time, and (3a) (3b) (3c) (3d) this is (4a) (4b) where is the nondimensional total receptor density, and is given by considering (5) and solving with, a polynomial in of degree. Clearly, the model is bistable only if (two stable nodes must be separated by an unstable node as the model is effectively one-dimensional in). We used as a measure of the apoptotic activation of a cluster. In principle, all open receptors contribute to apoptotic signaling, but is small, at least at steady state (since due to the assumed prevalence of the closed form [29]), and so can be neglected. While measures the coordination capacity of FasL and hence may be equated with its oligomeric order (e. g. , in the biological context), an appropriate value for, relating to the total receptor concentration, is somewhat more elusive. Therefore, we began our analysis by performing a simple receptor density estimate. Approximating the cell as a cube of linear dimension m, the associated volume of pL implies the correspondence nM molecules molecules/nm on restricting to the membrane, i. e. , by averaging over the surface area of m. Thus, for a conservative receptor concentration estimate of nM [7], [9], [12], [13], the number of Fas molecules in the neighborhood of each receptor is only, assuming a charateristic size of nm. We hence found that receptors may be very sparsely distributed. In this low density mode, high-order Fas interactions in the absence of ligand can be neglected (). Therefore, in this context, bistability is possible only if, and the trimerism of FasL thus demonstrates the lowest-order complexity required for bistability. From the form of, this bistability is reversible as a function of the FasL concentration since the governing polynomial for is of degree only at. This suggests that at the cluster level, the cell death decision can be reversed, which may have adverse effects on cellular and genomic integrity. However, irreversible bistability at higher receptor densities may also be achieved. Researchers have observed tendencies for death receptors both to pre-associate as dimers or trimers [30]–[32] and to selectively localize onto membrane lipid rafts [33]. The result of either of these processes may be to increase the local receptor concentration. In this high density mode, we set, as the preceeding approximation is no longer valid. Irreversible bistability then becomes attainable, representing a committed cell death decision. For the remainder of the study, we incorporated both the low and high receptor density regimes into a single model by setting, using as a continuous transition parameter. Furthermore, we set to correspond to observed biology. Calculation of the steady-state activation curves showed that the model indeed exhibits bistability (Figure 3) for reasonable parameter choices (Methods). Thus, we established the possiblity of a novel bistability mechanism in extrinsic apoptosis. The associated hysteresis enables threshold switching between well-separated low and high activation states. Biologically, these define local signals of life and death, which are integrated at the cell level to compute the overall apoptotic response. As per the previous analysis, reversibility of the bistability is dependent on, with irreversibility emerging for sufficiently high. This suggests a bivariate parameterization of, producing a multivalued steady-state surface over -space (Figure 4). The result is a cusp catastrophe, an elementary object of catastrophe theory, which studies how small perturbations in certain parameters can lead to large and sudden changes in the behavior of a nonlinear system [34]. A more instructive view of the dependence of the model' s qualitative structure on and is shown in Figure 5. We then focused on the activation and deactivation thresholds, respectively, defining the bistable regime. These are the points at which the steady state switches discontinuously from one branch to the other, and are given by the values of at which the hysteresis curve turns, i. e. , at (Figure 6). We performed a sensitivity analysis of by measuring the effects of perturbing the model parameters about baseline values (Methods). For each threshold-parameter pair, we computed a normalized sensitivity by linear regression. Strong effects of, , and were observed (Figure 7); for the corresponding Fas thresholds at, respectively, the parameters, , , and were emphasized. Thus, the bistability thresholds do not appear particularly robust. However, the data reveal that essentially all parameter sets sampled were bistable. This suggests a weaker form of robustness, namely, robustness of bistability, which nevertheless supports life and death decisions over a wide operating range. To probe this further, we sampled parameters with increasing spread about baseline values and computed the fraction of parameter sets that remained bistable (Methods). The results show that has an exponential form (Figure 8). Extrapolating to, the data suggest an asymptotic bistable fraction of. Hence, robustness of bistability remains substantial even under significant parameter variation. Thus far, we have considered only the apoptotic activation of an individual cluster. To obtain the more biologically relevant cell-level activation, we must integrate over all clusters. In principle, this integration should account for intercluster transport as well as any intrinsic differences between clusters, e. g. , as due to spatial inhomogeneities. Here, however, we provide as demonstration only a very simple integration scheme. Specifically, we assumed that clusters are identical (apart from their parameter values, which are drawn randomly) and independent, and that FasL is homogeneous over the cell membrane. Then we can express the normalized cell activation as (6) where the subscript denotes reference to cluster. A characteristic cell-level hysteresis curve is shown in Figure 9. As is immediately evident, such integration is a smoothing operator, averaging over the sharp thresholds of each cluster. Thus, the cell-level signal may be graded even though its constituents are not. Note, however, that the lack of a sudden switch from low to high Fas signaling does not necessarily imply the same at the level of the caspases which ultimately govern cell death, as downstream components may possess switching behaviors [7], [8], [10], [11], [25]. Finally, we sought to outline protocols to experimentally discriminate our model against the prevailing crosslinking model [26], which we considered in a slightly simplified form [35]. To be precise, the crosslinking model that we used has the reactions (7a) (7b) (7c) where is FasL, is Fas, and is the complex FasL∶Fas for. With (8a) (8b) (8c) (8d) (8e) (continuing the notational convention that lowercase letters denote the concentrations of their uppercase counterparts), the steady-state solution under mass-action dynamics is (9a) (9b) (9c) where (10a) (10b) and (11a) (11b) are the total ligand and receptor concentrations, respectively. In analogy with our proposed model, hereafter called the cluster model, we used (12) as a measure of the apoptotic signal.
In this work, we showed through analysis of a mathematical model that receptor clustering can support bistability and hysteresis in apoptosis through a higher-order analogue of biologically observed Fas pair-stabilization [27]. Hence we add to the signal processing activities in which receptor clustering has been suggested to participate [37]–[39]. This bistability plays an important functional role by enabling robust threshold switching between life and death states. Significantly, our results indicate potential key roles for ligand trimerism [17] and receptor pre-association [30]–[32] and localization onto membrane lipid rafts [33]. Thus, we provide novel interpretations for these phenomena within the unified context of bistability. Our model suggests an additional cell death decision, supplementing those that have been studied previously [7], [8], [10], [11], [25]. Critically, the proposed decision is implemented upstream at the very death receptors that initially detect the death signal encoded by FasL. This decision is therefore apical in that it precedes all others in the system. Consequently, it operates independently of all intracellular components and so offers a general mechanism for bistability, even in cell lines with, for example, only feedforward caspase-activation networks [7], [9], [13], [14]. Thus, receptor cluster-activation may explain how an effective apoptotic decision is implemented in such cells. Moreover, this suggests a novel target for induced cell termination in the treatment of disease [1]. We believe that our model provides an attractive theory for the observed biology. Although unlikely to be correct in mechanistic detail, the model may nevertheless reflect reality at a qualitative level. The significance of our work hence lies in its capacity to guide future research. We therefore readily invite experiment, which can reveal the true nature of the molecular mechanisms involved. Given their structural and functional homology, similar investigations on other members of the tumor necrosis factor receptor family may also prove fruitful. Such work serves to further our understanding of the formation and mode of action of complex signaling platforms such as the DISC, which in this view may be considered the macromolecular aggregates of active Fas.
The rationale for the choices is presented in the text; here, we further defend these by noting that no new behaviors are introduced with or. The remaining parameter values were guided by the following considerations. Specifically, we required and due to the assumed stabilities of the receptor species; all other parameters were assumed to be close to. Within these constraints, parameters were selected to ensure that is of the correct order of magnitude [7], [9], [12], [13]. The baseline parameter values used were, , , , , , and. To analyze the effects of variability in the model parameters, parameter values were sampled from a log-normal distribution, characterized by a variation coefficient, defined as the ratio of the standard deviation to the median of the distribution. For the sensitivity analysis, parameter sets were drawn at; for the robustness analysis, parameter sets were drawn over; and for the cell-level integration, parameter sets were drawn at. All parameters were drawn about baseline median values. For each threshold-parameter pair, linear regression was performed on the threshold data against the parameter data, each normalized by reference values. For parameters, the reference is the baseline (median) value; for thresholds, the reference is the threshold computed at baseline parameters. The normalized sensitivity was defined as the slope of the linear regression. The cluster invariant was derived by considering at steady state, i. e. , with, and identifying. Similarly, the crosslinking invariant was derived by considering at steady state and identifying. For the full forms of the invariants, see Protocol S1. For the model discrimination computation, parameter sets were drawn from log-normal distributions with median at. The active Fas concentration was calculated for each parameter set for each model at baseline parameters (for the crosslinking model); for the cluster model, if bistability was observed, one of the stable values of was chosen at random. The invariant error was minimized using SLSQP with a lower bound of for all parameters. All calculations were performed with Sage 4. 5 [40], using NumPy/SciPy [41] for numerical computation and matplotlib [42] for data visualization. The Sage worksheet containing all computations is provided in the Supporting Information (Protocol S1) and can also be downloaded from http: //www. sagenb. org/home/pub/1224/ or http: //www. courant. nyu. edu/~ho/.
Many prominent diseases, most notably cancer, arise from an imbalance between the rates of cell growth and death in the body. This is often due to mutations that disrupt a cell death program called apoptosis. Here, we focus on the extrinsic pathway of apoptotic activation which is initiated upon detection of an external death signal, encoded by a death ligand, by its corresponding death receptor. Through the tools of mathematical analysis, we find that a novel model of death ligand-receptor interactions based on recent experimental data possesses the capacity for bistability. Consequently, the model supports threshold-like switching between unambiguous life and death states; intuitively, the defining characteristic of an effective cell death mechanism. We thus highlight the role of death receptors, the first component along the apoptotic pathway, in deciding cell fate. Furthermore, the model suggests an explanation for various biologically observed phenomena, including the trimeric character of the death ligand and the tendency for death receptors to colocalize, in terms of bistability. Our work hence informs the molecular basis of the apoptotic point-of-no-return, and may influence future drug therapies against cancer and other diseases.
cell biology/cellular death and stress responses
cell biology/cell signaling
mathematics
biophysics/biomacromolecule-ligand interactions
computational biology/systems biology
Bistability in Apoptosis by Receptor Clustering