Spatio-temporal correlations in biochemical networks

Biochemical networks not only operate in time, but also in space. Arguably the best example is embryonic development where spatial gradients in the concentrations of morphogens – substances governing the patterning of tissue – encode positional information for differentiating cells. Also within cells, compartmentalization, scaffolding and localized interactions are actively exploited to enhance the regulatory function of biochemical networks. In this research line we aim to elucidate the spatial design principles of biochemical networks.

Gene expression

Diffusion can be a major source of noise
Experiments in recent years have beautifully demonstrated that gene expression is often highly stochastic. The theoretical analyses of these experiments typically employ the zero-dimensional chemical master equation. This approach takes into account the discrete character of the reactants and the probabilistic nature of chemical reactions. It does assume, however, that the cell is a ‘‘well-stirred’’ reactor, in which the particles are uniformly distributed in space at all times; the reaction rates only depend upon the global concentrations of the reactants and not upon the spatial positions of the reactant molecules. Yet, to react, reactants first have to move toward one another. They do so by diffusion, or in the case of eukaryotes, by a combination of passive and active transport. Both processes are stochastic in nature, and this could contribute to the noise in the network.

We studied by Green’s Function Reaction Dynamics [1,2] the effect of the diffusive motion of repressor molecules on the noise in mRNA and protein levels for a gene that is under the control of a repressor. We found that spatial fluctuations due to diffusion can drastically enhance the noise in gene expression [3]. After dissociation from the operator, a repressor can rapidly rebind to the DNA. Our results show that the rebinding trajectories are so short that, on this timescale, the RNA polymerase (RNAP) cannot effectively compete with the repressor for binding to the promoter. As a result, a dissociated repressor molecule will on average rebind many times, before it eventually diffuses away. These rebindings thus lower the effective dissociation rate, and this increases the noise in gene expression; indeed, the increase in the noise is due to slower promoter-state fluctuations (Fig.1). Another consequence of the timescale separation between repressor rebinding and RNAP association is that the effect of spatial fluctuations can be described by a well-stirred, zero-dimensional model by renormalizing the reaction rates for repressor-DNA (un) binding. Our results thus support the use of well-stirred, zerodimensional models for describing noise in gene expression. Yet, we stress that the rate constants can be renormalized so effectively, because the timescale for repressor rebinding is well separated from that of RNAP binding. This could be due to an artifact of our model, in which we describe the promoter as a sphere. We expect that in a more detailed model in which the repressors and the RNAP are allowed to also move along the DNA, this separation of timescales brakes down, leading to a much more intricate interplay between the motion of the repressor and the RNAP.

biochemical networks ten wolde SpatFluc

Fig.1 Noise in the expression of a gene that is under the control of a repressor protein. The figures on the left show time traces of the number of protein, NP, for different number of repressors, NR, and for a well-stirred (WS) model, in which it is assumed that the diffusion constant is infinitely large, and for a spatially resolved model, simulated with GFRD. The panel on the right shows the noise in the protein level as a function of the number of repressors NR; when NR is varied, the dissociation constant is adjusted such that the promoter occupancy remains unchanged; the nearly horizontal lines are the predictions from the well-stirred model, and the symbols are the results from the spatially resolved GFRD simulations. It is seen that when NR is large the WS and GFRD model predict similar protein concentration fluctuations; however, when NR is small, the WS model drastically underestimates the noise; now, diffusion of the repressor is the dominant source of noise in gene expression. In the spatially resolved model, a dissociated repressor molecule can rebind many times to the promoter, before it diffuses away into the cytosol. This lowers the effective dissociation time, leading to slower promoter-state fluctuations. These slow promoter-state fluctuations lead to bursts in gene expression as the time traces on the left show, especially when NR is small; for large NR, both the association and dissociation time are large (the dissociation constant is varied with NR so as to keep the promoter occupancy constant) and promoter-state fluctuations are fast and hence do not contribute to the noise in gene expression.

Spatial averaging during embryonic development

Diffusion can enhance robustness by washing out bursts in gene expression Embryonic development is driven by orderly, spatial patterns of gene expression that assign each cell in the embryo its particular fate. While gene expression is often highly stochastic, embryonic development is exceedingly precise. A vivid example is the Bicoid-Hunchback system in the early Drosophila embryo. Shortly after fertilization, the morphogen protein Bicoid (Bcd) forms an exponential concentration gradient along the anterior-posterior axis of the embryo, which provides positional information for the differentiating nuclei. One of the target genes of Bcd is hunchback (hb), which is expressed in the anterior half of the embryo. The posterior boundary of the hb expression domain is very sharp: by cell cycle 13, the position of the boundary varies only by about one nuclear spacing (Gregor et al., Cell, 2007). This precision is higher than the best achievable precision for a time-averaging based readout mechanism of the Bcd gradient. Intriguingly, the study of Gregor et al. revealed that the Hb concentrations in neighboring nuclei exhibit spatial correlations and the authors suggest that this implies a form of spatial averaging enhancing the precision of the posterior Hb boundary. However, the mechanism for spatial averaging remained unclear.

We have analytically and numerically studied the Bcd-Hb system [4]. Our analysis reveals a simple, yet powerful mechanism for spatial averaging, which is based on the diffusion of Hb itself. Hb diffusion between neighboring nuclei reduces the super-Poissonian part of the noise in its concentration, with a factor that depends on the diffusion length of Hb and the dimensionality of the system. In essence, diffusion reduces the noise by washing out bursts in gene expression. This mechanism is generic, and applies not only to any developmental system, but also to any biochemical network in general. For example, if a signaling protein is activated at one end of the cell and then has to diffuse to another place to activate another system, e.g., the messenger CheY in bacterial chemotaxis, then our results show that the non-Poissonian noise in the activation of the signaling protein is washed out by diffusion; for this reason it may be beneficial to spatially separate the in- and output of a signaling pathway. Our analysis also reveals that, while Hb diffusion reduces the noise, it also lessens the steepness of its expression boundary. The interplay between these two antagonistic effects leads to an optimal diffusion constant of D≅0.1 µm2s-1 that maximizes the precision of the hb expression domain.

biochemical networks ten wolde SpatAve

Fig.2 Spatial averaging during embryonic development. Shortly after fertilization,Bcd forms an exponential gradient along the anterior-posterior axis of the Drosophila embryo (solid line left panel). One of the targets of Bcd is hb, which is expressed in the anterior half of the embryo. Due to noise in gene expression, the Hb concentration profile fluctuates – the panel on the left shows a number of instantaneous Hb profiles. The panel on the right shows the width of the Hb-expression boundary as a function of the Hb diffusion constant. It is seen that even though hb expression is noisy, the hb mid-boundary is surprisingly sharp: for the optimal Hb diffusion constant, the width is about one nuclear spacing, in agreement with experiment (Gregor et al., Cell, 2007). The optimum arises from the interplay of two opposing factors, namely Hb diffusion making the Hb boundary less steep, and Hb diffusion washing out the effects of bursts in gene expression; diffusion reduces the non-Poissonian part of the noise. The solid line in the right panel is the result of an analytical calculation, where the boundary width Δx is calculated from the ratio of the noise σ in hb expression at the boundary and the steepness of the gradient at the boundary, |grad[Hb]|, Δx=σ/|grad[Hb]| (see [2]). It is seen that the agreement with the simulation results is remarkably good.

MAPK signaling

Diffusion can qualitatively change the response of biochemical networks
Multisite covalent modification of proteins is omnipresent in eukaryotic cells. A well-known example is the mitogen-activated protein kinase (MAPK) cascade where, in each layer of the cascade, a protein is phosphorylated at two sites. It has long been known that the response of a MAPK pathway strongly depends on whether the enzymes that modify the protein act processively or distributively. A distributive mechanism, in which the enzyme molecules have to release the substrate molecules in between the modification of the two sites, can generate an ultrasensitive response and lead to hysteresis and bistability. However, the theoretical studies performed so far have largely been based on mean-field calculations in which the spatio-temporal correlations of the particles are ignored.

We have studied by GFRD [5] a dual phosphorylation cycle in which the enzymes act according to a distributive mechanism (Fig.3). We find that the response of this network can differ dramatically from that predicted by a mean-field analysis based on the chemical rate equations. In particular, rapid rebindings of the enzyme molecules to the substrate molecules after modification of the first site can markedly speed up the response, with slower diffusion leading to a faster response (Fig.4A). Moreover, rapid enzyme-substrate rebindings can lead to loss of ultrasensitivity (Fig.4B) and bistability (Fig4C). In essence, rebindings can turn a distributive mechanism into a processive one. Our results also show that slow ADP release by the enzymes can protect the system against rapid enzyme-substrate rebindings, thus enabling ultrasensitivity and bistability.

Biological systems that exhibit macroscopic concentration gradients or spatio-temporal oscillations are typically considered to be reaction-diffusion problems. The MAPK system that we have studied is spatially uniform at cellular length scales. Yet, our analysis reveals that its macroscopic behavior, such as the presence of bistability, cannot be described by a mean-field analysis that ignores spatio-temporal fluctuations at molecular length scales. While it has long been known that spatio-temporal correlations at molecular length scales can renormalize rate constants, we believe that our results are the first to show that spatio-temporal correlations at molecular length scales can also qualitatively change the average, macroscopic behavior of a biological system that is uniform at cellular length scales. This also underscores the importance of particle-based modeling of biological systems in time and space: systems that appear to be well-stirred may have to be treated as reaction-diffusion systems.

biochmical networks ten wolde MAPK scheme

Fig.3 A cartoon of the MAPK pathway. The substrate X can be phosphorylated at two sites, which is controlled by the action of two antagonistic enzymes, Ea and Ed.  The enzymes operate according to a distributive mechanism, which means that they have to release the substrate in between the modification of the two sites.After the enzymes have modified a site, they are inactive for a reactivation time τrel. During this time they cannot bind a substrate


biochemical networks ten wolde MAPK Res




Fig.4 Enzyme-substrate rebindings can change the macroscopic response of a MAPK pathway. (A) The response time as a function of the diffusion constant D, for two different enzyme reactivation times (see Fig.3); blue and green lines are results from the GFRD simulations, while the black curves correspond to calculations based on the mean-field, deterministic chemical rate equations. It is seen that slower diffusion can lead to a faster response, especially when the enzyme-reactivation time is short (green line) (see Fig.3). (B) The steady-state input-output relation for different diffusion constants. When the diffusion constant is large, the response curve approaches that as predicted by the mean-field chemical rate equations (black solid line). For smaller diffusion constants, the response curve approaches that of a processive scheme (dashed blue line). (C) The concentration of doubly phosphorylated substrate Xpp as a function of the enzyme reactivation time τrel; the simulations were started from different initial conditions, and plotted are the final steady-state concentrations. When the enzyme reactivation time is large, the system is bistable, leading to two different steady-state concentrations of Xpp; however, for smaller enzyme-reactivation times, the bistability is lost. The results shown in these panels can be understood by noticing that rapid enzyme-substrate rebindings become more likely as the diffusion constant and/or the enzyme reactivation time decreases (see Fig.5).

biochemical networks ten wolde MAPK Preb

Fig.5 Enzyme-substrate rebindings. Plotted are the distributions of the second-association time, which is the time it takes for an enzyme to bind a substrate that has been modified once, thus the waiting time for the reaction Ea + Xp -> EaXp. The left panel corresponds to the scenario where the enzyme reactivation time τrel is short, while the right panel corresponds to the scenario where the enzyme reactivation time τrel is long (longer than the time τmol for molecules to cross a molecular distance). The peak at short times corresponds to rapid enzyme-substrate rebindings, where an enzyme and a substrate molecule that have just dissociated rapidly re-associate; the corner at long times corresponds to random encounters between enzyme and substrate molecules in the bulk. It is seen that rapid enzyme-substrate rebindings become more likely when the diffusion constant and/or the enzyme-reactivation time are decreased: the probability than an enzyme and a substrate molecule that have just dissociated from each other are still in close proximity by the time that the enzyme has become active again increases when this time becomes shorter (the molecules have had less time to diffuse away from each other) or the diffusion constant decreases (the molecules move slower and will thus cross a shorter distance during this time). Since these enzyme-substrate rebindings are much faster than random encounters, slower diffusion can lead to a faster response (left panel Fig.4). Moreover, since the probability of a rebinding event does not depend on the concentration, rebindings can turn a distributive scheme into a processive one (middle panel Fig.4). Whether a scheme is distributive or processive does not directly depend on the number of collisions required to fully phosphorylate the substrate, as commonly believed, but whether both modification sites are phosphorylated by the same enzyme molecule. While in a purely processive scheme the enzyme and substrate molecules remain physically connected to each other in between the modification of both sites, in a distributive scheme that acts processively the enzyme and substrate molecule do dissociate but then rapidly rebind.


  1. J. S. van Zon and P. R. ten Wolde , Simulating biochemical networks at the particle level and in time and space: Green’s function reaction dynamics , Physical Review Letters 94, 128103: 1- 4 (2005)
  2. J. S. van Zon and P. R.ten Wolde , Green’s-function reaction dynamics: a particle-based approach for simulating biochemical networks in time and space , Journal of Chemical Physics 123, 234910: 1-16 (2005)
  3. Jeroen S. van Zon , Marco J. Morelli , Sorin Tanase-Nicola and Pieter Rein ten Wolde , Diffusion of transcription factors can drastically enhance the noise in gene expression , Biophysical Journal 91, 4350- 4367 (2006)
  4. Thorsten Erdmann and Pieter Rein ten Wolde , On the role of spatial averaging for the precision of gene expression patterns , Physical Review Letters 103, 258101 1- 4 (2009)
  5. Koichi Takahashi , Sorin Tanase-Nicola and Pieter Rein ten Wolde , Spatio-temporal correlations can drastically change the response of a MAPK pathway , PNAS 107, 6:2473-2478 (2010).