|Home | About | Journals | Submit | Contact Us | Français|
Cell death plays an essential role in the development of tissues and organisms, the etiology of disease, and the responses of cells to therapeutic drugs. Here we review progress made over the last decade in using mathematical models and quantitative, often single-cell, data to study apoptosis. We discuss the delay that follows exposure of cells to prodeath stimuli, control of mitochondrial outer membrane permeabilization, switch-like activation of effector caspases, and variability in the timing and probability of death from one cell to the next. Finally, we discuss challenges facing the fields of biochemical modeling and systems pharmacology.
Apoptosis is a form of programmed cell death involving caspases, specialized cysteine proteases found in animal cells as inactive proenzymes (Fuentes-Prior and Salvesen, 2004). Dramatic progress has been made in recent years in identifying and determining the biochemical activities and cellular functions of biomolecules that regulate apoptosis and carry out its proteolytic program. However, current knowledge is largely qualitative and descriptive, and the complex circuits that integrate prosurvival and prodeath signals to control the fates of normal and diseased cells remain poorly understood. Successful creation of quantitative and predictive computational models of apoptosis would be significant from both basic research and clinical perspectives. From the standpoint of basic research, apoptosis is a stereotypical systems-level problem in which complex circuits involving graded and competing molecular signals determine binary life-death decisions at a single-cell level. Progress in modeling such decisions has had a significant impact on the small but growing field of mammalian systems biology. From a clinical perspective, diseases such as cancer involve disruption of the normal balance between cell proliferation and cell death, and anticancer drugs are thought to achieve their therapeutic effects by inducing apoptosis in cancer cells (Fadeel et al., 1999). However, it is difficult to anticipate whether a tumor cell will or will not be sensitive to a proapoptotic stimulus or drug based on general knowledge of apoptosis biochemistry because the importance of specific processes varies dramatically from one cell type to the next. Predictive, multifactorial, and context-sensitive computational models relevant to disease states will impact drug discovery and clinical care.
Apoptosis can be triggered by intrinsic and extrinsic stimuli. In intrinsic apoptosis, the death-inducing stimulus involves cellular damage or malfunction brought about by stress, ultraviolet (UV) or ionizing radiation, oncogene activation, toxin exposure, etc. (Kaufmann and Earnshaw, 2000). Extrinsic apoptosis is triggered by binding of extracellular ligands to specific transmembrane receptors, primarily members of the tumor necrosis factor receptor (TNFR) family (Kaufmann and Earnshaw, 2000). Receptor binding by TNF family ligands activates caspase-dependent pathways that are quite well understood in molecular terms. In general, extrinsic apoptosis has received more attention than intrinsic apoptosis from investigators seeking to develop mathematical models, but extrinsic and intrinsic apoptosis share many components and regulatory mechanisms.
The best studied inducers of extrinsic apoptosis are TNF-α, Fas ligand (FasL, also known as Apo-1/CD95 ligand), and TRAIL (TNF-related apoptosis-inducing ligand, also known as Apo2L; Figure 1A). Binding of these ligands to trimers of cognate receptors causes a conformational change that promotes assembly of death-inducing signaling complexes (DISCs) on receptor cytoplasmic tails (Gonzalvez and Ashkenazi, 2010). DISCs contain multiple adaptor proteins, such as TRADD and FADD, which recruit and promote the activation of initiator procaspases. The composition of the DISC differs from one type of death receptor to the next and also changes upon receptor internalization (Schutze et al., 2008). A remarkable feature of TNF-family receptors is that they activate both proapoptotic and prosurvival signaling cascades and the extent of cell death is determined in part by the balance between these competing signals. Prodeath processes are triggered by activation of initiator procaspases-8 and -10 at the DISC, a process that can be modulated by the catalytically inactive procaspase-8 homolog FLIP (Fuentes-Prior and Salvesen, 2004). Prosurvival processes are generally ascribed to activation of the NF-κB transcription factor, but other less well-understood processes are also involved, such as induction of the mitogen-activated protein kinase (MAPK) and Akt (protein kinase B) cascades (Falschlehner et al., 2007).
Initiator caspases recruited to the DISC directly cleave effector procaspases-3 and -7 generating active proteases (Fuentes-Prior and Salvesen, 2004). Effector caspases cleave essential structural proteins such as cytokeratins and nuclear lamins and also inhibitor of caspase-activated DNase (iCAD), which liberates the DNase (CAD) to digest chromosomal DNA and cause cell death. So-called “type I” apoptosis, which comprises a direct pathway of receptor→initiator caspases→effector caspases→death, is thought to be sufficient for death in certain cell types, but in most cell types apoptosis occurs by a “type II” pathway in which mitochondrial outer membrane permeabilization (MOMP) is a necessary precursor to effector caspase activation (Scaffidi et al., 1998). MOMP is triggered by the formation of pores in the mitochondrial membrane. Pore formation is controlled by the ~20 members of the Bcl-2 protein family, which can be roughly divided into four types: the “effectors” Bax and Bak whose oligomerization creates pores; “inhibitors” of Bax and Bak association such as Bcl-2, Mcl1, and BclxL; “activators” of Bax and Bak such as Bid and Bim; and “sensitizers” such as Bad, Bik, and Noxa that antagonize antiapoptotic Bcl-2-like proteins (Letai, 2008). In extrinsic apoptosis, initiator caspases that have been activated at the DISC cleave Bid into tBid, which in turn promotes a conformational change in Bax and Bak leading to oligomerization. Bax or Bak oligomers create pores in the mitochondrial outer membrane and promote cytoplasmic translocation of critical apoptosis regulators such as cytochrome c and Smac/Diablo, which normally reside in the space between the outer and inner mitochondrial membranes. MOMP does not occur until proapoptotic pore-forming proteins overwhelm antiapoptotic Bcl-2-like proteins (the so-called rheostat model) (Korsmeyer et al., 1993). Under most circumstances, MOMP is a sudden process that lasts a few minutes and marks the point of no return in the commitment to cell death (Chipuk et al., 2006; Tait et al., 2010). Once translocated to the cytosol, cytochrome c combines with Apaf-1 and caspase-9 to form the apoptosome, which cleaves and activates effector procaspases (Fuentes-Prior and Salvesen, 2004). XIAP associates with the catalytic pocket of active effector caspases-3 and -7 blocking protease activity and promoting their ubiquitin-dependent degradation. Binding of Smac to XIAP relieves this inhibition, allowing effector caspases to cleave their substrates and cause cell death (Fuentes-Prior and Salvesen, 2004).
In this Review, we describe how combining theoretical and computational approaches with live-cell imaging and quantitative biochemical analysis has provided new insight into mechanisms controlling the dynamics of extrinsic apoptosis. We start with a brief description of modeling concepts and methods relevant to apoptosis research. Next, we survey the recent literature. Modeling apoptosis, like quantitative analysis of mammalian signal transduction in general, is a field in its infancy fraught with many technical and conceptual challenges. Thus, only a subset of the known biochemistry of extrinsic apoptosis has been subjected to computational analysis, and this analysis has been performed only in a few cell lines. Key questions, such as differences between normal and transformed cells, have not yet been addressed in terms amenable to modeling. This Review, therefore, focuses on the subset of questions for which modeling has provided new insight (Figure 2). These include: (1) How is all-or-none control over effector caspase activity achieved? (2) How are activated effector caspases inhibited during the pre-MOMP delay while initiator caspase activity rises? (3) How do prosurvival and prodeath signals interact to determine if and when MOMP occurs? (4) What causes cell-to-cell variation in the timing and probability of apoptosis? We close this Review with an evaluation of current and emerging methods and future prospects. Readers interested in a more thorough description of the biology of extrinsic apoptosis are referred to several excellent reviews (Fuentes-Prior and Salvesen, 2004; Gonzalvez and Ashkenazi, 2010; Hengartner, 2000) and to Douglas Green's new book Means to an End: Apoptosis and Other Cell Death Mechanisms (Green, 2011).
The term “model” is used in a variety of fields in the natural and applied sciences to describe a mathematical or computational representation of a physical system. In molecular biology, the term usually refers to a “word model” or narrative description accompanied by a diagram, although it can also refer to a cell line or genetically engineered mouse that recapitulates aspects of a human disease. In this Review, we restrict use of the term “model” to describe an executable set of rules or equations in mathematical form. We are primarily interested in models that are built and tested using detailed cellular or biochemical experiments. Models of cellular biochemistry can be based on different mathematical formalisms, from Boolean logic to differential equations, depending on the degree of detail and the scope of the modeling effort. Most models of apoptosis have been encoded using ordinary differential equations (ODEs), which describe the evolution of a system in continuous time. ODEs are the mathematical representation of mass action kinetics, the familiar biochemical approximation in which rates of reaction are proportional to the concentrations of reactants (Figure 1B) (Chen et al., 2010). Diffusion, spatial gradients, or transport can be modeled explicitly using partial differential equations (PDEs), which represent biochemical systems in continuous time and space. For example, Rehm et al. (2009) used PDEs to model the spread of mitochondrial permeabilization through a cell following an initial, localized MOMP event. Using sets of differential equations it is possible to encode a complex network of interacting biochemical reactions and then study network dynamics under the assumption that protein concentrations and reaction rates can be estimated from experimental data. Differential equation models often increase rapidly in complexity as species are added, as each new protein can give rise to a large number of model species differing in location, binding state, and degree of posttranslational modification. This problem has effectively limited data-dependent ODE/PDE models to fewer than ~20 gene products (and on the order of 50–100 model species), although efforts are underway to increase this limit.
In addition to differential equations, several other formalisms have been used to model apoptosis. Stochastic models make it possible to represent reactions as processes that are discrete and random, rather than continuous and deterministic. Stochastic models are advantageous when the number of individual reactants of any species is small (typically fewer than ~100) or reaction rates very slow (Zheng and Ross, 1991). In these cases, a Monte Carlo procedure is used to represent the probabilistic nature of collisions and reactions among individual molecules (Gillespie, 1977). For example, stochastic cellular automata have been used to model the movement of molecules on the mitochondrial outer membrane (Chen et al., 2007). When sufficient time-resolved quantitative data are lacking, a less precise modeling framework is usually advantageous, and logic-based models have proven particularly popular. Boolean models, for example, are discrete two-state logical models in which each node in a network is represented as a simple on/off switch. Boolean models have been used to represent the interplay among survival, necrosis, and apoptosis pathways and to predict the likelihood that each phenotype would result following changes in the levels of regulatory proteins (Calzone et al., 2010). However, the more qualitative and phenomenological the modeling framework, the less mechanistic the insight.
Regardless of modeling framework, a trade-off exists between model tractability and model detail or scope. The inclusion of more species makes it possible to analyze biochemical processes in greater detail or to represent the operation of large networks involving many gene products, but larger models are more difficult to constrain with experimental data, and excess detail can mask underlying regulatory mechanisms. A Jorge Luis Borges story comes to mind in which the art of cartography achieved such a perfection of detail that cartographers built a map of their empire with 1:1 correspondence to the empire itself, rendering the map useless (Borges and Hurley, 1999). On the other hand, although small models have the advantage of relative simplicity and even analytical tractability (i.e., capable of being solved exactly without simulation), they run the risk of grossly simplifying the underlying biochemistry and of including an insufficient number of regulatory processes. As yet, no clear principles exist to guide decisions about model scope and complexity, and most studies remain constrained by the relative immaturity of modeling software and a paucity of experimental data.
Estimating values for rate constants and initial protein concentrations (the parameters in differential equation models) remains extremely challenging both computationally and experimentally. Each reaction in an ODE model is associated with one or more “initial conditions” (the concentrations of reactants at time zero) and rate constants, usually a forward and reverse rate constant. Some of these parameters are available in the literature, typically from in vitro biochemical experiments, and these values may hold true in the context of a cell. In many other cases, however, no estimates of rate constants are available and parameters must be estimated directly from experiments (Chen et al., 2010). In addition, protein concentrations vary from cell type to cell type and should be measured directly in the cell type under investigation, although this is often not done because it is time consuming. The estimation of unknown parameter values based on data (typically, time-dependent changes in the abundance or localization of proteins in the model) is called model calibration, model training, or model fitting. Almost all realistic models of biological systems are too large for all parameters to be fully constrained by experimental data, and the models are therefore “nonidentifiable.” Thus far, the process of model calibration has been approached rather informally, but more rigorous approaches are in development (e.g., Kim et al., 2010). Careful analysis is expected to confirm the common-sense view that solid conclusions can be reached even in the case of partial knowledge.
Modeling biological processes requires the collection and analysis of quantitative experimental data. An ODE model, which assumes that each compartment is well mixed, necessarily represents a single cell, and calibrating and testing ODE models therefore require collecting data on single cells over time. However, live-cell imaging experiments usually rely on genetically modified cell lines carrying fluorescent reporters. Creating these lines is relatively time-consuming, and the extent of multiplexing is limited by phototoxicity and the availability of noninterfering fluorophores. It is not always clear that an engineered reporter correctly represents the activity or state of modification of endogenous proteins (see, for example, discrepancies regarding initiator caspase activity reporters, discussed below; Albeck et al., 2008a; Hellwig et al., 2008; Hellwig et al., 2010). Flow cytometry, immunofluorescence, and single-cell PCR are also effective means to assay single cells, and biochemical experiments (immunoblotting or ELISAs for example) performed on populations of cells remain essential for quantitative biology. Although rarely addressed, effective integration of data arising from multiple measurement methods is an area in which computational models are likely to play a key role (Albeck et al., 2006).
The construction and parameterization of even a well designed model do not lead directly to a better understanding of the system—model analysis is required. The dependence of the system on parameter values is of particular interest and can be approached using sensitivity analysis. Sensitivity analysis involves systematically varying parameters (initial conditions or rate constants) while monitoring the consequences for model output (the time at which a cell undergoes apoptosis, for example). Sensitivity analysis reveals which outputs are sensitive to variation in which parameters and can be viewed as the computational equivalent of experiments that knock down or overexpress proteins while monitoring phenotype. For example, Hua et al. (2005) created an ODE model of Fas signaling and performed sensitivity analysis by varying the initial concentration of each protein species 10- or 100-fold above or below a baseline value. Using the half-time of caspase-3 activation as an output, they predicted (and confirmed experimentally) that increases but not decreases in Bcl-2 levels would alter sensitivity to FasL. From a practical perspective, sensitive parameters must be estimated with particular care if a model is to be reliable, but from a biological perspective, they represent possible means of regulation. Points in a network that exhibit extreme sensitivity to small perturbations are often referred to as “fragile” (the converse of “robust”), and considerable interest exists in the idea that fragility analysis, a concept borrowed from control theory, might be applied to biological pathways. In this view, fragile points might identify processes frequently mutated in disease or potentially modifiable using therapeutic drugs (Luan et al., 2007).
Stability analysis is another commonly used method of model analysis. Some models of biochemical networks have the interesting property of converging at equilibrium to a small set of stable states known as fixed points, where the rate of change in the concentrations of all model species is zero. Identification and characterization of fixed points can provide valuable insight into the dynamics of a system, its responses to perturbation, and the nature of regulatory mechanisms. Of particular interest in biology is bistability, a property in which a system of equations has two stable fixed points separated by an unstable fixed point. Bistability has obvious appeal in the case of apoptosis, in which cells are either alive or dead, and has been proposed to underlie a variety of binary fate decisions such as maturation of Xenopus oocytes (Ferrell and Machleder, 1998) and lactose utilization in E. coli (Ozbudak et al., 2004). From the perspective of control, many bistable systems have two valuable properties: (1) they are insensitive to minor perturbations because the system is “attracted” to the nearest stable state (in apoptosis, a bistable system would be resistant to spontaneous activation of proapoptotic proteins, for example), and (2) they exhibit “all-or-none” transitions from one stable state to another in response to small changes in the level of a key regulatory input (a property known in biochemistry as “ultrasensitivity”). Bistable processes often exhibit hysteresis (path dependence): once in the on state, they do not readily slip back to off. It is often assumed that the regulatory machinery for apoptosis must be bistable in the mathematical sense with one equilibrium state corresponding to caspases off and “alive” and the other to caspases on and “dead” (Figure 3A). Although bistability remains the favorite framework for thinking about the switch between life and death, bistability is not strictly necessary for a switch-like transition between two distinct states (Albeck et al., 2008b). A monostable system in which the landscape changes through time can create a temporal switch between two states; in this case, the change in the landscape involves the creation, destruction, or translocation of precisely those proteins (caspases, cytochrome c, etc.) that are known to regulate apoptosis. In this regard, it should be noted that the “sharpness” of a switch in a conventional bistable system refers to the steepness of the dose-response curve (to a change in the concentration of a regulatory protein, for example), not necessarily sharpness in time. In contrast, the “all-or-nothing” switch observed by time-lapse microscopy of cells undergoing apoptosis refers to a switch from alive to dead that is sharp in a temporal sense. These considerations do not imply that the biochemical pathways controlling apoptosis are not bistable systems, but rather that bistability is not necessary a priori.
The first model of extrinsic apoptosis was published a decade ago and set the stage for subsequent work in the field. Fussenegger et al. (2000) used emerging understanding of MOMP and caspase activation by death receptors to assemble a simple ODE model. By increasing or decreasing the levels of pairs of proteins in the model, the authors determined which combinations promoted or blocked effector caspase activation, thereby providing insight into ratiometric control over cell death by caspase-3 and XIAP (Fussenegger et al., 2000). At the same time, the development of fluorescent reporters for MOMP and caspase substrate cleavage allowed several groups to collect data on the dynamics of apoptosis in single cells. These data showed that following exposure to inducers of either intrinsic or extrinsic apoptosis (UV light, actinomycin D, staurosporine, or TNF), cells wait for several hours before initiating a rapid chain of events that triggers MOMP and activates effector caspases (Goldstein et al., 2000a, 2000b; Tyas et al., 2000). This contrasts with data obtained by western blotting and other population-average biochemical assays that suggested that MOMP and caspase activation occur gradually over a period of several hours. The two types of data can be reconciled by noting that apoptosis is sudden and switch like in individual cells, but that it takes place at different times in different cells (Figure 3B) (Goldstein et al., 2000a; Goldstein et al., 2000b; Tyas et al., 2000).
Goldstein et al. (2000b) used time-lapse imaging of cytochrome c translocation to obtain the first data on the kinetics of MOMP. They observed the time between proapoptotic insult and MOMP to vary depending on the type and strength of the stimulus (ranging from 4–20 hr following exposure to the pan-specific kinase inhibitor staurosporine and 9–17 hr following exposure to UV light), but the rate and extent of cytochrome c release were constant, taking ~5 min to reach completion. Further understanding of the link between MOMP and caspase activation was made possible by the development of intramolecular Förster resonance energy transfer (FRET) reporters for caspase-mediated proteolysis. The first FRET reporters for monitoring caspase activity by time-lapse microscopy linked cyan fluorescent protein (CFP) to yellow fluorescent protein (YFP) using a polypeptide linker containing the amino acid sequence DEVD, a substrate for caspase-3 (CFP-DEVD-YFP) (Rehm et al., 2002; Tyas et al., 2000). Prior to reporter cleavage, CFP lies in close proximity to YFP, causing FRET between the two fluorescent proteins and reducing CFP emission. Following cleavage of the DEVD-containing linker, the efficiency of FRET drops dramatically, increasing the CFP to YFP fluorescence ratio. Time-lapse imaging of cells expressing CFP-DEVD-YFP revealed that caspase-3 is also activated rapidly, taking <15 min to reach completion (Rehm et al., 2002; Tyas et al., 2000).
Rehm et al. (2002) asked whether the cleavage kinetics of effector caspase substrates depended on the identity or strength of the apoptotic stimulus. Like Goldstein et al. (2000b), Rehm and colleagues observed the delay between exposure to a prodeath stimulus and the onset of effector caspase activation to vary from cell to cell. They also noted that the average delay varied with the dose and identity of the proapoptotic stimulus (3 ng/ml or 200 ng/ml TNF, 3 μM staurosporine, 10 μM etoposide), but the kinetics of reporter cleavage did not. The authors developed a quantitative description of these data, showing that caspase activation in individual cells fits a sigmoidal Boltzmann equation in which the lag time is dose and stimulus dependent, but cleavage kinetics are dose invariant (Figures 3C and 3D) (Rehm et al., 2002). Subsequent multiplex imaging of MOMP and effector caspase reporters in single cells showed that MOMP precedes effector caspase activation by ~10 min (Rehm et al., 2003). In electromechanical terminology, the regulation of MOMP and effector caspase activity constitutes a variable delay, snap-action switch.
Intrigued by the idea that a switch is central to the regulation of apoptosis, several groups have attempted to understand how such a switch might arise, based on models in which bistability is assumed as a design principle. Eissing et al. (2004) created an 8-equation ODE model of apoptosis in a type I cell that included activation of caspase-8, consequent cleavage and activation of caspase-3, inhibition and degradation of activated caspase-3 by XIAP, and activation of residual caspase-8 by activated caspase-3 in a feedback loop. The small size of the model made it possible to identify stable states analytically, and the authors found that adding a mechanism to inhibit active caspase-8 (via the protein Bar) was necessary to ensure bistability at the level of effector caspase activity (Eissing et al., 2004). A subsequent modeling study that examined how cells would resist spontaneous procaspase-8 activation argued against a major role for Bar, however (Wurstle et al., 2010). Legewie et al. (2006) created a 13-ODE model that described activation of caspase-9 by Apaf-1, consequent activation of caspase-3, and inhibition of caspases by XIAP. The authors identified an “implicit” or hidden positive feedback loop as a key contributor to bistability; in this loop, caspase-3 promotes its own activation by sequestering XIAP away from caspase-9, allowing caspase-9 to cleave additional procaspase-3 (Legewie et al., 2006). Bagci et al. (2006) built models of varying complexity (the largest being 31 ODEs) centered on apoptosome formation and caspase-3 activation and concluded that cooperativity in the formation of the apoptosome was a key element for ensuring bistability (Bagci et al., 2006). Chen et al. (2007) constructed both an ODE model and a stochastic cellular automaton model to examine the potential for interactions among Bcl-2-family members to generate bistability at MOMP. These models included activation of Bax by an activator such as tBid, inhibition of the activator and Bax by Bcl-2, and displacement of the activator in the activator-Bcl-2 complex by Bax. This description of Bax and Bcl-2 also encoded an implicit positive feedback loop in which freed activator could bind more Bax, leading to bistability in pore formation. Addition of cooperativity in Bax multimerization resulted in a model encoding a one-way (as opposed to bidirectional or hysteretic) switch (Chen et al., 2007). In a corroborating study that used flow cytometry, an antibody against activated Bax (clone 6A7) revealed a bimodal distribution in the staining of HeLa cells treated with 400 nM staurosporine for 6 hr (Sun et al., 2010). However, antibody staining is unreliable in dying/dead cells, so proving the point will require showing bimodality in Bax activation in cells that have not yet undergone effector caspase activation. Most recently, Ho and Harrington (2010) built a small ODE model in which FasL acts as a clustering agent for Fas receptors. The reactions described spontaneous receptor opening and closing, constitutive destabilization of open receptors, and ligand-independent and -dependent stabilization of receptor clusters. Analytical methods showed the system to exhibit reversible bistability (hysteresis) at low receptor concentrations but irreversible bistability at higher local receptor densities (Ho and Harrington, 2010). In summary, this set of papers reveals that almost every point in the apoptosis pathway has the potential to generate bistability in the mathematical sense. However, many of the papers were written in an era in which it was not yet common for mathematical modeling to be combined with quantitative experimentation in a single manuscript. The results of simulation were compared to data from the literature, but proposed regulatory mechanisms were not confirmed using RNAi or other perturbation-based experiments.
Whereas the first generation of apoptosis models focused on specific steps in the process of cell death (MOMP, apoptosome formation, etc.), Albeck et al. (2008b) built a model that spanned the entire pathway of extrinsic apoptosis from ligand binding to cleavage of effector caspase substrates, albeit in simplified form. A model comprising 58 differential equations was sufficient to capture the essence of TRAIL-receptor binding, cleavage of initiator and effector caspases, initiation of MOMP, release of Smac and cytochrome c into the cytoplasm, and finally caspase-3 activation and substrate cleavage. The model was trained against experimental data that included live-cell microscopy, immunoblotting, and flow cytometry in wild-type HeLa cells or cells perturbed by protein overexpression or RNAi-mediated protein depletion. Model analysis and experiments confirmed earlier evidence that MOMP is the point in receptor-mediated apoptosis at which upstream signals are transformed into an all-or-none snap-action signal (Goldstein et al., 2005, 2000b; Madesh et al., 2002; Rehm et al., 2003; von Ahsen et al., 2000). To understand how this switch might work in molecular terms, Albeck et al. (2008b) analyzed a series of models of increasing complexity and biochemical realism that linked tBid cleavage by initiator caspases to Smac/cytochrome c release. The performance of each model was analyzed for its ability to create a variable-delay, snap-action switch. A useful insight was that the “rheostat model” (Korsmeyer et al., 1993), in which it was postulated that MOMP is triggered when levels of active Bax/Bak exceed those of Bcl-2/BclxL, only functioned in its simplest form as a switch if Bax and Bcl-2 were assumed to associate irreversibly at a rate faster than diffusion. In contrast, snap-action switching emerged naturally from the biochemistry of Bax and Bcl-2 if more complex reaction topologies were assumed; these included slow activation of Bax by tBid, partitioning of reactants into cytosolic and mitochondrial compartments, and a requirement for Bax multimerization. Rapid and complete translocation of Smac and cytochrome c was ensured in part by the favorable kinetics of moving proteins down a steep concentration gradient from the mitochondrion (where they are abundant) to the cytosol (where they are initially absent). Despite the apparent success of the Albeck et al. (2008b) model, it is important to realize that it involves a simple picture of DISC formation as well as a simplified version of MOMP that lacks the multiplicity of positive- and negative-acting Bcl-2-family members present in real cells.
Better understanding of caspase substrate specificity (Luo et al., 2003; Stennicke et al., 2000; Thornberry et al., 1997) along with a direct comparison of CFP-DEVD-YFP cleavage kinetics with those of endogenous substrates (Albeck et al., 2008a) made clear that the CFP-DEVD-YFP biosensor is processed by both effector and initiator caspases. Changing the biosensor linker to DEVDR made it 20-fold more selective for caspase-3 relative to caspase-8 (Albeck et al., 2008a), and changing the cleavage recognition site to IETD resulted in a FRET reporter selective for initiator caspases (Luo et al., 2003). Combining this selective effector caspase reporter with a MOMP reporter showed that effector caspase activity is negligible during the pre-MOMP delay (Albeck et al., 2008a); this had correctly been assumed to be true by Rehm et al. (2002), despite the use of a less specific CFP-DEVD-YFP reporter. In contrast, initiator caspases are active throughout the pre-MOMP delay (Albeck et al., 2008a; Hellwig et al., 2008), and their substrates Bid and procaspase-3 accumulate in cleaved form. Caspase-3 is a very potent enzyme, and model-based simulation and experiments suggest that ~400 active molecules are sufficient to cleave 106–107 molecules of cellular substrate within several hours (Albeck et al., 2008a). However, during the pre-MOMP delay, no processing of effector caspase substrates can be detected using live-cell FRET reporters or flow cytometry (c.f. Figures 2A and 2B). This raises the interesting question: How are effector caspases maintained in an off state despite being continually processed by initiator caspases from a zymogen into a cleaved and potentially active form?
One mechanism for keeping processed effector caspases “off” is binding of XIAP to the catalytic cleft of caspase-3. This tight interaction (~1 nM) might seem sufficient to hold caspase-3 in check, but modeling shows that a >100-fold molar excess of XIAP over caspase-3 would be required to ensure effective inhibition of caspase-3 proteolytic activity over the course of a typical 2–6 hr pre-MOMP delay. The requirement for such a large excess of XIAP over caspase-3 arises because competitive inhibition is reversible whereas substrate cleavage is not and because substrates, which are abundant, are in competition with XIAP for access to the caspase catalytic site. As XIAP and caspase-3 are present at roughly equal concentrations in HeLa cells, simple competitive inhibition cannot be the sole inhibitory mechanism. XIAP is an E3 ligase able to promote ubiquitination and degradation of caspase-3, and simulation suggests that a combination of competitive inhibition and caspase degradation would constitute an effective means of regulation (Albeck et al., 2008a). Confirming these predictions, depletion of XIAP by RNAi or pharmacological inhibition of the proteasome was observed to cause effector caspase activation prior to MOMP (Albeck et al., 2008a). Deletion of XIAP in the mouse or truncation of the ubiquitination-promoting RING domain also caused elevated caspase-3 activity and sensitivity to apoptosis (Schile et al., 2008), demonstrating a critical role for XIAP-mediated ubiquitination of caspase-3 in vivo. The pre-MOMP delay evidently constitutes a “latent” death state in which effector procaspases are actively processed by initiators but are held in check by XIAP until Smac is released during MOMP. The reasoning that led to this conclusion illustrates the value of making models explicit and analyzing them computationally: a biochemical mechanism that seems adequate on its face (Bcl-2-Bax binding in the rheostat model or competitive inhibition of C3 by XIAP during the pre-MOMP delay) proves insufficient when actual protein levels and rates of reaction are taken into account. In this sense, quantitative analysis can fundamentally change our qualitative understanding of a regulatory mechanism. It should be noted, however, that current models of receptor-mediated apoptosis in type II cells cannot completely restrain pre-MOMP caspase-3 activity when experimentally measured procaspase-3 and XIAP concentrations are used. Although XIAP-mediated degradation of active caspase-3 is necessary, raising this degradation rate too much compromises the switch-like activation of effector caspase substrate cleavage post-MOMP. Reconciliation of all experimental observations awaits the development of more sophisticated and complete models.
If XIAP is partially depleted by RNAi and MOMP is blocked by overexpression of Bcl-2, a sublethal level of effector caspase activity is generated and effector caspase substrates are only partially processed; moreover incomplete cleavage of caspase-3 substrates does not necessarily cause cell death (at least in HeLa cells) (Albeck et al., 2008a). Modeling and experiments with XIAP overexpression suggest three possible outcomes depending on XIAP levels: with [XIAP] < 0.15 μM, effector caspase substrate cleavage is complete; at [XIAP] > 0.30 μM, cleavage is fully inhibited; and at intermediate XIAP concentrations, slow submaximal effector caspase substrate cleavage occurs (Figure 3E) (Rehm et al., 2006). Thus, alteration of XIAP levels disrupts normal switch-like control over effector caspase activation and interferes with the normal link between caspase activation and cell killing. Activation of CAD in the absence of cell death is expected to be particularly problematic since it has the potential to cause genomic instability (Lovric and Hawkins, 2010) and has been proposed to be the trigger of the chromosomal translocations observed in some leukemias (Betti et al., 2005; Vaughan et al., 2002, 2005; Villalobos et al., 2006).
The role of XIAP in restraining caspase-3 in the absence of MOMP makes it a central factor in controlling type I versus type II apoptosis. Jost et al. (2009) observed that inhibition of XIAP function by gene targeting or a Smac mimetic drug caused type II cells to adopt a type I phenotype. Bid deficiency protected hepatocytes and pancreatic β cells from FasL-induced apoptosis (fulfilling the definition of mitochondria-dependent type II death), but concomitant loss of XIAP (in Bid−/− XIAP−/− mice) restored FasL sensitivity, thereby demonstrating a switch to type I behavior (Jost et al., 2009). Type I cells are defined as not requiring MOMP for apoptosis, but blockade of the mitochondrial pathway via Bid depletion or Bcl-2 overexpression in type I cells has been observed to reduce effector caspase activity and to increase the number of cells surviving TRAIL exposure (Maas et al., 2010). Both type I and type II pathways, therefore, appear to depend to a greater or lesser extent on the mitochondrial pathway, either for regulating XIAP and activating effector caspases or for killing cells by disrupting essential mitochondrial functions.
Apoptosis proceeds at different rates in different cells, even among members of a clonal population. Some cells die within 45 min of exposure to FasL or TRAIL, whereas other cells in the same dish wait 12 hr or more. A simple way to conceptualize control over the timing of apoptosis in single cells is that the level of active receptor determines the amount of active caspase-8/10, which sets the rate of tBid cleavage and, thus, the rate of approach to a threshold that must be overcome for MOMP to occur (Figure 2C). The height of this threshold is set by the relative levels of competing pro- and antiapoptotic Bcl-2-family proteins (Chipuk and Green, 2008). We discuss below recent advances in our understanding of the MOMP threshold and return later to the determinants of the rate of approach to the threshold. Using fluorescent measurements in a purified in vitro system, Lovell et al. (2008) simultaneously measured the rates of three reactions leading to pore formation and determined the following order of events. First tBid binds rapidly to mitochondrial membranes where tBid and Bax interact, promoting insertion of Bax into the membrane, a rate-limiting step. Bax then oligomerizes to form pores, and membranes become permeable. In vitro, Bax oligomerization continues even after membranes are permeabilized (Lovell et al., 2008). In cell culture, Bax multimerization is first detected immediately prior to MOMP and then continues for at least 30 min, ultimately generating many more Bax puncta or pores than the number required for MOMP (Albeck et al., 2008b; Dussmann et al., 2010). Formation of the first observable Bax (or Bak) puncta correlates temporally and spatially with the first subset of mitochondria to undergo MOMP. Pore formation and MOMP then spread through the cell as a wave with a velocity of ~0.6 μm/s, a process that has been modeled using a PDE network (Rehm et al., 2009). The process of pore formation proceeds more rapidly at higher doses of TRAIL, presumably due to an increased rate of procaspase-8 activation (Rehm et al., 2009). However, it has recently been observed that in a subset of HeLa cells, MCF-7 cells, and murine embryonic fibroblasts, some mitochondria fail to undergo MOMP in response to diverse proapoptotic stimuli (actinomycin D, UV, staurosporine, or TNF). The subset of mitochondria that remain intact fail to accumulate GFP-Bax puncta but undergo complete MOMP when treated with the Bcl-2 antagonist and investigational therapeutic ABT-737, suggesting that resistance of mitochondria to MOMP lies at the point of Bax/Bak activation (Tait et al., 2010). These findings suggest that mitochondria in a single cell differ from each other with respect to their sensitivities to proapoptotic stimuli and that MOMP might not always be an all-or-none event at the single-cell level (Tait et al., 2010).
Time-lapse imaging of initiator caspase and MOMP reporters shows that the height of the MOMP threshold varies from cell to cell. MOMP is triggered following cleavage of ~10% of a reporter carrying one IETD recognition site (Hellwig et al., 2008, 2010) or 30%–60% of a reporter carrying two recognition sites (Albeck et al., 2008a, 2008b; Spencer et al., 2009). Variation in the height of the MOMP threshold from cell to cell (presumably arising from variation in the levels of Bcl-2-family proteins, see below) can most easily be resolved using the sensitized dual-IETD reporter (Figure 2C) and contributes ~20% of the total variability in the time of death among HeLa cells in clonal population exposed to 10 ng/ml TRAIL (Spencer et al., 2009). The remaining 80% of the variability appears to reflect differences in the rate of Bid cleavage, although these percentages are expected to change with stimulus and cell type. However, the precise dynamics of Bid cleavage have recently been thrown into some doubt: a FRET reporter containing full-length Bid rather than an artificial IETD caspase recognition site exhibits minimal cleavage prior to MOMP (Hellwig et al., 2010). One explanation for this discrepancy is that IETD-only reporters might be overly sensitive and not reflect the kinetics of endogenous substrate cleavage. In this view, cleavage of Bid by initiator caspases is subject to additional forms of regulation so that tBid does not accumulate until just before MOMP (Hellwig et al., 2010). Conversely, the Bid-containing FRET reporter might simply be insufficiently sensitive, and levels of tBid required for MOMP (estimated to be <3% of the total Bid pool) might be below the level of detection. In this view, IETD-only reporters conveniently amplify a signal that would otherwise be undetectable. Resolving this question will require careful comparison of reporter constructs with endogenous proteins, which will itself depend on the availability of antibodies that can distinguish different caspase-8 substrates. It seems likely that carefully calibrated models will also help with data integration.
The rate at which initiator caspase substrates are cleaved varies from cell to cell (Figure 2C). The current view is that the strength of receptor signaling and the amount of active DISC control the rate of initiator caspase substrate cleavage and thus the rate of approach to the MOMP threshold, with lower levels of prodeath stimulus leading to slower Bid cleavage and slower onset of apoptosis. Models of DISC formation in FasL-treated cells have questioned whether apoptosis simply slows down with decreasing ligand concentrations (a continuous decrease), or whether there is a minimum ligand:receptor ratio needed for induction of apoptosis (a threshold; Bentele et al., 2004). Modeling predicted that below a critical ligand:receptor ratio of 1:100, apoptosis is completely blocked due to the presence of the inhibitory DISC component c-FLIP. Above the critical threshold, c-FLIP is insufficient to block all DISC activity prior to the formation of active caspase-8. A follow-up study refined this view by showing that active DISC is formed at concentrations of a receptor crosslinking antibody (anti-APO-1, which activates Fas receptors) below a critical threshold. However, because c-FLIP has a higher affinity than procaspase-8 to the few DISCs that are formed, activation of caspase-8 is effectively inhibited (Lavrik et al., 2007). Continuing this line of reasoning, Fricker et al. (2010) used modeling, biochemical assays, and live-cell imaging to explore how levels of c-FLIP iso-forms determine sensitivity to Fas signaling. Although c-FLIPS/R is well established as an inhibitor of Fas-mediated apoptosis, the role of c-FLIPL has been controversial because it plays both pro- and antiapoptotic roles depending on expression level. Model analysis suggested that the effects of c-FLIPL on activation of procaspase-8 vary with FasL levels: relative to cells having endogenous levels of c-FLIPL, 20-fold overexpression of c-FLIPL blocks cell death when FasL levels are low but accelerates death when FasL levels are high. However, even at high FasL levels, a further increase in c-FLIPL concentration inhibits procaspase-8 processing and decreases the extent of cell death. Models can be quite helpful in exploring these sorts of quantitative relationships. One explanation supported by model analysis involves the fact that c-FLIPL has higher affinity for DISCs than procaspase-8 but that procaspase-8 is present in cells in substantial molar excess. At low FasL levels, few DISCs are formed, relative affinities dominate, and the ratio of c-FLIPL to caspase-8 at DISCs is high. At high levels of FasL, many DISCs are formed, and the small number of c-FLIPL molecules (~300 per HeLa cell) is exhausted, allowing DISC-bound procaspase-8 to overwhelm c-FLIPL. Thus, subtle changes in the levels of c-FLIPL and FasL can change the timing and probability of death in nonlinear ways that can be understood only if the concentrations of interacting proteins are taken into account (Fricker et al., 2010).
An additional factor affecting the life-or-death fate of a cell exposed to death ligand is the interplay between prosurvival and proapoptotic pathways. The relative strength of these competing regulatory processes is also thought to be controlled by the composition of the DISC. Induced survival signaling has been largely attributed to NF-κB, and many negative regulators of apoptosis are known to be induced by NF-kB, including c-FLIP, BclxL, and members of the IAP family (Gonzalvez and Ashkenazi, 2010). It is not yet clear which of these factors is most important nor whether NF-κB-independent processes, such as MAPK signaling, also play important prosurvival functions. Whereas Bentele et al. (2004) and Fricker et al. (2010) focused on the presence or absence of prodeath signaling at the DISC, Lavrik et al. (2007) and Neumann et al. (2010) explicitly focused on the balance between prodeath versus prosurvival activities. Lavrik et al. (2007) demonstrated that Erk kinase is activated in response to anti-APO-1 over a wide range of doses, even in the presence of a pan-caspase inhibitor. This work showed that survival signaling occurred in parallel with death signaling, but it was not clear how the survival signal was initiated. Neumann et al. (2010) built and tested an ODE model of Fas-mediated apoptosis with a postulated link between apoptosis and survival pathways in which c-FLIPL is cleaved by caspase-8 into p43-FLIP, which binds and activates IκB kinase (IKK). IKK then phosphorylates and inhibits IκB, a negative regulator of NF-κB, leading to induction of NF-κB-mediated transcription. Simulation and experiments suggested that both proapoptotic (caspase-8 dependent) and prosurvival (NF-κB-dependent) pathways are activated in parallel and that a subtle balance between c-FLIPL and initiator caspase levels determines which one predominates. By visualizing the levels of c-FLIPL and procaspase-8 on a parameter landscape (Figure 3F), the authors showed that c-FLIPL can disable, promote, or inhibit NF-κB activation depending on whether the level of c-FLIPL is low, intermediate, or high. This effect arises because high levels of c-FLIPL prevent caspase-8-mediated processing of c-FLIPL into the IKK-binding p43-FLIP form (Neumann et al., 2010).
Future work on the topic of induced survival signaling would benefit from single-cell measurements combining reporters for NF-κB target gene expression (Nelson et al., 2004) and reporters of initiator caspase activity so as to capture feedback. Both the apoptosis and the NF-κB fields have a tradition of utilizing live-cell imaging and mathematical modeling (e.g., Ashall et al., 2009; Hoffmann et al., 2002; Lee et al., 2009), and it would be valuable to combine models of both processes. This would lead to better understanding of competing prosurvival and pro-death death processes in different cell types.
Individual cells differ widely in their responses to apoptotic stimuli (Figure 2D). Potential sources of cell-to-cell variability in the timing and probability of apoptosis include genetic or epigenetic differences, differences in cell-cycle phase, stochastic fluctuations in biochemical reactions, and natural variation in protein concentrations. To distinguish among these possibilities, three independent groups followed dividing cells using time-lapse microscopy and compared the timing and probability of apoptosis in sister cells and in randomly selected pairs of cells (Bhola and Simon, 2009; Rehm et al., 2009; Spencer et al., 2009). At a dose of TRAIL sufficient to induce apoptosis in half of the cells, the probability of death was observed to be highly correlated between sisters, as was the time at which cells died. Correlation in death time among sister cells has been observed in a variety of cell types (HeLa, MCF-10A, NIH 3T3, HT1080, and murine embryonic fibroblasts) following exposure to a variety of apoptosis-inducing agents (TRAIL, TNF-α, staurosporine, and etoposide). In contrast, randomly selected cells were found to be uncorrelated, and no obvious correlation with cell-cycle phase or with position in the dish could be detected (Bhola and Simon, 2009; Spencer et al., 2009), although the way in which these experiments were performed does not rule out some contribution from cell-cycle state (Rehm et al., 2009). Importantly, the degree of similarity between sisters fell as the time since cell division increased so that within one to two generations, sisters were no more correlated than randomly chosen pairs of cells. This transient heritability in timing of death argues against a genetic or epigenetic explanation for cell-to-cell variability in apoptosis, as genetic and epigenetic differences tend to be stable over much longer timescales. The initial correlation between sister cells also rules out a significant role for stochasticity in the reactions that regulate caspase activation, a conclusion supported by simulation (Eissing et al., 2005). Transient heritability in the timing and probability of death data are most consistent with an explanation rooted in natural cell-to-cell variation in the levels or activities of proteins among genetically identical cells. Sister cells are known to inherit similar levels of relatively abundant biomolecules during cytokinesis, but levels then diverge due to random fluctuations in protein synthesis and degradation (Sigal et al., 2006). In support of this, experiments with cycloheximide show that the rate of sister cell decorrelation is highly sensitive to the rate of protein synthesis (Spencer et al., 2009).
Over the last decade, modeling and experimentation in bacteria, yeast, and more recently in mammalian cells, have provided a mechanistic framework for understanding stochastic variation (“noise”) in rates of transcription and translation (Raj and van Oudenaarden, 2008). The number of transcriptional initiation complexes on any single gene is small (potentially as small as 1–2), and the probability that a transcript will be created in any time interval is therefore highly stochastic. Fluctuations in mRNA levels result in fluctuating rates of protein synthesis. With short-lived or low-copy-number proteins, this can cause large fluctuations in protein levels, whereas with relatively abundant proteins, such as those controlling apoptosis, the most significant effect is that different cells contain different concentrations of each protein, and thus unique proteomes. Current models predict that the distribution of concentrations across a population of cells should be long-tailed, following a log-normal or gamma distribution (Friedman et al., 2006; Krishna et al., 2005). In the case of proteins controlling apoptosis, flow cytometry reveals a nearly log-normal distribution with a coefficient of variation (CV; a unit-less measure of variability equal to the standard deviation divided by the mean) ranging from 0.2 to 0.3 (Spencer et al., 2009). Such a spread results in cells in the top 5th percentile having >2.5× higher protein expression compared to cells in the bottom 5th percentile (Niepel et al., 2009). The question then arises of whether such modest variation in protein levels is sufficient to explain the observed variation in the timing of cell death. Model-based simulation suggested that it is: when the distribution of cell death times was computed for TRAIL-induced apoptosis assuming log-normally distributed protein concentrations, a close match was observed between the variability in simulation and experiment (Spencer et al., 2009). In the absence of any simple experimental test, the match between simulation and measurement increases our confidence in the hypothesis that natural variation in the levels of apoptotic regulators is responsible for variability in the time and probability of cell death.
Is it possible to establish a direct link between the levels of any single protein and the probability and timing of apoptosis? In principle such a measurement could be made by fluorescently tagging proteins of interest at the endogenous locus and then relating their levels to time of death using live-cell microscopy. However, mathematical modeling suggests that achieving reasonable predictability over cell fate would require single-cell measurement of many protein levels (as well as some posttranslational modifications), a difficult task. Alternatively, simulation suggests that predictability can be achieved by measuring the rates of critical reactions, such as the processing of caspase-8 substrates. Because this rate depends on the levels of multiple upstream proteins, measuring it is much more informative than simply knowing protein levels (Spencer et al., 2009). This conclusion implies a fundamental limit to our ability to predict cell fate based on single-cell proteomics.
Key goals for a combined model- and experiment-driven analysis of apoptosis are to understand how multiple cooperating and competing signals are integrated to effectively execute a binary death-survival decision, to determine why some processes and proteins are important in one cell type and not in another, and to predict the responses of cells to death ligands and chemotherapy drugs. A review of the literature thus far suggests that these goals remain largely unfulfilled. Skeptics will argue that quantitative analysis can only add details to existing conceptual frameworks or that mathematical models are too theoretical and too dependent on assumptions to be useful (although drawing a pathway diagram may involve just as many assumptions). A more generous and realistic assessment would be that mechanistic modeling of apoptosis has had an impact in motivating the collection and analysis of quantitative single-cell data, critically evaluating potential regulatory mechanisms, and investigating the origins of cell-to-cell variability.
Addressing the long-term goals of quantitative, model-driven biology will require major conceptual and technical advances. Most computational tools currently in use have been adapted from other fields, but understanding a biological system is nothing like fixing a radio. Cells are not well-mixed systems as encountered in chemistry, nor are they easily understood in terms of fundamental physical laws or obviously subject to the design principles (such as modularity) encountered in engineered systems. They resemble all of these to some extent, but systems biology is currently immersed in the uncharted process of working out which concepts from chemistry, physics, and engineering are most useful in understanding cells and tissues.
It is already evident that different research groups will continue to build models differing in scope and level of detail and customized to the biological questions being addressed. Current approaches to model building typically involve de novo creation of complex sets of equations in each paper. A lack of transparency in the underlying assumptions makes it difficult for practitioners, nevermind the general research community, to understand how models differ from each other. Fortunately, “rules-based” modeling methods now in development promise to address the issue of model reusability and intelligibility (Faeder et al., 2009; Hlavacek et al., 2006). More rigorous means for linking models to experimental data and for understanding which aspects of a model are supported by data are required. Progress in this area is slow, but the basic principles are understood in the context of engineering and the physical sciences (Jaqaman and Danuser, 2006). Finally, we must work to ensure basic familiarity with dynamical systems among trainees. It is widely accepted that a working knowledge of statistical methods such as clustering and regression is essential in contemporary biomedicine, but it is unfortunate that few students are taught that familiar Michaelis-Menten equations are simply approximations to a mass-action formalism written as networks of ODEs (Chen et al., 2010).
Cancer pharmacology is the area of translational medicine in which models of apoptosis are most obviously of value. Critical questions in the development of rational and personalized treatment of cancer involve understanding precisely how anticancer drugs induce apoptosis, why the extent of cell killing varies so dramatically from one tumor to the next, and how we can predict response to chemotherapy, both “targeted” and cytotoxic. As yet no quantitative, model-based studies of these issues have been reported, but it seems almost certain that sensitivity and resistance will be controlled in a multifactorial manner. Genes and proteins that are important in one cellular setting will not be significant in another. In the case of TRAIL, for example, conventional molecular approaches have implicated the levels of O-glycosylation enzymes (GALNT3, GALNT14; (Wagner et al., 2007), TRAIL decoy receptors (DcR1, DcR2, and osteoprotegerin), c-FLIP, BclxL, and inhibitor of apoptosis proteins (IAPs) in TRAIL resistance in different cell lines (reviewed in Zhang and Fang, 2005). It is likely that all of these explanations are correct to some degree, and the key task therefore becomes understanding the role of context. This is precisely where models hold great promise, as they are able to quantify and weigh the contributions of multiple factors. Such context sensitivity could be implemented by using a model in which the topology and rate constants remain the same for all cell types but protein concentrations (initial conditions) are altered to match experimentally measured protein levels.
Ultimately, we need to understand the regulation of apoptosis in the context of real human tissues and tumors. Because mechanistic modeling is dependent on quantitative, multiplex data, this will not be straightforward, even in model organisms. New in vivo caspase activity probes (Edgington et al., 2009) and high-resolution intravital microscopy (Condeelis and Weissleder, 2010) will play an important role in data acquisition in vivo, but it also seems probable that the development of mechanistic models able to store, simulate, and rationalize results obtained across a panel of cancer cell lines will be essential. Such context-sensitive modeling might uncover a multifactorial measurement that could be made on real human tumors. Expression profiling and cancer genome sequencing also aspire to personalize cancer therapy, but the framework we envision is complementary in focusing on biochemical mechanism. A multiplex measurement method (BH3 profiling) already exists to estimate the propensity of cells to undergo apoptosis; it involves permeabilizing cells and then monitoring their responses to diverse BH3-only peptides (Deng et al., 2007). BH3 profiling can predict sensitivity to conventional chemotherapies and to the Bcl-2/BclxL antagonist ABT-737 (Deng et al., 2007). It would be valuable to construct a predictive mathematical framework for BH3 profiling and thereby generate precise mechanistic understanding of drug sensitivity and resistance that could be translated clinically.
Single-cell analysis of cellular responses to FasL and TRAIL has highlighted the dramatic impact of cell-to-cell variability in determining the timing and probability of response. That cells surviving exposure to a death ligand or cytotoxic drug can resume normal proliferation is a testament to the “stiff trigger,” “all-or-nothing” nature of the apoptotic switch. Cells that cross the threshold for MOMP are normally fully committed to die, whereas cells that remain below it can recover and continue to proliferate. In the case of receptor-mediated apoptosis, the presence of a dose-dependent variable delay preceding MOMP followed by a dose-independent and nearly invariant post-MOMP period likely reflects the evolutionary advantages of such a system. Variability in the timing and probability of apoptosis makes it possible for a uniform population of cells to respond to a prodeath stimulus in a graded manner, even though the response is binary at the single-cell level. In contrast, by undergoing MOMP and effector caspase activation in a rapid and invariant way, cells avoid the highly deleterious effects of initiating but not completing apoptosis; these effects include formation of “undead” cells with damaged genomes.
Variability in response appears to be universal across diverse cell lines and proapoptotic stimuli (Cohen et al., 2008; Gascoigne and Taylor, 2008; Geva-Zatorsky et al., 2006; Orth et al., 2008; Sharma et al., 2010; Spencer et al., 2009; Huang et al., 2010). For example, Gascoigne and Taylor (2008) characterized the response of 15 cell lines to three different classes of antimitotic drugs and found significant inter- and intra-cell line variation, with cells exhibiting multiple distinct phenotypes in response to the same treatment. Cohen et al. (2008) correlated variability in the levels of two proteins with the life-or-death response to the cancer drug camptothecin. Most recently, Sharma et al. (2010) detected a small subpopulation of reversibly “drug-tolerant” cells following treatment with cisplatin or the epidermal growth factor receptor inhibitor erlotinib. The significance of these findings is that cancer therapy is beset by the problem of fractional, or incomplete, killing of tumor cells. Multiple explanations have been proposed for fractional killing, including drug insensitivity during certain phases of the cell cycle, genetic heterogeneity, incomplete access of tumor to drug (Chabner and Longo, 2006; Skeel, 2003), and the existence of drug-resistant cancer stem cells (Reya et al., 2001). Single-cell imaging and computational modeling of apoptosis have added to this list cell-to-cell variability in protein levels arising from stochasticity in protein expression (Spencer et al., 2009). A critical task for the future will be to ascertain the relative importance of these processes in determining the extent of fractional killing with real tumors and therapeutic protocols. Because a wide variety of biochemical processes are involved, all operating on different timescales, developing an appropriate quantitative framework will be a key step to better understanding.
The authors thank J. Albeck, J. Bachman, D. Flusberg, S. Gaudet, T. Letai, and C. Lopez for their help and acknowledge NIH grants GM68762 and CA139980 for support.