Skip Navigation


Bioinformatics Advance Access originally published online on December 8, 2006
Bioinformatics 2007 23(3):358-364; doi:10.1093/bioinformatics/btl627
This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (Print PDF) Freely available
Right arrow Supplementary Data
Right arrow All Versions of this Article:
23/3/358    most recent
btl627v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrow Search for citing articles in:
ISI Web of Science (10)
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by Bagheri, N.
Right arrow Articles by Doyle, F. J.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Bagheri, N.
Right arrow Articles by Doyle, F. J., III
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?

© The Author 2006. Published by Oxford University Press. All rights reserved. For Permissions, please email: journals.permissions@oxfordjournals.org

Quantitative performance metrics for robustness in circadian rhythms

Neda Bagheri 1, Jörg Stelling 2,3 and Francis J. Doyle, III 1,4,*

1 Department of Electrical and Computer Engineering, University of California in Santa Barbara CA 93106-9560, USA
2 Institute of Computational Science 8092 Zurich, Switzerland
3 Swiss Institute of Bioinformatics, ETH Zurich 8092 Zurich, Switzerland
4 Department of Chemical Engineering, University of California in Santa Barbara CA 93106-5080, USA

*To whom correspondence should be addressed.


    ABSTRACT
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 

Motivation: Sensitivity analysis provides key measures that aid in unraveling the design principles responsible for the robust performance of biological networks. Such metrics allow researchers to investigate comprehensively model performance, to develop more realistic models, and to design informative experiments. However, sensitivity analysis of oscillatory systems focuses on period and amplitude characteristics, while biologically relevant effects on phase are neglected.

Results: Here, we introduce a novel set of phase-based sensitivity metrics for performance: period, phase, corrected phase and relative phase. Both state- and phase-based tools are applied to free-running Drosophila melanogaster and Mus musculus circadian models. Each metric produces unique sensitivity values used to rank parameters from least to most sensitive. Similarities among the resulting rank distributions strongly suggest a conservation of sensitivity with respect to parameter function and type. A consistent result, for instance, is that model performance of biological oscillators is more sensitive to global parameters than local (i.e. circadian specific) parameters. Discrepancies among these distributions highlight the individual metrics' definition of performance as specific parametric sensitivity values depend on the defined metric, or output.

Availability: An implementation of the algorithm in MATLAB (Mathworks, Inc.) is available from the authors.

Contact: frank.doyle{at}icb.ucsb.edu

Supplementary information: Supplementary Data are available at Bioinformatics online.


    1 INTRODUCTION
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 
Oscillatory processes are omnipresent in nature, comprising the cell cycle, neuron firing, ecological cycles and others; they govern many organisms' behaviors. A well-studied example of a biological oscillator is the circadian clock. The term circa- (about) diem (a day) describes a biological event that repeats approximately every 24 h. Circadian rhythms are observed at all cellular levels since oscillations in enzymes and hormones affect cell function, cell division and cell growth (Edery, 2000). They serve to impose internal alignments between different biochemical and physiological oscillations. Their ability to anticipate environmental changes enables organisms to organize their physiology and behavior such that they occur at biologically advantageous times during the day (Edery, 2000): visual and mental acuity fluctuate, for instance, affecting complex behaviors.

The concept of robustness relates to how a system maintains desired performance or functionality, despite internal or environmental perturbations (Stelling et al., 2004b). For oscillatory systems, specific functions include amplitude, period or phase (the relative timing with respect to a reference) (Fig. 1). The combination of robustness and sensitivity allows biological systems to adapt to their environments. For instance, jet lag occurs as a result of robust circadian properties that maintain the synchrony between biological clocks and their previous environment. Jet lag diminishes, however, due to the system's sensitivity toward light cues; a change in light pattern resets the circadian phase. The ability to maintain robust performance (i.e. circadian rhythms) in the face of perturbations and uncertainty (e.g. a displacement in longitudinal location), is a long-recognized and critical property of living systems.


Figure 1
View larger version (21K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Fig. 1 Characteristics of oscillators. A 15% parametric change in an asymptotically stable limit cycle oscillator forces the amplitude and period of oscillation to either increase or decrease relative to nominal (upper panel), and also affects the phase (lower panel). If the peak of the perturbed trajectory leads (lags) that of the nominal, there is a phase advance (delay).

 
Local robustness can be assessed through sensitivity analysis (Stelling et al., 2004b; Ma and Iglesias, 2002), by measuring the degree to which parametric perturbations dictate specific output dynamics. To date, performance metrics have only been formulated as a function of state, period or amplitude behavior. However, phase is a unique property of many biological oscillators because these oscillators are able to synchronize (or entrain) their phase to that of a forcing signal (e.g. light cues), while amplitude or period need not be affected (Pittendrigh and Daan, 1976a). Therefore, phase-based sensitivity analysis may provide a more discerning assessment of the system as it capitalizes on one of the more prominent features in biology, phase resetting (Winfree, 2001).

Although most methodologies of sensitivity analysis investigate a system's state (Varma et al., 1999), phase has been analyzed as a key performance attribute in few studies. Kramer et al. use isochrons (collections of points that evolve to the same position on the limit cycle) to measures phase advances/delays with respect to parametric perturbation (Kramer et al., 1984), reflecting information contained in existing phase response curves (Daan and Pittendrigh, 1976; Winfree, 2001). Similarly, Rand et al. make use of infinitesimal response curves to investigate phase dynamics as a function of independent finite length parametric perturbations (Rand et al., 2004). In each case, the respective output, or performance measure, is coupled: by quantifying phase as a single independent measure, the methods fail to isolate it from period dynamics. Parametric perturbations, however, often influence amplitude, period and phase simultaneously (Fig. 1), which requires new phase-based performance metrics.


    2 APPROACH
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 
Oscillatory systems exhibit a variety of characteristic output behavior, some of which include period, shape and phase. We investigated these output measures, or performance attributes, via parametric sensitivity analysis. A Drosophila melanogaster 10-state mathematical model (Leloup and Goldbeter, 1998) served to compare sensitivity distributions between performance metrics (Fig. 2). This moderately complex system consists of two coupled negative feedback loops that model the transcription, translation, phosphorylation and effective delays associated with period and timeless genes, and their protein counterparts (per and tim, and PER and TIM, respectively). While the model does not account for the complexity of the real network that, for instance, includes additional positive feedback loops (Smolen et al., 2001; Cyran et al., 2003), it has been experimentally validated (Leloup and Goldbeter, 1998) and is widely employed as a reference model (Stelling et al., 2004b; Gonze et al., 2002). Here, we focused on the development of phase-based methodologies to establish new metrics for investigating robust and sensitive properties that highlight specific network components used in maintaining system behavior. To confirm that the results can be generalized to other models, we included a 16-state Mus musculus model (Leloup and Goldbeter, 2003) in our analysis.


Figure 2
View larger version (16K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Fig. 2 The 10-state circadian model (adapted from Leloup and Goldbeter, 1998). Auto-inhibition of per and tim gene expression occurs via the nuclear PER/TIM complex. per and tim genes are transcribed in the nucleus, after which their mRNAs are transported into the cytosol where they undergo protein synthesis. The newly formed PER and TIM proteins are phosphorylated. The doubly phosphorylated proteins form a PER/TIM complex that enters the nucleus and closes the feedback loop.

 

    3 METHODS
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 
In order to generalize the application of sensitivity analysis, we treat the circadian system as a set of nonlinear ordinary differential equations (ODEs) such that x(t) defines the m-length state vector, {rho} defines the n-length time-invariant parameter vector and f (x(t), {rho}) defines the m-length system dynamics (see Supplementary material):


Formula 1

(1)
Sensitivity measures describe a change in output with respect to a change in input. Since biological circuits must function under a range of different parameters (Wagner, 2005), we regard the input as a parameter and the output as system performance.

We investigate performance by means of characteristic oscillatory output measures (Fig. 1) (Dunlap et al., 2004). The period of oscillation reflects the time required for a reference point in successive waves to pass a fixed point. The shape of an oscillatory system captures dynamic state behavior values typically observed by analyzing perturbation effects on the limit cycle. Shape sensitivity measures include amplitude sensitivity as it relates perturbed state maximum/minimum values to their nominal values. Phase dynamics are calculated via an angle (radian) or time (hour) framework. As an angular measure, certain state dynamics are projected onto a 2D space and investigated by means of the system's limit cycle (Fig. 3a). While the real-time state vector cycles around the asymptotically stable limit cycle, an angular relationship with respect to a fixed reference vector is established. These angular measures may define either an advance or delay of phase relative to the nominal trajectory (Fig. 1), with a delay defined as a shifting of the cycle to a later time.


Figure 3
View larger version (8K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Fig. 3 Time-dependent phase dynamics. (a) Illustrates how real-time phase dynamics are captured in a 2D limit cycle. Phase is the radian-based measure {theta}(t, {rho}) that describes the angular relationship between the state vector, x(t, {rho}), and some predetermined reference, r. Phase measures are recorded as a function of time under various parameter sets. (b) Illustrates the coupling of period and phase as trajectories diverge in time. (c) Depicts a decoupled measure of phase behavior as circadian phase dynamics are normalized with respect to their perturbation-induced periods.

 
3.1 State-based sensitivity metrics
Parametric state sensitivity is captured in an m by n matrix consisting of individual state performance values with respect to isolated parametric perturbations. In the study we employ the direct method (Varma et al., 1999; Khalil, 2002) as a means to determine exact parametric state sensitivity measures, Formula . This approach relies on the continuity of f(x(t), {rho}) with respect to the parameter vector, {rho}. Applying the chain rule results in partial derivatives of the function with respect to states and parameters producing an ordinary differential equation of sensitivity dynamics:


Formula 2

(2)
The initial conditions for the sensitivities are zero unless they rely on system parameters.

Coefficients of raw state sensitivity, Formula , rely on multiple coupled outputs (such as period and phase) and grow unbounded in time for parameters whose period sensitivities are non-zero (Larter, 1983; Tomovic and Vukobratovic, 1972; Zak et al., 2005). The secular term is due to computation of sensitivities involving a non-uniformly valid expansion of a periodic system (Larter, 1983). Tomovic and Vukobratovic (1972) demonstrate that state sensitivity, evaluated at the constant nominal period, {tau}, provides a specific measure corresponding to limit cycle behavior. This measure is unbiased to changes in the period/frequency of oscillation. The resulting m by n matrix is referred to as the cleaned-out or shape sensitivity measure, Formula . Sc(t) is periodic in time and describes how parametric perturbations affect the shape of state trajectories (Larter, 1983). In its raw form, state sensitivity for oscillatory systems may be decomposed into a combination of shape and period sensitivity measures (Tomovic and Vukobratovic, 1972):


Formula 3

(3)

This decomposition of raw state sensitivity highlights its linear time growth, Formula , while isolating period sensitivity, Formula . The decomposition is generally accomplished by calculating state and period sensitivities, and then solving for the cleaned-out shape sensitivity matrix.

We employ the method proposed by Zak et al. (2005) to determine period sensitivity, S{tau}, by making use of the decomposition (Equation (3)) and evaluating state sensitivities at a large time, t1 >> {tau}. At this time, the second term on the right-hand side of Equation (3) dominates the cleaned-out sensitivity matrix. Singular value decomposition of state sensitivity produces S(t1) = U{Sigma}VT, where {Sigma} is an m by n diagonal matrix of non-negative singular values, {sigma}, and matrices U and V contain the eigenvectors of SST and STS, respectively. Hence, period sensitivity of any given state may be approximated by:


Formula 4

(4)
where {sigma}1 is the largest singular value in {Sigma}, ||f (x(t1), {rho})|| is the vector norm of the state matrix evaluated at time t1, and v1 is the first column vector of V. This approximation holds true at any large times t >> {tau}, when the system is non-zero, and when the period of oscillation is sensitive to at least one parameter (Zak et al., 2005).

Amplitude sensitivity, S{Lambda}, describes how maximum state values vary due to independent parametric perturbations. As discussed, cleaned-out sensitivity defines how parametric perturbations affect the shape of state trajectories at every time in the cycle; at peak concentration times this measure relates directly to the state's maximum concentrations. Thus, we calculated amplitude sensitivity values from shape sensitivity by evaluating the time-dependent measure at peak concentrations times, tpeak:


Formula 5

(5)

3.2 Phase-based sensitivity metrics
In this work, we examine both standard (raw) phase, and decoupled (corrected) phase as performance measures, in addition to examining phase-based period and relative phase dynamics. Numerical methods are to approximate the exact solution of phase. We extract radian-based phase angles, {theta}(t, {rho}), from the system's limit cycle (Fig. 3a) using the cosine rule (Strang, 1988). Resulting phase dynamics reflect the oscillator's real-time position, or concentration, with respect to a static reference, r:


Formula 6

(6)
Phase measures are recorded under varying parameter sets, {theta} (t, Formula ) where Formula indicates measurements with respect to a perturbation of magnitude {delta} affecting a single parameter, {rho}j. This perturbation changes the nominal parameter vector from {rho} to a perturbed parameter vector, Formula :


Formula 7

(7)

The collected phase trajectories capture the oscillator's position and map them onto nominal time (Fig. 3b). If the perturbation strength and length are the same, the information contained in these trajectories is analogous to that contained in phase transition curves (Johnson, 1999).

Direct evaluation of phase trajectories (Fig. 3b) yields two types of sensitivity measures: period and phase because the change in period (or period sensitivity) is merely an accumulation in the change in phase (or phase sensitivity) evaluated over an entire cycle. In the case of period sensitivity, phase trajectories are evaluated at a 2{pi} · L-radian interval where the integer interval, L, is chosen arbitrarily. The difference between the perturbed 2{pi} · L -radian crossing time, Formula , and the nominal 2{pi} · L-radian crossing time, Formula , yields the system's periodic performance. A normalized time difference between perturbed and nominal phase trajectories defines the system's period sensitivity, Formula : a series of measurements denoting the quantitative change in period with respect to a change in the jth parameter,


Formula 8

(8)

In the case of strict phase sensitivity, we evaluate phase trajectories at specific times, t = tk. The radian phase difference between perturbed and nominal trajectories describes the raw phase sensitivity, S{theta}: a series of measures reflecting the induced change in phase with respect to a change in the jth parameter tk-hours after the perturbation,


Formula 9

(9)

Just as state sensitivity diverges over time, phase sensitivity grows unbounded due to the non-uniform expansion of a periodic system. To correct for the integrated response of perturbation effects (in this case, the coupling of period and phase), each phase trajectory is normalized with respect to the period of the system after parametric perturbation, Formula : dividing each time series by Formula decouples the system's phase from its period. As a result, normalized datasets begin and end at the same relative time points (0 and 100% of their respective cycles). These modified datasets allow for a comparison between nominal and perturbed phase at every point in the cycle (Fig. 3c). This corrected phase sensitivity assumes a linear scaling of raw phase measures resulting in a time-dependent performance quantity that identifies specific points of the cycle most susceptible to uncertainty:


Formula 10

(10)

Phase-based period, phase and corrected phase sensitivity analysis examine the biological network relative to a static reference. In some cases, a relative analysis that studies relationships within the perturbed network may be more useful. The timing effects relating transcription, translation, phosphorylation, and transport are governed by global cellular processes. Variation of these specific time intervals as a result of parameter manipulation indicates a degree of sensitivity. Relative sensitivity, S{varphi}, investigates the time interval relating the hour-difference between the occurrence of particular events; for instance, the time interval between peak mRNA concentrations and their corresponding protein concentrations. This time interval, {varphi}(x(t), {rho}), is a function of the system's state and parameter vectors. It explores how a system's individual components change relative to one another due to parametric perturbation:


Formula 11

(11)

3.3 Principles of comparison
To ensure an accurate comparison of these sensitivity distributions, several standards were applied. Time-varying metrics—state, shape, phase and corrected phase—were observed at 6, 12, 18 and 24 h of their respective nominal cycles (see Supplementary material in excel spreadsheet). Unless otherwise noted, we present the 18 h dataset. For the Drosophila model, phase-based analysis used the per mRNA and PER/TIM nuclear protein complex concentrations, respectively to capture circadian phase dynamics; for other combinations of states, and for the mammalian model, we refer to the Supplementary material. To further facilitate metric comparison, state, shape and amplitude performance measures refer only to per mRNA concentrations. As a result, their respective sensitivity distributions were biased toward per gene dynamics.


    4 RESULTS
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 
Parametric sensitivity analysis assigns absolute performance measures to each parameter. The greater the absolute sensitivity value, the more susceptible the system is to an isolated parametric perturbation. Plotting these sensitivities against one another enables the assessment of metric properties via correlation plots. Two metrics yield strongly similar performance distributions if their sensitivity measures (suitably normalized) align along a 45 degree line. State- and phase-based period sensitivity metrics produced such high-correlation (Fig. 4a), validating the accuracy of the proposed numerical phase methods.


Figure 4
View larger version (34K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Fig. 4 Metric evaluation. (a–c) Parametric sensitivity metrics ordered from least to greatest absolute value for the Drosophila model: (a) state- and phase-based period sensitivity, (b) corrected phase (decoupled angular phase trajectories) and relative phase (time interval between peak per mRNA and nuclear PER/TIM protein concentrations) sensitivity ranks and (c) amplitude-based and state-based metrics. Spearman rank correlation coefficients for these pairs of metrics are 1.00, 0.64 and 0.85, respectively. Legends describe the parameters' particular biological processes as the shading of each data point describes their type: global (open), mixed (gray) and local (black) parameters. Data points outlined in red reflect per gene dynamics, and those in green reflect tim gene dynamics. (d and e) Color-coded average sensitivities (values are scaled between 10–4 and 1) among parameter function (upper subplots) and parameter type (lower subplots) for each metric in the fly (d) and the mammalian (e) model. See the Supplementary material for alternate correlation diagrams and ranking plots. A colour version of this figure is available as supplementary data.

 
In an earlier publication, we introduced parameter ranking as a means of comparing results between different networks and/or analyses (Stelling et al., 2004b) where parameters are assigned integer numbers, or ranks, reflecting their absolute sensitivity distribution. Parameters are then plotted against one another from least to most sensitive; those nearest the origin are least sensitive. A correlation diagram of the two period sensitivity rank distributions shows minor deviations due to the nature of parameter ranking (refer to Supplementary material). Symmetry of the Drosophila model forces similar parameters (such as per and tim transcription, vs P,T) to have equal sensitivity values. Meanwhile, parameters with effectively equal sensitivities are assigned unique integer ranks, causing an artificial discrepancy between the two parameters. Note that metric assessment based on parameter ordering is not absolute, as no single rank, state or time point captures the complete system dynamics.

4.1 Parameter classification
Sensitivity measures reveal a system's susceptibility to changes in a particular set of parameters, as mild perturbations may significantly alter performance. Such findings highlight biochemical processes that impact performance. To further explore this phenomenon, we grouped model parameters according to their biochemical process—auto-inhibiting gene expression, transcription, translation, degradation, association/disassociation, phosphorylation/dephosphorylation and transport (for details on the classification of model parameters, see the Supplementary material of Stelling et al., 2004b). For Drosophila, the analysis identified phosphorylation and dephosphorylation rates as insensitive parameters when compared to rates of degradation, transport, translation and transcription (Figs 4a–c). Parameters associated with per gene and protein dynamics were highlighted in red, while those associated with tim gene and protein dynamics were highlighted in green. Each metric maintained a similar distribution: degradation and transport were consistently more sensitive from one performance measure to the next, while (de)phosphorylation rates were consistently less sensitive. There exists a conservation of sensitivity throughout the network regardless of the employed performance criterion.

Figures 4d and e (upper subplots) provide a more complete description of parameter sensitivity with respect to function for the fly and the mammalian model, respectively. The color of each cell reflects the average sensitivity rank associated with the parametric function for each performance measure. The color bar was indexed according to the range of output sensitivities within each plot, providing distinct shades or measures of sensitivity. The black horizontal lines serve to visually separate functional groups, while the vertical height of each group reflects the abundance of parameters in a group relative to the total number of model parameters. In Drosophila there were 8 parameters associated with phosphorylation, for example, 2 associated with transcription and only 1 associated with unspecified degradation. Phosphorylation and dephosphorylation encompassed much of the parameter space as they represented 16 of the 37 parameters.

4.2 Global versus local parameters
The investigated performance metrics depict classes of sensitivity that associate with parameter function, proving a certain conservation of robustness for specific biochemical processes. A similar sensitivity distribution was found when parameters were separated into three types: global, mixed and local (Stelling et al., 2004b). Global parameters are involved in core cellular reactions non-specific to the circadian rhythm; they encompass properties consistent with the entirety of the cellular network. Local parameters are primarily attributed to the circadian system; their processes and/or elements are not shared with many other cellular circuits. Global parameters include transcription rates, various mRNA and protein degradation rates, and translation rates, all of which were consistently found in the first (upper right) quadrant of Figures 4ac.

Local parameters, including protein phosphorylation and dephosphorylation, were consistently in the lower left (less sensitive) corner (see Stelling et al., 2004b) for details of parameter classification). The clustering of parameter sensitivities into global versus local groups was coherent in all sensitivity metrics. Performance was consistently more sensitive to perturbations in global and mixed parameters than local parameters. Figures 4d and e (lower subplots) emphasize this conservation of sensitivity apparent in both models by assigning a color reflective of the average sensitivity ranking within each parameter type.

Our results demonstrated that every defined performance metric was more sensitive to perturbations involving global and mixed parameters than it was to perturbations involving local parameters. Grouping parametric sensitivity based on parameter type provided a more consistent and distinct distribution of sensitivity measures among various metrics than the grouping of sensitivity by function. This outcome agreed with a previous study suggesting that circadian performance is greatly affected by changes in global parameters and less susceptible to changes in local parameters (Stelling et al., 2004a,b).

4.3 Specific performance correlations
Similarities among correlation diagrams related overall performance to parameter function and type. Differences between correlation diagrams related individual performance to specific biochemical processes (Figs 4ac). Note, that the significance of these relations cannot be stated in a rigorous statistical sense because for each combination of metric and parameter only one data point are available. Experimental observations, however, support the bias of certain output measures toward particular networks or biochemical processes. For instance, Allada et al. discussed the relevance of the Clk gene in Drosophila (not included in our model) in maintaining circadian functionality since its phase and amplitude perturbation (the timing and concentration levels of Clk oscillation) did not cause large changes in circadian periodicity (Allada et al., 2003).

Periodicity was consistently more sensitive to PER and TIM protein degradation rates (vd P,T), and less sensitive to their auto-inhibiting of gene expression (KI P,T) and transcription rates (vs P,T) when compared to other metric rank distributions (see Supplementary material). Uncorrected phase rejected changes in protein degradation rates (vd P,T) as compared to other metrics and was more receptive to disturbances in translation (ks P,T) and PER/TIM complex degradation rates (kd C,N). States were more responsive to changes in protein degradation (Kd P,T) and mRNA degradation rates (vm P,T), and less responsive to translation rates (ks P,T) relative to other performance measures.

Figure 4b establishes corrected phase performance as being more susceptible to parametric changes involving transcription and auto-inhibition rates, and less susceptible to parametric changes involving (dis)association and mixed degradation rates. The corrected phase-based method ranked these parameters with lower order (or sensitivity) when compared to the relative phase method. These conclusions, however, were not apparent when comparing corrected phase results to other metrics, or comparing them at different times. Therefore, time-dependent metrics, such as shape and corrected phase sensitivity may need to be evaluated over an entire cycle to justify such generalizations.

The relative phase between peak concentrations of per mRNA and nuclear PER/TIM protein complex appeared to be most sensitive to mixed protein degradation (Kd P,T) and protein complex degradation rates (kd C,N); and least sensitive to TIM transcription (vsT), tim mRNA degradation (vmT), and auto-inhibition (KI P,T) rates. By definition, the metric was biased toward per gene dynamics. Changes in the tim gene did not readily affect the phase function. A change in PER or TIM protein, however, directly affected PER/TIM association, thereby amplifying relative phase susceptibility to protein degradation rates.

Studies by Meyer et al. (2006) support the need for relative phase analysis: evidence supporting the disassociation of PER/TIM in the cytoplasm and their independent travel into the nucleus highlights a key interval of time. In any given cell, nuclear accumulation rates were different and independent for each protein. Therefore, PER and TIM proteins appeared to act as constituents of an intracellular interval timer. Similarly, transcriptional regulation of wc-1 in Neurospora was shown to be responsible for the phase of the clock (Kaldi et al., 2006). Relative phase sensitivity analysis is directly applicable to the analysis of such phenomena.

4.4 Similarity of performance metrics
Interpretation of greater performance based on metric similarity may provide a top-down assessment of cellular topology as similar performance criteria may be structurally related. We expect certain performance criteria to be strongly correlated, and others, (such as amplitude and relative phase) to be more distant. A dendrogram using the Spearman rank correlation coefficient as a distance measure condenses the information captured in the correlation diagrams (Fig. 5). For both models, we find a clear separation of raw or corrected phase and period sensitivity, respectively (see also Figs 4de), which is robust even when only subsets of parameters are considered (see Supplementary material for details). This conservation of the distance across models underlines the significance of analyzing phase sensitivity for full characterization of the performance, which is enabled by the methods presented in this paper. Consequently, experimental design should consider consistently distant measures, such as period and phase sensitivity, for the comprehensive characterization of biological oscillators because closely related measures do not provide new information.


Figure 5
View larger version (18K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Fig. 5 Metric similarity analysis. Similarity is defined by the Spearman rank correlation coefficient for all pairs of rank distributions, where the length of each link represents the normalized distance between two metrics (solid lines = Drosophila, dashed lines = Mus musculus). For both models, correlation coefficients rij varied between 0.58 and 1.00; the corresponding distances were defined as dij = 1 – rij and normalized by the maximal dij. Complete linkage was employed to construct the dendrogram.

 
For some states, relative phase and period dynamics are very similar. Relative phase examines the time interval between different cellular components while period measures the time interval between events of the same component in successive intervals. In theory, a parametric perturbation may not change the periodicity of the system while it may advance an event in one state and delay an event in another. Our results, however, suggest that the relative timing between certain states is critical for periodicity (see below and Supplementary material). Other correlations, such as the similarity and dissimilarity between corrected phase and phase (fly versus mammalian, respectively), may be specific to a biological model or due to bias/errors. Based on our analysis of other models for Drosophila (Stelling et al., 2004b) and for mammalian circadian clocks (unpublished data), we favor the former possibility. In principle, such specificity is not surprising since each of the metrics is related to a particular definition of robustness. Structural features point to the possibility of unraveling nature's design principles as each biological output is not tied to a specific metric, but to the performance of the organism.

Many biological circuits are multi-functional, but their analysis often focuses on few functions. As we prove in this study, the robustness of one function (or performance metric) does not imply a robustness of other functions, though it may allow the network to acquire greater functionality (Wagner, 2005). State-based metrics examine how parameters affect concentration values while phase-based metrics investigate the frequency and speed of concentrations relative to one another as they travel about their asymptotically stable limit cycle. Consequently, there are few instances in which any two metrics are strongly correlated. In fact, weak correlation among metrics is expected since each metric provides a unique measure of performance that relies on the defined output.

These different analyses of clock functions reveal details on the nature of the intracellular timer. For instance, in Drosophila, relative phase of per mRNA and peak PER/TIM nuclear protein concentrations is more sensitive to changes in protein degradation than to changes in tim transcription. Interestingly, relative phase correlates with period sensitivity for those pairs of states involved in establishing time delays (i.e. cascaded (de)phosphorylation) in the negative feedback loops required for oscillations (Supplementary Fig. 4). This conclusion derives from a similarity analysis of three definitions of relative phase: the time interval between peak PER protein and phosphorylated PER protein concentrations, peak per mRNA and nuclear PER/TIM protein concentrations and peak doubly phosphorylated PER protein and PER/TIM protein concentrations. The investigation of relative phase between PER and TIM nuclear accumulation may even better address the existence and performance of such an intracellular timer (Meyer et al., 2006). Additionally, it may better focus on entrainment phenomena associated with non-parametric models through use of skeletal photoperiodism, by examining the time interval between the onset of light and circadian activity (Pittendrigh and Daan, 1976b).

Performance metrics may also be used to highlight bias toward particular network functions. Two distinct perceptions of network performance arise from enhancing the Drosophila model through addition of genes, such as Clock. Allada et al. investigate performance with respect to individual feedback loops: they show that the PER/TIM protein complex loop regulates period and phase dynamics, while the CLK/CYC protein complex loop (that is not included in our model) regulates amplitude, but not periodicity.


    5 CONCLUSION
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 
The increasing complexity of biological models makes it difficult to untangle the roles of different mechanisms in determining rhythmic period, phase and amplitude (Allada et al., 2003). With molecular detail emerging on connected oscillators, such as the cell cycle and the circadian clock (Pregueiro et al., 2006), a quantitative and systematic framework for the analysis of increasingly complex mathematical models becomes mandatory. Here, we present such a framework by introducing four unique metrics that use phase dynamics as their measure of performance: period, phase corrected phase and relative phase sensitivity. These metrics and their calculation are general and not confined to the models investigated in this study; however, they rely on the existence of regular limit cycles and monotonic phase dynamics (Fig. 3b).

State- and phase-based analysis yield metric-specific sensitivity distributions, suggesting that robust performance is a function of the output measure for a given metric. However, there exists a conservation of sensitivity among parameter function and type; this consistency supports general theories relating performance more to system structure than to fine-tuning of parameters (Stelling et al., 2004b). Such findings motivate further investigation of phase-based sensitivity to better understand the delegation of performance to specific network components. As our methods are generic and provide parametric sensitivity values for a variety of performance measures, the identification of similar architectures in different biological contexts could provide a basis to assume similar functional properties without the need to characterize the system experimentally (Guantes and Poyatos, 2006). Hence, through use of performance metrics, we are able to identify critical control mechanisms and to better address nature's design principles.


    Acknowledgments
 
The authors thank the anonymous referees for their suggestions that helped improve the paper. This work was supported by the Institute for Collaborative Biotechnologies through grant DAAD19-03-D-0004 from the ARO; IGERT NSF grant DGE02-21715; DARPA BioSpice Program; and the Research Participation Program between the US DOE and AFRL/HEP.

Conflict of Interest: none declared.


    FOOTNOTES
 
Associate Editor: Martin Bishop

Received on August 15, 2006; revised on December 5, 2006; accepted on December 5, 2006

    REFERENCES
 TOP
 ABSTRACT
 1 INTRODUCTION
 2 APPROACH
 3 METHODS
 4 RESULTS
 5 CONCLUSION
 REFERENCES
 

    Allada, R., et al. (2003) A recessive mutant of Drosophila Clock reveals a role in circadian rhythm amplitude. EMBO J, . 22, 3367–3375[CrossRef][Web of Science][Medline].

    Cyran, S., et al. (2003) vrille, Pdp1, and dClock form a second feedback loop in the Drosophila circadian clock. Cell, 112, 329–341[CrossRef][Web of Science][Medline].

    Daan, S. and Pittendrigh, C.S. (1976) A functional analysis of circadian pacemakers in nocturnal rodents. II. The variability of phase response curves. J. Comp. Physiol, . 106, 253–266[CrossRef].

    In Dunlap, J.C., Loros, J.J., DeCoursey, P.J. (Eds.). Chronobiology: Biological Timekeeping, (2004) , Sunderland, MA Sinauer Associates, Inc.

    Edery, I. (2000) Circadian rhythms in a nutshell. Physiol. Genomics, 3, 59–74[Abstract/Free Full Text].

    Gonze, D., et al. (2002) Robustness of circadian rhythms with respect to molecular noise. Proc. Natl Acad. Sci. USA, 99, 673–678[Abstract/Free Full Text].

    Guantes, R. and Poyatos, J.F. (2006) Dynamical principles of two-component genetic oscillators. PLoS Comput. Biol, . 2, 0188–0197.

    Johnson, C.H. (1999) Forty years of PRCs–what have we learned? Chronobiol. Int, . 16, 711–743[Web of Science][Medline].

    Kaldi, K., et al. (2006) Transcriptional regulation of the Neurospora circadian clock gene wc-1 affects the phase of circadian output. EMBO Rep, . 7, 199–204[CrossRef][Web of Science][Medline].

    Khalil, H.K. Nonlinear Systems, (2002) , Upper Saddle River, NJ Prentice hall.

    Kramer, M.A., et al. (1984) Sensitivity analysis of oscillatory systems. Appl. Math. Modelling, 8, 328–340.

    Larter, R. (1983) Sensitivity analysis of autonomous oscillators: separation of secular terms and determination of structural stability. J. Phys. Chem, . 87, 3114–3121[CrossRef].

    Leloup, J.-C. and Goldbeter, A. (1998) A model for circadian rhythms in Drosophila incorporating the formation of a complex between the PER and TIM proteins. J. Biol. Rhy, . 13, 70–87.

    Leloup, J.-C. and Goldbeter, A. (2003) Toward a detailed computational model for the mammalian circadian clock. Proc. Natl Acad. Sci. USA, 100, 7051–7056[Abstract/Free Full Text].

    Ma, L. and Iglesias, P.A. (2002) Quantifying robustness of biochemical network models. BMC Bioinformatics, 3, 38[Medline].

    Meyer, P., et al. (2006) PER-TIM interaction in living Drosophila cells: an interval timer for the circadian clock. Science, 311, 226–229[Abstract/Free Full Text].

    Pittendrigh, C.S. and Daan, S. (1976a) A functional analysis of circadian pacemakers in nocturnal rodents. I. The stability and lability of spontaneous frequency. J. Comp. Physiol, . 106, 223–252[CrossRef].

    Pittendrigh, C.S. and Daan, S. (1976b) A functional analysis of circadian pacemakers in nocturnal rodents. IV. Entrainment: Pacemaker as clock. J. Comp. Physiol, . 106, 291–331[CrossRef].

    Pregueiro, A.M., et al. (2006) The Neurospora checkpoint kinase 2: a regulatory link between the circadian and cell cycles. Science, 313, 644–649[Abstract/Free Full Text].

    Rand, D.A., et al. (2004) Design principles underlying circadian clocks. J. R. Soc. Interface, 1, 119–130[Medline].

    Smolen, P., et al. (2001) Modeling circadian oscillations with interlocking positive and negative feedback loops. J. Neurosci, . 21, 6644–6656[Abstract/Free Full Text].

    Stelling, J., et al. (2004a) Robustness of cellular functions. Cell, 118, 675–685[CrossRef][Web of Science][Medline].

    Stelling, J., et al. (2004b) Robustness properties of circadian clock architectures. Proc. Natl Acad. Sci. USA, 101, 13210–13215[Abstract/Free Full Text].

    Strang, G. Linear Algebra and its Applications, (1988) , New York, NY Saunders College Publishing.

    Tomovic, R. and Vukobratovic, M. General Sensitivity Theory, (1972) , New York, NY Elsevier.

    Varma, A., Morbidelli, M., Wu, H. Parametric Sensitivity in Chemical Systems, (1999) , New York, NY Oxford University Press.

    Wagner, A. (2005) Circuit topology and the evolution of robustness in two-gene circadian oscillators. Proc. Natl Acad. Sci. USA, 102, 11775–11780[Abstract/Free Full Text].

    Winfree, A. The Geometry of Biological Time, (2001) , New York, NY Springer.

    Zak, D.E., et al. (2005) Sensitivity analysis of oscillatory (bio)chemical systems. Comput. Chem. Eng, . 29, 663–673[CrossRef].


Add to CiteULike CiteULike   Add to Connotea Connotea   Add to Del.icio.us Del.icio.us    What's this?


This article has been cited by other articles:


Home page
J Biol RhythmsHome page
N. Bagheri, M. J. Lawson, J. Stelling, and F. J. Doyle
Modeling the Drosophila melanogaster Circadian Oscillator via Phase Optimization
J Biol Rhythms, December 1, 2008; 23(6): 525 - 537.
[Abstract] [PDF]


Home page
J Biol RhythmsHome page
R. Gunawan and F. J. Doyle III
Phase Sensitivity Analysis of Circadian Rhythm Entrainment
J Biol Rhythms, April 1, 2007; 22(2): 180 - 194.
[Abstract] [PDF]


Home page
Cold Spring Harb Symp Quant BiolHome page
J. Rougemont and F. Naef
Stochastic Phase Oscillators and Circadian Bioluminescence Recordings
Cold Spring Harb Symp Quant Biol, January 1, 2007; 72(0): 405 - 411.
[Abstract] [PDF]


This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (Print PDF) Freely available
Right arrow Supplementary Data
Right arrow All Versions of this Article:
23/3/358    most recent
btl627v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrow Search for citing articles in:
ISI Web of Science (10)
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by Bagheri, N.
Right arrow Articles by Doyle, F. J.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Bagheri, N.
Right arrow Articles by Doyle, F. J., III
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?