Abstract
Background
Sustained stimulation with tumour necrosis factor alpha (TNFalpha) induces substantial oscillations—observed at both the single cell and population levels—in the nuclear factor kappa B (NFkappa B) system. Although the mechanism has not yet been elucidated fully, a core system has been identified consisting of a negative feedback loop involving NFkappa B (RelA:p50 heterodimer) and its inhibitor Ikappa Balpha. Many authors have suggested that this core oscillator should couple to other oscillatory pathways.
Results
First we analyse singlecell data from experiments in which the NFkappa B system is forced by short trains of strong pulses of TNFalpha. Power spectra of the ratio of nucleartocytoplasmic concentration of NFkappa B suggest that the cells' responses are entrained by the pulsing frequency. Using a recent model of the NFkappa B system due to Caroline Horton, we carried out extensive numerical simulations to analyze the response frequencies induced by trains of pulses of TNFalpha stimulation having a wide range of frequencies and amplitudes. These studies suggest that for sufficiently weak stimulation, various nonlinear resonances should be observable. To explore further the possibility of probing alternative feedback mechanisms, we also coupled the model to sinusoidal signals with a wide range of strengths and frequencies. Our results show that, at least in simulation, frequencies other than those of the forcing and the main NFkappa B oscillator can be excited via sub and superharmonic resonance, producing quasiperiodic and even chaotic dynamics.
Conclusions
Our numerical results suggest that the entrainment phenomena observed in pulsestimulated experiments is a consequence of the high intensity of the stimulation. Computational studies based on current models suggest that resonant interactions between periodic pulsatile forcing and the system's natural frequencies may become evident for sufficiently weak stimulation. Further simulations suggest that the nonlinearities of the NFkappa B feedback oscillator mean that even sinusoidally modulated forcing can induce a rich variety of nonlinear interactions.
Background
Nuclear factor kappa B (NFκB) transcription factors are critical to the control of response to cellular stress and are also involved in the regulation of cellcycle/growth, survival, apoptosis, inflammation and immunity [15]. They are dimeric molecules, composed of either homoor heterodimers, with the most common form being the RelA:p50 heterodimer.
In resting cells, NFκB dimers are sequestered in the cytoplasm by members of a family of molecules, the inhibitors of κB or IκBs. When the cells are stimulated with tumour necrosis factor alpha (TNFα) certain kinases, the IKKs, are activated and phosphorylate IκBs, marking them for degradation. The degraded IκBs then release NFκB dimers which are free to translocate into the nucleus, where they bind to specific sequences in the promoter or enhancer regions of target genes, including those for IκBα and the zinc finger protein A20 [6,7]. Newly synthesized IκBα migrates to the nucleus, binds to NFκB dimers and removes them from the nucleus, while A20 protein stays in cytoplasm and represses the activity of TNFα receptors [8]. Hence the NFκB system includes at least two negative feedback loops, one involving cytoplasmic sequestration mediated by IκBα and another involving A20.
Coupled negative feedback loops are known to have the potential to support oscillations [9,10]. Biological evidence suggesting NFκB oscillations was reported by Hoffmann et al. [11], who did populationlevel studies of IκBα/embryonic fibroblasts using electromobility shift assays, and also observed by Nelson et al. [12] in the cellular concentration and nuclear: cytoplasmic localisation of IκBα and RelAfluorescent fusion proteins in single SKNAS cells. Single cell timelapse imaging data showed persistent cycling of NFκB localisation between the cytoplasm and nucleus of SKNAS cells in response to continual TNFα stimulation [12,13].
The expression of genes regulated by NFκB is tightly coordinated with the activities of many other signalling and transcriptionfactor pathways [1418] including the p53 signalling pathway. Though the canonical NFκB signalling pathway has been studied extensively, the existence and mechanisms of the interactions between the NFκB pathway and other signalling pathways are unclear.
Ashall et al. [19] reported that different patterns of pulsatile stimulation of the NFκB system lead to different patterns of NFκB dependent gene expression, supporting the view that the frequency of the oscillations may have a functional role. It is therefore important to understand how the system's response frequencies are influenced by interactions with other oscillatory pathways. If both the NFκB network and its couplings to other other oscillatory pathways were purely linear, then it would be straightforward to use the machinery of transfer functions to characterise their interactions. In particular, one would expect the power spectrum of a periodically stimulated system to have its power concentrated at the forcing signal's frequency and its harmonics. But the NFκB network is a highly nonlinear system of coupled chemical reactions and hence, as we will demonstrate in modelling studies below, its response to periodic stimulation can include complex interactions between the intrinsic negative feedback oscillator and the stimulus. In what follows we work with a deterministic model first described in [19] and examine the power spectral densities of the time courses of NFκB localisation when the model is subjected to two types of periodic signals: trains of rectangular pulses and sinusoidal signals. Both sorts of stimulus are intended as a proxies for the influence of other oscillatory pathways on the core NFκB feedback loop and we find that such interactions can produce a rich variety of nonlinear dynamical behaviour. Our results on sinusoidal forcing are in broad agreement with those of Fonslet et al. [20], who applied a related family of sinusoidal stimulation protocols to a reduced model of the core NFκB oscillator originally developed in [21].
Results and Discussion
Experimental phenomena and time series analysis
We begin by analyzing the time courses of the localisation of fluorescently tagged NFκB in single cells subjected to TNFα stimulation (see [22] for details). The panels in the top row of Figure 1 show the ratio of nucleartocytoplasmic fluorescence for three patterns of stimulation. In all cases the cells received three strong (10 ng/ml) pulses of TNFα, each of five minute duration. These pulses were separated by various interpulse intervals—55 minutes for panel (a),95 minutes for panel (b) and 195 minutes for panel (c)—resulting in net stimulus periods of 60, 100 and 200 minutes. The first of these is rather shorter than the observed period of the oscillations induced by constant TNFα stimulation, which is close to 100 minutes.
Figure 1. Power spectra of singlecell data. Panels (a), (b) and (c) show time courses of the ratio of nucleartocytoplasmic fluorescence recorded from single cells while panels (d), (e) and (f) show the corresponding power spectral densities, normalized to have unit area beneath the curve. Note that the frequencies in panels (d)(f) are measured in units of the stimulus frequency and so the horizontal scale varies from panel to panel.
The bottom row of Figure 1 shows the corresponding power spectral densities and makes the point that the intrinsic oscillation—the one induced by constant stimulation—does not appear to have been excited by these pulsed stimulation protocols: in all cases the peaks in the power spectral density appear at multiples of the stimulus frequency. In mathematical terms, the cells respond as though they have been perturbed away from a stable resting state and do not offer evidence of any more complex internal dynamics. Nonetheless, these data provided crucial input to the development of the models in [19]. In particular, note that the amplitude of the response to the second and third pulses of stimulation varies as a function of the interpulse period, being somewhat reduced in panels (a) and (b), but not in (c): this constrains the rate at which the model should relax toward equilibrium.
One possibility is that the pulsing signal is too strong and suppresses the response at the natural frequency. In the following section we report in numerical experiments, based on the deterministic model in [19], through which we investigate this possibility by simulating the system's response to a periodic train of square pulses with strength varying from 0 to 10 ng/ml and interpulse intervals ranging from 20 to 295 minutes.
Model Introduction and Bifurcation Analysis
In this section, we introduce a deterministic mathematical model developed in [19]. It includes two coupled negative feedback loops—one for NFκBIκBα interactions and another that models the dynamics of A20—as is illustrated in Figure 2. Here A20 acts by modifying the activity of IκB kinase (IKK), which transduces the TNFα signal by phosphorylating IκBα, thus initiating the production of free NFκB. IKK is assumed to exist in one of three states: neutral IKK (denoted IKKn, a state ready for activation), active IKK (the state capable of phosphorylating IκBα) and inactive IKK (IKKi, a state unable to phosphorylate IκBα and also incapable of being activated by the TNFα signal). TNFα activates IKK by transforming neutral IKK into active IKK, which then acts to liberate NFκB. Active IKK is assumed to convert to its inactive form spontaneously, with linear kinetics. In the absence of A20, this inactive IKK would then transform back into neutral IKK, also with linear kinetics. A20 downregulates NFκB activity indirectly, by retarding the transformation of inactive IKK to the neutral form: Equations (1) provide a simple model for these dynamics:
Figure 2. The deterministic, twoloop model of [19]. Proteins are shown as coloured ovals, while mRNA's are shown roundcornered rectangles. In addition to the direct negative feedback through which nuclear NFκB leads to the production of IκBα, and thus to its own inactivation, the model also includes an indirect negative feedback mediated by A20's influence on the signal transduction machinery.
Here we have suppressed the time dependence of the concentrations of the various forms of IKK and of the A20 protein. The k_{* }are fixed parameters and T_{R }is a parameter indicating the presence or absence TNFα stimulation: T_{R }= 1 when the system is being stimulated with 10 ng/ml TNFα and T_{R }= 0 otherwise. This is a simplification: TNFα does not act directly on IKK, but rather produces its effect through a cascade of chemical reactions that begins when TNFα binds its receptor at the cell surface.
Following Ashall et al., we treat only the final stages of this chain. Thus when we model strong stimulation by setting T_{R }= 1, we abstract away a great deal of transduction machinery. In the studies that follow we will want to explore the consequences of weaker stimulation (TNFα concentrations of 10 ng/ml are orders of magnitude higher than physiological levels: Matalka et al. [23] report measurements from various tissues in healthy mice and found concentrations in the range 15 pg/ml while Prabha et al. [24] studied the plasma of healthy human subjects and found TNFα concentrations on the order of 100 pg/ml.) and it is natural to generalize the role of T_{R}, allowing it to vary across the interval 0 ≤ T_{R }≤ 1. In light of the modelling assumptions discussed above, one should not imagine that T_{R }depends linearly on the concentration of TNFα, but only that T_{R }increases monotonically with dose.
The model of Ashall et al. captures two important dynamical features observed in single cell data [12,19]: (i) continuous stimulation with high doses of TNFα leads to sustained oscillations with a period of around 100 minutes and (ii) pulsatile stimulation with the same high concentration of TNFα leads—as illustrated in Figure 1—to entrainment of the response by the pulsing signal. Given that continuous, strong stimulation induces sustained oscillations, one is prompted to ask how the existence of these oscillations depends on the strength of stimulation. We addressed this question by doing a bifurcation analysis using T_{R }as the parameter: Figure 3 is the resulting diagram. The model has a Hopf bifurcation (HB) at T_{R }= T_{* }≈ 0.366, which means that for T_{R }> T_{* }the ratio of nucleartocytoplasmic NFκB concentration exhibits sustained oscillations, but for T_{R }< T_{* }only damped oscillations. A forthcoming paper, Wang et al. [25], will present a comprehensive survey of the bifurcation structure of this model, but here we will concentrate on spectral analysis of the responses.
Figure 3. Qualitative behaviour as a function of T_{R}. At left, a bifurcation diagram showing how the oscillatory response depends on T_{R}. The solid black curve in the region 0 ≤ T_{R }≤ T_{* }represents a branch of stable equilibriasteady, nonoscillating solutionswhile the dashed black curve for T_{R }>T_{* }indicates a branch of unstable steady states. The red curves show the limitsthe peak and trough valuesof the stable oscillatory responses that exists for these values of T_{R}. The panel at right shows period of the oscillation as a function of T_{R}.
Numerical experiments
All the numerical experiments reported here begin from the same initial condition, whose preparation is described in the Methods section. We simulated the consequences of two TNFα stimulation protocols: pulselike stimulation similar to that used in Nelson's experiments and sinusoidally modulated stimulation.
Pulsed Stimulation
From the experiments and simulations in [19], we know that stimulation with a finite train of three strong pulses tends to entrain the cell's response. In this section, we first show that when subjected to a long periodic train of strong pulses, the model's response frequencies are also entrained, but that when the model is driven by sufficiently weak periodic pulse trains, various resonance phenomena connected to the underlying natural oscillations become observable.
Figure 4 illustrates the continuouslypulsed analogue of the experiments from Figure 1. It shows the spectral content of the model's response to a long periodic pulse train in which five minute periods of strong (T_{R }= 1) stimulation alternate with periods during which T_{R }= 0. The panel at left shows a typical response to this sort of strong, pulsed, periodic forcing: after a brief transient the response is periodic with the same period as the forcing pulse train. The right panel, which is a heat map showing the power spectral density of the response as a function of pulsing frequency, shows that similar entrainment occurs over a wide range of frequencies: the bright lines—which correspond to peaks in the power spectral density—have integer slopes, indicating that the power in the response is concentrated at harmonics of the forcing frequency.
Figure 4. Periodic stimulation with strong pulses. The panel at left shows the model's response to sustained periodic stimulation in which strong (T_{R }= 1) pulses of 5 minute duration alternate with stimulusfree intervals of 55 minute duration. The panel at right shows the power spectral density (PSD) of the response as a function of pulsing frequency: a vertical strip cut from this panel provides a colourcoded version of the PSD of the models's response to pulsed forcing at the corresponding frequency.
However, when the stimulating pulses are weaker the nonlinearity of the system leads to more complex patterns of resonance. In particular, both subharmonic resonance—when periodic forcing excites a response at a rationallyrelated lower frequency—and superharmonic resonance—in which a pure sinusoid excites responses containing higher harmonics—are possible. Note that these terms are defined with reference to the forcing frequency, a convention used in, for example, [26]. Superharmonic resonance is harder to identify when, as with the rectangular pulses used here, the periodic forcing already has power at higher harmonics, but subharmonic resonance occurs when 0.01 < T_{R }< 0.2, as is evident in the power spectral densities summarized in Figure 5. In addition to the bright lines with integer slopes, lines with slopes and also appear. The first of these provides evidence that subharmonic resonance occurs over a wide range of pulsing frequencies.
Figure 5. Response to periodic trains of weak pulses. A heat map summarising the power spectral density of NFκB localisation when the model is subjected to pulsed forcing with T_{R }= 0.1. For forcing frequencies in the range 1.53.0 × 10^{4}. Hz. the lowestlying line segment, which has slope ≈ 1/2, provides evidence of subharmonic resonance.
When T_{R }≥ 0.2 the power spectral density of the response shows power only at integer multiples of the pulsing frequency, indicating that the response is fully entrained by the forcing. Given that the pulsestrength used in the experiments corresponds to T_{R }= 1, these modelling results are in qualitative agreement with the experimental data in Figure 1: the stimulation was so strong that we should have expected it to have entrained the responses completely.
Sinusoidally modulated stimulation
In this section we study the response of Horton's model to sinusoidally modulated stimulation of the form
Here 0 ≤ ε ≤ 1/(1 + η) is the time average of the stimulus strength while 0 ≤ η ≤ 1 and ν > 0 are, respectively, the relative amplitude and the frequency of the sinusoidal modulation. With this parameterization T_{R }(t) has period τ = 1/ν and range
The roles of the parameters are illustrated in Figure 6.
Figure 6. Sinusoidal forcing function. The function T_{R }(t) for the parameter values ν = 1, η = 0.5 and ε = 0.6.
The case with η = 1 is the straightforward substitution of pulse trains with sinusoidal waves having peaktotrough amplitude 2ε. Our motivations for this formulation are twofold: firstly, provided η < 1 we have T_{R}(t) > 0 at all times and so expect the system more readily to exhibit oscillations similar to those induced by constant stimulation. Secondly, periodic forcing of the form (2) converts our model into a sinusoidally forced nonlinear oscillator, a class of systems that has been studied very extensively (see, for example, Pikovsky et al. [27], Wiggins [26] or Nayfeh and Mook [28]).
When η = 0 the forcing (2) reduces to a constant stimulus with T_{R }(t) = ε and the bifurcation analysis illustrated in Figure 3 leads us to expect stable oscillations when ε > 0.366 and a stable steadystate otherwise. But for η > 0 the system's response depends delicately on the relationship between the forcing frequency ν and the natural frequency ν_{0}. When the forcing frequency ν and natural frequency ν_{0 }are close (the requisite degree of closeness depends on η) the response will become mode locked, or synchronized with the forcing: it will then have the same frequency ν as the forcing and its power spectral density will be concentrated around the harmonics of ν. When the difference (ν ν_{0}) is larger—when it lies just outside of a critical interval around zero whose size depends on η —the response becomes quasiperiodic and its power spectrum has features at frequencies f of the form
where p and q are integers.
If pν and qν_{0 }are close, the system's nonlinearity will permit what is called synchronization of order p/q: there will be periodic responses in which the intrinsic oscillator goes through p cycles for every q periods of the forcing, so that the response has frequency
This idea allows us to give precise definitions for the terms suband superharmonic resonance used above: the former corresponds to resonances where p = 1 and q > 1 in (4), while the latter corresponds to q = 1 and p > 1. Finally, when both η and the difference (ν ν_{0}) are large, the system's response can become chaotic, so that the features in its power spectral density bear no simple relationship to the frequencies ν and ν_{0}.
Figure 7 which is a heat map of the power spectral density of the response generated by relatively strong sinusoidal forcing of the from (2), illustrates many of the behaviours described above. Consider first the vertical strip with ν/ν_{0 }≈ 1. In this region the sinusoidal modulation entrains the NFκB system essentially completely and so the heat map resembles the corresponding region in the left panel Figure 4: the power in the response is concentrated along lines corresponding to harmonics of the forcing frequency. Away from the region ν/ν_{0 }≈ 1 the power spectrum of the response is considerably more complex: the strong horizontal bands in Figure 7 which occur at integer multiples of ν_{0}, show that there is substantial power at the NFκB system's natural frequency and its harmonics. Additionally, the nonlinearity of the system means that the response has power at frequencies given by
Figure 7. Nonlinear resonances to sinusoidal forcing. A heat map showing the power spectral density of the response of Horton's model to forcing of the form (2) with ε = η = 0.5 and frequencies in the range 0 ≤ ν ≤ 4.5 × ν_{0}.
which gives rise to the network of lines with slope of ±1. Finally, Figure 7 also exhibits subharmonic resonance: vertical strips for which ν ≈ qν_{0 }(with q a whole number) show power concentrated at the forcing frequency, but also at frequencies f = ν/p, where p is a whole number. This is perhaps clearest in the strip ν ≈ 2ν_{0}, where the strongest spectral feature lies along the line f = ν/2 ≈ ν_{0}.
The two panels of Figure 8 are analogues of Figure 7 but with especially strong (η = 1, left panel) or weak (η = 0.1, right panel) modulation. The qualitative features are much the same, though it is interesting to note that in the limit of very strong modulation—when η = 1 and so T_{R }(t) vanishes once per forcing period—the system is very strongly entrained by the forcing and does not show much power at its natural frequency ν_{0 }or its harmonics until the forcing frequency ν > 1.5 ν_{0}. By contrast, when the modulation is weak evidence of modal interaction is also weak, with very narrow modelocking regions near forcing frequencies of the form ν ≈ mν_{0}.
Figure 8. Weak and strong modulations. Heat maps showing the power spectral density of the responses to forcing of the form (2) with ε = 0.5 and η = 1 (left) or η = 0.1 (right).
As Figures 7 and 8 illustrate, the susceptibility of our model NFκB system to resonance with sinusoidal modulations depends strongly on the amplitude of the modulation: Figure 9 provides a quantitative survey of this phenomenon. The shaded regions are examples of Arnol'd tongues: their precise shapes can be calculated with the methods outlined in the Appendix. When the modulation frequency and amplitude (ν, η) lie inside these tongues, the response of the system will be periodic, with a period τ that is an integer multiple of the modulation period 1/ν and that also lies close to an integer multiple of the natural Period 1/ν_{0}:
Figure 9. Arnol'd tongues. The qualitative behaviour of Horton's NF κB system when subjected to forcing of the from 2 with ε = 0.5: when the modulation frequency ν and amplitude η lie inside the shaded regions the response is periodic with a period that is related, by a relation like (5), to both the modulation's period and the system's natural period.
where p and q are whole numbers. Generally speaking these periodic responses—which are a nonlinear generalization of the familiar phenomenon of resonance in linear systems—are easiest to excite when the numbers p and q are small: the tongues in Figure 9 are labelled by ratios q:p where p and q are as in (5). Although Figure 9 illustrates this story for the specific family of modulations with ε = 0.5, the qualitative picture is essentially the same for all values in the range of 0.366 < ε < 0.5: in all cases there will be narrow tongues of parameter combinations (ν, η) for which the response is periodic with a period given by a resonance relation like (5). For modulations whose parameters lie outside the tongues the qualitative behaviour of the response will be more complex, having a power spectrum qualitatively similar to those in Figure 7, with features at harmonics of the modulation frequency, at harmonics of the natural frequency ν_{0 }(which varies with ε) and at nearresonant combinations (3) of the two.
When T_{R }< 0.366, the unforced system is a damped oscillator and so, when subjected to periodic forcing, may exhibit resonant phenomena. By using the same numerical simulation and power spectral density analyses as we did for the case ε = 0.5, we find that the responses of the pathway to forcing with with T_{R }< 0.366 can be divided into two groups: those in which harmonic, suband superharmonic resonances can be observed (in the region 0 < ε ≤ 0.01), and those for which only harmonic and subharmonic resonance can be observed.
We have focussed here on periodic and quasiperiodic oscillatory interactions amenable to powerspectral analyses of the sort illustrated in Figures 4, 5, 6, 7, 8 and, for singlecell recordings, in the lower panels of Figure 1. But there is every reason to expect a much richer range of dynamical behaviour: Fonslet et al. [20], who applied forcing of the form (2) with η = 1 to a simplified NFκB model, saw evidence of a perioddoubling cascade as well as chaos and strange attractors. Complete analysis of these more complex dynamical regimes requires tools beyond the power spectrum, and so we defer their exploration to a future paper.
Conclusions
Given that intrinsically nonlinear chemical kinetics underpin cellsignalling networks, one shouldn't expect these systems to be linear and any analysis of experimental time series, whether in the time or frequency domain, must take this intrinsic nonlinearity into account. Our modelling studies suggest that coupling even the simplest, sinusoidal signal into the NFκB network can give rise to a host nonlinear phenomena, including harmonic, subharmonic and superharmonic resonances as well as quasiperiodic and even chaotic behaviour.
The simulation studies reported here used an externallyimposed oscillatory forcing as a proxy for interactions between the core NFκB feedback loop and other oscillatory networks. Our results suggest that interactions between the NFκB oscillator and other oscillatory pathways can give rise to extremely rich temporal signalling programs and so, perhaps, to many distinct patterns of expression for target genes.
Methods
Differential Equations
The model studied in this paper is specified by the following ODEs. External forcing via TNFα stimulation is represented by the function T_{r }(t), whose value ranges between 0 and 1, indicating the strength of stimulation as a fraction of the maximum possible.
Initial conditions
All the numerical experiments reported above began with the initial concentrations
while all others start at zero. We then integrated the ODEs through 2000 minutes of simulated time and used the result as an initial condition for the various TNFαstimulation studies.
Numerics
The bifurcation diagram in Figure 3 and the plot of the Arnol'd tongues in Figure 9 were prepared using Bard Ermentrout's XPPAUT [29], a front end for the powerful numerical bifurcation package AUTO [30]. All our other numerical work used MATLAB^{© }to integrate the ODEs, compute the power spectra and plot the figures.
Authors' contributions
YW did the numerical experiments reported here, basing her work on code originally written by CAH, who also developed the underlying model. PP and CAH provided biologicallyinformed advice on modelling, especially on how best to couple the model to an external signal. MRHW, DBK and DSB provided critical readings of the manuscript, as well scientific guidance throughout the project. MRM worked closely with YW on the computational side and, with YW, wrote the first draft and prepared the figures. All authors reviewed and approved the final draft.
Appendix: circle maps
Here we discuss a standard mathematical tool, the circle map, used to study the response of a nonlinear oscillator subjected to periodic forcing. The idea is to consider sufficiently weak forcing that the response is close to that of the unforced system, and then define a phase angle ϕ that is close to the phase of the corresponding unforced oscillator.
If ω_{0 }is the frequency of the unperturbed oscillator and ω is the forcing frequency then, after one period of the forcing, ϕ will have advanced by
where τ = 2π/ω is the period of the forcing, ε is a measure of its amplitude and G(ϕ) is a phasedependent function that characterizes the response.
The simplest circle map, studied here by way of illustration, is given by a function F : [0, 2π) → [0, 2π) of the form of
where α = (ω_{0}/ω), and 0 < ε ≪ 2π. The dynamics of a circle map can be characterized with a single parameter called the rotation number. For a given initial point ϕ_{0}, the rotation number is defined as the longterm average phase shift per one iteration (that is, per one period of the forcing):
where
One can show that the limit in the definition of the rotation number (7) exists and does not depend on the initial point ϕ_{0 }[31], p.102. The qualitative dynamics of the forced system are thus of two types: motions with rational rotation numbers and those with irrational ones. Further, the rotation number is retionalsay , with p ∈ ℕ and q ∈ ℤ^{+ }—if and only if there exists some ϕ_{0 }∈ [0, 2π) such that
That is, rational rotation numbers correspond to periodic dynamics: the orbit of ϕ_{0 }is called a (p, q) cycle. On the other hand, irrational rotation numbers correspond to quasiperiodic dynamics [27,31,32].
According to a theorem of Denjoy [33], if the rotation number is irrational there is a continuous, invertible change of coordinates h: [0, 2π) → [0, 2π), say, h(ϕ) = θ, such that
Where . This implies that forced oscillators whose corresponding circle maps have irrational rotation number never repeat periodically.
Results of Arnold [31] show that the qualitative behaviour of the circle map (6) is stable against arbitrary small perturbations if and only if the rotation number is rational. We can thus expect that if F has a rational rotation number , there exists a region of parameter values (α, ε) such that the all the forced systems whose parameters lie in this region share the same rational rotation number: that is, they all have the same sorts of periodic response. Such regions of parameter values are called Arnol'd tongues.
A similar story holds for the sinusoidally forced NFκB system: when the timeaverage of the stimulation is sufficiently strong that the system supports periodic oscillations with frequency ω_{0}, and when the amplitude of the forcing is sufficiently weak that it is sensible to measure a phase ϕ with respect to that of the intrinsic oscillator, then one can construct numerically a circle map of the form
where T = 2π/ω is the period of the forcing. Although the function G' in (9) is not as simple as the the G (ϕ) in (6), the map F' still has (p, q) periodic responses and a corresponding system of Arnol'd tongues: Figure 9 shows examples.
We'll conclude this Appendix by explaining why, at least for the simple circle map (6), the Arnol'd tongues are wedgeshaped. Consider the case where F has rational rotation number . Then, as mentioned above, there exists a periodic point ϕ_{0 }satisfying (8). If we expand F^{q }(ϕ_{0}) in powers of ε we find
Consider the case where α is very close to a rational number: with β ≪ 1/q. Then
Where
It is not hard to show that G_{P/q }is bounded and continuous when regarded as a function of ϕ_{0 }and so attains its maximum and minimum values on [0, 2π]. Thus, for each fixed ε, there is an interval of β on which G_{P/q }= 0 for some ϕ ∈ [0, 2π]. Now consider the way in which the end points of this interval depend on ε. The case q = 1 and p = 0 is especially clear: Eqn. (10) becomes
which has some solutions ϕ ∈ [0, 2π] provided that 2πβ ≤ ε or β ≤ ε/2π.
The values q = 1 and p = 1 correspond to the case where the periodic response has the same frequency as the forced system and thus undergoes one complete phaserotation per period of the forcing. The analysis sketched above shows that this particular form of modelocking will pertain for (ε, α) in a wedgeshaped region whose width increases linearly with forcing amplitude ε. For more general modelocking resonances, where q ≥ 2, Arnol'd showed in [31] that the width of the tongue satisfies:
Acknowledgements
CAH and YW were supported by BBSRC grant BBD0088081. PP and MRHW were supported by MRC Grant G0500346. PP, MRM, DSB and MRHW were supported by BBSRC grants BBF005381 & BBF00561X1.
References

Hayden M, Ghosh S: Signaling to NFκB. [http://dx.doi.org/10.1101/gad.1228704] webcite

Bonizzi G, Karin M: The two NFκB activation pathways and their role in innate and adaptive immunity.
Trends in Immunology 2004, 25(6):280288. PubMed Abstract  Publisher Full Text

Pasparakis M, Luedde T, SchmidtSupprian M: Dissection of the NFκB signalling cascade in transgenic and knockout mice.

Pahl HL: Activators and target genes of Rel/NFκB transcription factors. [http://www.nature.com/onc/journal/v18/n49/abs/1203239a.html] webcite
Oncogene 1999, 18(49):68536866. PubMed Abstract  Publisher Full Text

Ghosh S, May M, Kopp E: NFκB and rel proteins: Evolutionarily conserved mediators of immune responses.
Annual Review of Immunology 1998, 16:225260. PubMed Abstract  Publisher Full Text

Scott ML, Fujita T, Liou HC, Nolan GP, Baltimore D: The p65 subunit of NFkappa B regulates I kappa B by two distinct mechanisms. [http://genesdev.cshlp.org/content/7/7a/1266.abstract] webcite

Song HY, Rothe M, Goeddel DV: The tumor necrosis factorinducible zinc finger protein A20 interacts with TRAF1/TRAF2 and inhibits NFkappaB activation.
Proceedings of the National Academy of Sciences of the United States of America 1996, 93(13):67216725. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Wertz IE, O'Rourke KM, Zhou H, Eby M, Aravind L, Seshagiri S, Wu P, Wiesmann C, Baker R, Boone DL, Ma A, Koonin EV, Dixit VM: Deubiquitination and ubiquitin ligase domains of A20 downregulate NFκB signalling.
Nature 2004, 430:694699. PubMed Abstract  Publisher Full Text

Pigolotti S, Krishna S, Jensen MH: Oscillation patterns in negative feedback loops.
Proceedings of the National Academy of Sciences 2007, 104(16):65336537. Publisher Full Text

Thomas R: Logical analysis of systems comprising feedback loops.
Journal of Theoretical Biology 1978, 73(4):631656. PubMed Abstract  Publisher Full Text

Hoffmann A, Levchenko A, Scott M, Baltimore D: The I kappa BNFkappa B signaling module: Temporal control and selective gene activation.
Science 2002, 298(5596):12411245. PubMed Abstract  Publisher Full Text

Nelson DE, Ihekwaba AEC, Elliott M, Johnson JR, Gibney CA, Foreman BE, Nelson G, See V, Horton CA, Spiller DG, Edwards SW, McDowell HP, Unitt JF, Sullivan E, Grimley R, Benson N, Broomhead DS, Kell DB, White MRH: Oscillations in NFkappa B signaling control the dynamics of gene expression.
Science 2004, 306(5696):704708. PubMed Abstract  Publisher Full Text

Horton CA: Computational Modelling of Cell Signalling Pathway Dynamics. PhD thesis. School of Biological Sciences, University of Liverpool; 2006.

Garside H, Stevens A, Farrow S, Normand C, Houle B, Berry A, Maschera B, Ray D: Glucocorticoid ligands specify different interactions with NFkappa B by allosteric effects on the glucocorticoid receptor DNA binding domain.
Journal of Biological Chemistry 2004, 279(48):5005050059. PubMed Abstract  Publisher Full Text

Ikeda A, Sun X, Li Y, Zhang Y, Eckner R, Doi T, Takahashi T, Obata Y, Yoshioka K, Yamamoto K: p300/CBPdependent and independent transcriptional interference between NFkappa B RelA and p53.
Biochemical and Biophysical Research Communications 2000, 272(2):375379. PubMed Abstract  Publisher Full Text

Perkins ND: Integrating cellsignalling pathways with NFκkB and IKK function.
Nature Reviews Molecular Cell Biology 2007, 8:4962. PubMed Abstract  Publisher Full Text

Salminen A, Ojala J, Huuskonen J, Kauppinen A, Suuronen T, Kaarniranta K: Interaction of agingassociated signaling cascades: Inhibition of NFkappa B signaling by longevity factors FoxOs and SIRT1.
Cellular And Molecular Life Sciences 2008, 65(78):10491058. PubMed Abstract  Publisher Full Text

Wadgaonkar R, Phelps K, Haque Z, Williams A, Silverman E, Collins T: CREBbinding protein is a nuclear integrator of nuclear factorkappa B and p53 signaling.
Journal of Biological Chemistry 1999, 274(4):18791882. PubMed Abstract  Publisher Full Text

Ashall L, Horton CA, Nelson DE, Paszek P, Harper CV, Sillitoe K, Ryan S, Spiller DG, Unitt JF, Broomhead DS, Kell DB, Rand DA, See V, White MRH: Pulsatile stimulation determines timing and specificity of NFkappaBdependent transcription.
Science 2009, 324(5924):242246. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Fonslet J, RudPetersen K, Krishna S, Jensen MH: Pulses and Chaos: Dynamical Response in a Simple Genetic Oscillator.
International Journal of Modern Physics B 2007, 21(2324):40834090. Publisher Full Text

Krishna S, Jensen MH, Sneppen K: Minimal model of spiky oscillations in NFκB signaling.
PNAS 2006, 103(29):1084010845. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Nelson DE: Charaterisation and Manipulation of NFκB Signal Oscillations in Single Living Cells. PhD thesis. School of Biological Sciences, University of Liverpool; 2005.

Matalka KZ, Tutunji MF, AbuBaker M, Abu Baker Y: Measurement of protein cytokines in tissue extracts by enzymelinked immunosorbent assays: Application to lipopolysaccharideinduced differential milieu of cytokines.
Neuroendocrinology Letters 2005, 26(3):231236. PubMed Abstract

Prabha C, Jalapathy KV, Matsa RP, Das SD: Role of TNFalpha in host immune response in tuberculous pleuritis.

Wang Y, Paszek P, Horton CA, Hong Y, White MRH, Kell DB, Muldoon MR, Broomhead DS: A systematic survey of the response of a model NFκB signalling pathway to TNFα stimulation. [http://eprints.ma.man.ac.uk/1404/] webcite
MIMS EPrint 2009.104 Manchester Institute for Mathematical Sciences, The University of Manchester; 2010.

Wiggins S: Introduction to Applied Nonlinear Dynamical Systems and Chaos. New York; London: Springer Verlag; 1990.

Pikovsky A, Rosenblum M, Kurths J: Synchronization : A Universal Concept in Nonlinear Sciences. Cambridge, UK: Cambridge University Press; 2003.

Nayfeh AH, Mook DT: Nonlinear Oscillators. New York: John Wiley & Sons; 1979.

Ermentrout B: Simulating, Analyzing, and Animating Dynamical Systems: A Guide to XPPAUT for Researchers and Students. Philadelphia: SIAM; 2002.

Doedel EJ, Paffenroth RC, Champneys AR, Fairgrieve TF, Kuznetsov YA, Oldeman BE, Sandstede B, Wang XJ: AUTO2000 : Software for continuation and bifurcation problems in ordinary differential equations. In Tech rep. California Institute of Technology, Department of Applied and Computational Mathematics, Pasadena CA 91125; 2000.

Arnold VI: Geometrical Methods in the Theory of Ordinary Differential Equations. 2nd edition. New York; London: Springer Verlag; 1988.

Arrowsmith DK, Place CM: An introduction to Dynamical Systems. Cambridge, UK: Cambridge University Press; 1990.

Denjoy A: Sur les courbes definies par les equations differentielles a la surface du tore.