Discussion of Models
Emergent Properties of Networks of Biological Signaling Pathways
Discussion of Models: Accuracy, Robustness, and Diffusion
Source: Science. 1999 Jan;283(5400):381-7.
We address three major issues regarding accuracy of these models.
Accuracy
Numerical accuracy
The exponential Euler formulation was used for integration. This is an explicit method based on the Forward Euler method, but it assumes exponential decays of variables. This makes it well suited for handling kinetic simulations. The accuracy improves with shorter simulation time-steps.
Implicit simulation techniques with variable time-steps were also examined, but these did not offer significant improvement in simulation times. The simpler Exponential Euler method was therefore used in all the model calculations.
The numerical accuracy of the computation was verified in three ways:
- Simple kinetic schemes which could be calculated analytically were modeled and simulated results compared to the analytical results.
- The law of mass conservation and microscopic reversibility principles were used as tests of accuracy in complex reaction schemes.
- The same model was run at different time steps and the resulting simulated values were compared. The calculations were shown to be stable and to provide convergent solutions for the range of time steps used in the study.
In each case the simulations met accuracy criteria of 5% or better. We did not seek to improve accuracy further, as the experimental uncertainty was considerably larger. Computations were carried out using the well stirred assumption, i.e. that all molecules had equal access to all others. This means that diffusion and compartmentalization were neglected (see discussion below).
Experimental accuracy and limitations of data
In carrying out simulations of natural processes, there are two major issues that must be considered with regards to the experiments from which the constants are obtained. The first is whether all the relevant experimental variables have been considered. The second is whether the experimental data is accurate enough.
In our models, missing experimental variables would correspond to signaling interactions which we have omitted. It is likel that some if not many interactions may be missing in our models. This is partly due to our inability to model every reported interaction, and partly because in all likelihood further interactions remain to be discovered. We would, however argue, that we have probably accounted for a large proportion of the important interactions. Current experimental strategy is well suited to defining major interactions first, and then homing in on details. The analysis of robustness of emergent behaviour also indicates that newly discovered interactions or regulation within the individual pathways may not have a major effect on the overall properties of the systems. Additional positive or negative connections between pathways may strengthen or weaken the system as a whole.These interactions as they are discovered should allow for further refinement of these models.
The second issue, experimental accuracy, is one which we have taken pains to estimate (Figure 2F in the paper). We feel that barring the conceptual limitations mentioned below, our simulations are likely to behave in qualitatively and quantitatively very similar ways within wide ranges of experimental uncertainty.
One of the important functions of a theoretical analysis of experimental results is to indicate where more rigorous experimental data would be required to fully understand a system. Biochemical reports are often usually good for modeling (Bower JM and Koch, C, Trends Neur. Sci. 15:11 1992.) in that a well-defined set of parameters such as Vmax and Km are reported with a complete estimates of error. However, a few recurring gaps limit accurate modelling, even though the raw data should be available from the same experiments. The most prominent is lack of information about time-courses, when reporting Kd for an equilibrium reaction. Direct time-course values such as profiles of reaction progress were available in most cases. In some cases the time-course had to be inferred from reports of the time of onset of a response, or by time-courses of down-stream effects. Additional data to estimate the ratio k2:k3 for Michaelis-Menten enzymes would also be very useful in constraining the models. In our simulations the ratio 4 was chosen as it results in a low percentage of enzyme sequestered in a complex with the substrate over most of the physiological range of enzyme and substrate concentrations. We varied this ratio over an order of magnitude in either direction and observed very little change in simulation results.
Conceptual limitations
Conceptual limitations to our models include the following issues which are topics of active research in many laboratories:
Cytoskeletal interactions and reaction channeling.
Protein scaffolds have recently been identified as important determinants of reaction specificity in several yeast MAP Kinase cascades. Potential scaffolding proteins for the MAP-kinase1,2 exist in mammalian systems as well (Schaeffer H. J et al Science 281:1668, 1998). Such anchors for active components can result in an assembly-line organization of signaling reactions that may be quite distinct from solution interactions both in terms of specificity and rates. We anticipate that such interactions may significantly change our understanding of signaling processes, and will need considerably more experimental data and possibly new modeling approaches.
Single-molecule interactions.
At physiological concentrations, there may be only a few individual molecules interacting in a dendritic spine, vesicle release site or other recation centers. In such situations the common approximation of molecular concentration as a continuous variable breaks down, and stochastic Markov processes are required to describe chemical interactions.
Robustness
The modelling approach outlined here is limited by the reliability of experimental information. Experimental error and variability between preparations introduces uncertainty in parameter estimates. Nevertheless, many aspects of our models and probably signalling pathways in general, are relatively insensitive to the exact parameter values (Figure 2F). The concentration of protein kinase C (PKC) was one of the most sensitive parameters, allowing 25% variation above or below the reference value while retaining the property of bistability. Most other parameters were robust within a factor of 2 in either direction. Work in preparation indicates that the bistable property of the system exists over a remarkably large region of parameter space. It is noteworthy that effective concentration of the signaling components appear to be critical parameters in determining bistable properties of the systems. This could well be a useful property of a biological system where the first stimulus stimulus triggers the synthesis of a crucial component to reach critical mass wherby subsequent stimuli are capable of tripping the switch. Such a system remains to be experimentally described. Current data do not definitely establish that bistability exists in systems which has a wide range of experimental values and large errors. Nevertheless, we feel that the available data suggest that the system we describe exists in a region where bistability is highly likely.
Diffusion and Compartmentalization
All reactions were modeled as point reactions, using the "well-stirred cell" assumption. Spatial aspects of calcium fluxes were examined in a limited manner, using explicit compartmentalization in the calcium model. We represented three calcium compartments: extracellular, intracellular and intracellular stores. We also carried out simple one-dimensional calcium diffusion calculations in the dendritic spine model. Several of our enzymes were also represented as transiting between cytosolic and membrane-bound forms. In our models this did not involve explicit spatial compartments. Instead we limited the reactions which the enzymes could participate in. This treatment clearly oversimplifies the enormous implications of spatial segregation and intracellular messenger diffusion in real cells. Additionally interactions between reactants that function in two and three dimensions may need to appropriately analyzed. Nevertheless, a point model such as ours is an essential first step towards building a more comprehensive model which incorporates diffusive and anchoring and compartmental effects.