Stochastic switching between multistable oscillation patterns of the Min-system

For the particle-based simulation, we use the reaction scheme of the minimal model for cooperative attachment [39, 47, 49, 52]. The model uses the following interactions between Min-proteins and the inner bacterial membrane (schematically shown in figure 1(A)). Freely diffusing cytoplasmic ${\mathrm{MinD}}_{\mathrm{ATP}}$ can bind to the membrane with a rate constant ${\sigma }_{D}$

$$\begin{eqnarray}&&{\mathrm{MinD}}_{\mathrm{ATP}}\mathop{\longrightarrow }\limits^{{\sigma }_{D}}{\mathrm{MinD}}_{\mathrm{bound}}.\end{eqnarray} \tag{ 1a }$$

${\mathrm{MinD}}_{\mathrm{ATP}}$ binds preferably to high density ${\mathrm{MinD}}_{\mathrm{bound}}$ regions (cooperative MinD binding)

$$\begin{eqnarray}&&{\mathrm{MinD}}_{\mathrm{ATP}}+{\mathrm{MinD}}_{\mathrm{bound}}\mathop{\longrightarrow }\limits^{{\sigma }_{{dD}}}2\ {\mathrm{MinD}}_{\mathrm{bound}}.\end{eqnarray} \tag{ 1b }$$

Membrane-bound MinD also recruits cytoplasmic MinE to the membrane with rate ${\sigma }_{E}$, creating a ${\mathrm{MinDE}}_{\mathrm{bound}}$ complex

$$\begin{eqnarray}&&{\mathrm{MinD}}_{\mathrm{bound}}+\mathrm{MinE}\mathop{\longrightarrow }\limits^{{\sigma }_{E}}{\mathrm{MinDE}}_{\mathrm{bound}}.\end{eqnarray} \tag{ 1c }$$

All membrane-bound proteins diffuse in the plane of the membrane, but with a much smaller diffusion constant than in the bulk.

MinE attachment activates ATP hydrolysis of MinD. The hydrolysis of ATP to ADP breaks up the membrane-bound complex and releases ${\mathrm{MinD}}_{\mathrm{ADP}}$ and MinE back into the cytoplasmic bulk with rate constant ${\sigma }_{\mathrm{off}}$

$$\begin{eqnarray}&&{\mathrm{MinDE}}_{\mathrm{bound}}\mathop{\longrightarrow }\limits^{{\sigma }_{\mathrm{off}}}\mathrm{MinE}+{\mathrm{MinD}}_{\mathrm{ADP}}.\end{eqnarray} \tag{ 1d }$$

Finally ${\mathrm{MinD}}_{\mathrm{ADP}}$ exchanges ADP by another ATP molecule (nucleotide exchange) with the rate λ

$$\begin{eqnarray}&&{\mathrm{MinD}}_{\mathrm{ADP}}\mathop{\longrightarrow }\limits^{\lambda }{\mathrm{MinD}}_{\mathrm{ATP}}.\end{eqnarray} \tag{ 1e }$$

This completes the reaction cycle. In table 1 we list our parameter values as set A. For comparison, we also list parameter values used in other studies (set B in [30] and set C in [47, 49]).

Table 1.  Parameters sets used for simulations of the Min-system. Set A is the main parameter set used here following [39, 52]. Parameter set B is taken from [30] and C from [47, 49].

Parameter	Set A	Set B	Set C	Unit	Description (reaction type)
DD	2.5	16	16	μm2 s−1	Bulk diffusion coefficient of MinD
DE	2.5	10	10	μm2 s−1	Bulk diffusion coefficient of MinE
${D}_{\mathrm{bound}}$	0.01	0.013	0.013	μm2 s−1	Membrane diffusion coefficient
λ	0.5	1	6	${{\rm{s}}}^{-1}$	First order, unimolecular
${\sigma }_{D}$	0.025	0.075	0.1	μm s−1	First order, membrane attachment
${\sigma }_{{dD}}$	0.0149	0.054	0.1	μm3 s−1	Second order, bimolecular
${\sigma }_{E}$	0.093	0.254	0.435	μm3 s−1	Second order, bimolecular
${\sigma }_{\mathrm{off}}$	0.7	0.33	0.5	${{\rm{s}}}^{-1}$	First order, unimolecular
cD	0.797	0.85	1.0	μM	Total MinD concentration
cE	0.207	0.31	0.5	μM	Total MinE concentration

⁴In [30] the unit for these bimolecular reactions is specified as $\mu {{\rm{m}}}^{2}\,{{\rm{s}}}^{-1}$.

We use custom-written code to simulate the stochastic dynamics of the Min-system with very good statistics. For all simulations we use a fixed discrete time step of ${\rm{\Delta }}t={10}^{-4}\ {\rm{s}}$. During every time step each particle is first propagated in space. Thereafter every particle can react according to the previously introduced Min reaction scheme (1a)–(1e). The movement of both free and membrane-bound particles is realized through Brownian dynamics. Individual molecules are treated as point-like particles without orientation. Therefore we can monitor the propagation separately for each Cartesian coordinate. During a simulation step of ${\rm{\Delta }}t$ the displacements of the diffusing particles with diffusion constant D are drawn from a Gaussian distribution with standard deviation ${\sigma }_{x}=\sqrt{2D{\rm{\Delta }}t}$ [53] such that

$$\begin{eqnarray}&&x(t+{\rm{\Delta }}t)=x(t)+{X}_{{\rm{G}}},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{p}_{{X}_{{\rm{G}}}}(x)=\displaystyle \frac{1}{\sqrt{4\pi D{\rm{\Delta }}t}}\exp \left(-\displaystyle \frac{{x}^{2}}{4D{\rm{\Delta }}t}\right),\end{eqnarray} \tag{ 3 }$$

where ${p}_{{X}_{{\rm{G}}}}$ is the probability distribution of X_G. The same update step is used for the y and z direction. Free particles in the bulk of the simulated volume undergo three-dimensional diffusion with reflective boundary conditions at the borders of the simulation volume. The membrane-bound particles perform a two-dimensional diffusion on the membrane with a much smaller diffusion constant ${D}_{\mathrm{bound}}$ (compare table 1). Membrane-bound particles are allowed to diffuse between different membrane areas that are in contact with each other.

The different reactions in the Min reaction scheme (1a)–(1e) can be classified into three different types (more details on the corresponding implementations are given in the supplementary information). The first type considered here are first order reactions without explicit spatial dependence. The conversion of ${\mathrm{MinD}}_{\mathrm{ADP}}$ to ${\mathrm{MinD}}_{\mathrm{ATP}}$ and the unbinding of the ${\mathrm{MinDE}}_{\mathrm{bound}}$ complex from the membrane are of this type. Such reactions are treated as a simple Poisson process. For a reaction rate κ, the probability to react during a time step ${\rm{\Delta }}t$ is given by

$$\begin{eqnarray}&&{p}_{\kappa }=1-\exp (-\kappa {\rm{\Delta }}t).\end{eqnarray} \tag{ 4 }$$

The second type is also a first order reaction, but with confinement to a reactive area at a border of the simulated volume. The membrane attachment of ${\mathrm{MinD}}_{\mathrm{ATP}}$ proteins is a reaction of this type. For a given reaction rate σ, we implement these reactions by allowing particles that are closer to the membrane than $d=0.02\,\mu {\rm{m}}$ to attempt membrane attachment with a Poisson rate $\kappa =\sigma /d$. This results in a reaction probability of

$$\begin{eqnarray}&&{p}_{\sigma }=1-\exp \left(-\displaystyle \frac{\sigma }{d}{\rm{\Delta }}t\right).\end{eqnarray} \tag{ 5 }$$

The last reaction type is a second order reaction between free and membrane-bound particles. The cooperative recruitment of cytosolic ${\mathrm{MinD}}_{\mathrm{ATP}}$ and MinE to membrane-bound MinD are reactions of this third type. In our simulation we adopt the algorithm implemented in the software package Smoldyn, which has been used earlier to simulate the Min-system [52]. This algorithm is based on the Smoluchowski framework in which two particles react upon collision [54]. However, the classical treatment by Smoluchowski only considers diffusion-limited reactions and therefore assumes instantaneous reactions upon collision. In order to take finite reaction rates into account, one imposes a radiation boundary condition [55–57]. From the diffusion constant D, the reaction rate σ and the simulation time step ${\rm{\Delta }}t$, a reaction radius r_σ is calculated [58]. Whenever a freely diffusing particle comes within the distance of r_σ to a membrane-bound particle, the free particle reacts. For intermediate values of ${\rm{\Delta }}t$ (such as the time step of ${10}^{-4}\ {\rm{s}}$ that we use for the Min-system) the value of ${r}_{\sigma }(D,\sigma ,{\rm{\Delta }}t)$ is obtained numerically [58]. Those numerical values are taken from the Smoldyn software. For example, for parameter set A the reaction radius for the rate ${\sigma }_{{dD}}$ is ${r}_{{\sigma }_{{dD}}}=0.0091\,\mu {\rm{m}}$, and for ${\sigma }_{E}$ it is ${r}_{{\sigma }_{E}}=0.0179\,\mu {\rm{m}}$.

In our simulations we use rectangular reaction compartments. We considered three different membrane setups as illustrated in figure 1(B). To mimic in vitro experiments, where substrates or open compartments are functionalized with a membrane layer [26–28], we place the reactive membrane at the bottom (geometry A) or at the side walls and the bottom of the simulation compartment (geometry B). To simulate rectangular shaped E. coli cells, inspired by the cell sculpting approach from [30, 47], fully membrane-covered volumes are used (geometry C). We refer to the long side of the lateral extension as the major or the x-axis, and the smaller side as minor or y-axis, and accordingly consider the compartment height to extend in the z-direction, aligning the rectangular geometry perpendicular with the coordinate frame.

In our simulations we investigate a wide range of compartment dimensions. For a simulation box with a length of 10 μm, width of 5 μm and height of 0.5 μm, we use 6003 ${\mathrm{MinD}}_{\mathrm{ATP}}$ particles, 6003 ${\mathrm{MinD}}_{\mathrm{ADP}}$ particles and 3124 MinE particles as initial condition [52]. These particle numbers amount to a total MinD concentration of $0.797\,\mu {\rm{M}}$ and a MinE concentration of $0.207\,\mu {\rm{M}}$. For other simulation compartment sizes we scale the particle numbers linear with the volume, since in experiments E. coli bacteria typically have a constant Min-protein concentration [30].

In our simulations of the Min-system different oscillation patterns emerge along the major or minor axis of the simulation compartment. In order to analyze the frequency of different modes and the stability of the oscillations in the large amount of simulation data, an oscillation mode recognition algorithm is needed. Therefore we monitor the MinD protein densities at the poles of the different axes over time. To determine the axis along which the oscillation takes place, we compare the Fourier transformation of the normalized densities over time (${\rho }_{{t}_{i}}$ where i denotes the discretized time resolution). If an oscillation takes place, there is a dominant peak in the Fourier spectrum and the overall maximal amplitude of the Fourier spectrum is significantly higher than the one from the non-oscillating axis. The same Fourier spectrum is also used to determine the oscillation period T. To differentiate between pole-to-pole oscillations and striped oscillations of a given axis in the system, we extract the phase difference between the density oscillation at the poles of the cell.

When identifying switching events, one has to be careful because they have to be identified in local time windows and stochastic fluctuations might lead to temporal changes that might be mistaken to be mode switches. We therefore smoothen the data. In detail, we calculate the convolution C_i between the densities over time and a Gaussian time window G_i

$$\begin{eqnarray}&&{C}_{i}=\displaystyle \sum _{j}{\rho }_{{t}_{j}}{G}_{i-j}\,\quad \mathrm{where}\quad {G}_{i}=\displaystyle \frac{1}{\sqrt{2\pi }\omega }\exp \left(-\displaystyle \frac{{({i}\tau )}^{2}}{2{\omega }^{2}}\right).\end{eqnarray} \tag{ 6 }$$

Here τ is the time between successive density measurements and we set $\omega =100\,{\rm{s}}$ as width of the time window. The current oscillation mode is now determined from the convoluted densities C_i and assigned to the time ${i}\tau$. In this way, only switches are identified that persist for a sufficiently long time.

First we investigated the oscillations that emerge in geometry A with parameter set A using a rectangular simulation volume with dimensions $10\,\mu {\rm{m}}\times 5\,\mu {\rm{m}}\times 0.5\,\mu {\rm{m}}$ ($x,y,z$). With this particular choice the width of the system approximately matches the typical length of wild-type E. coli cells and the length of the system corresponds to the length of a grown E. coli cell which can roughly double in length before septum formation and division. As shown by the kymographs in figure 2 and in agreement with the results of Hoffmann and Schwarz [52], in our simulations two different oscillation modes occur (compare also supplemental movie S1). Note from the color legend that dark and light colors correspond to low and high concentrations, respectively, as used throughout this work. In the first mode the Min-proteins oscillate along the minor y-axis from one pole to the other (pole-to-pole oscillation). In the second mode the proteins oscillate along the major x-axis between the poles and the middle of the compartment (striped oscillation). The system stochastically switches between the two modes, sometimes via a short oscillation along the diagonal of the compartment. The mode switching behavior of the Min-system in large volumes is in agreement with the experimental results of  Wu et al [30] and cannot be analyzed completely with conventional PDE-models of the Min-oscillations because they do not account for the noise in the system leading to the stochastic switch.

Zoom In Zoom Out Reset image size

Detailed analyses of the pole-to-pole and the striped oscillations from figure 2 are shown in figures 3(A) and (B), respectively (compare also supplemental movies S2 and S3, respectively). First, we consider the pole-to-pole oscillations in figure 3(A). In the kymographs of the bound MinD and MinE proteins (top row figures) we see clusters of bound proteins that detach from the membrane beginning in the middle and from there move towards one of the poles of the compartment in an alternating way. The shapes of the bound MinD protein density clusters in the kymographs have a triangular form, in contrast to the line-like structures of the bound MinE proteins. Those density lines in the bound MinE kymograph show that the MinE proteins form a high density cluster in the middle of the cell which propagates to one of the compartment poles. This behavior is similar to the experimentally observed ringlike structures of MinE proteins in E. coli bacteria that travel from the middle to the poles of the cell, leading to the dissociation of MinD proteins from the membrane. The kymographs of the free particles (middle row in figure 3(A)) are averaged over all heights and have the inverse shape of the corresponding kymographs of the bound particles (top row in figure 3(A)). Where the density of bound particles is high, the density of free particles is low and vice versa. During the simulations the spatial density differences of both MinD and MinE are higher for membrane-bound particles than for free particles in the bulk. Therefore in the two bottom kymographs in figure 3(A), which are showing total particle densities, the pattern of membrane-bound particles is dominant.

Zoom In Zoom Out Reset image size

The kymographs of the striped oscillation in figure 3(B) have the same structure as the ones of the pole-to-pole oscillations. However, the edges of the bound Min-protein clusters in the kymographs, that indicate the traveling Min-waves, are curved, in contrast to the straight lines that we see for the pole-to-pole oscillations as shown in figure 3(A). The time-averaged density profiles of MinD proteins for the pole-to-pole and striped oscillations are shown in figures 3(C) and (D), respectively. As expected the density of the proteins is minimal between the oscillation nodes of the emerging standing wave patterns.

It is highly instructive to compare the time evolution of the MinD and MinE protein densities. In figure 4(A) we plot the time evolution of the particle densities of the transverse pole-to-pole oscillation at a fixed position $y=4.9\,\mu {\rm{m}}$, which is at the edge of the minor axis along which the oscillation takes place. The shape of the transient density profiles is similar to experimentally observed density profiles of traveling Min-protein waves on flat membrane surfaces [24, 25]. The period of both oscillations modes was $T=(33.8\pm 0.1){\rm{s}}$, which is in agreement with the results of Huang et al [39].

Zoom In Zoom Out Reset image size

To analyze the influence of the bulk volume on the oscillations, we have also monitored how the MinD and MinE densities changes at different heights above the membrane. For this we again use geometry A, but now with a bottom surface of only $6\,\mu {\rm{m}}\times 3\,\mu {\rm{m}}$, which robustly produces longitudinal pole-to-pole oscillations along the major axis. For a compartment height of 4 μm the volume was divided in three layers, and for each layer the mean density on one side of the major axis was plotted over time as shown in figure 4(B). We see that the largest density changes take place in the layer directly above the membrane, however, the oscillations persist up to the top layer even for the highest compartment with $z=4\,\mu {\rm{m}}$. Strikingly, the MinD-density is always much higher than the one of MinE. Furthermore we reduced the compartment height to 0.2 μm. Interestingly, the pole-to-pole oscillations were still present. This implies that the bulk of the simulation volume has only a mild effect on the oscillations in geometry A, despite the fact that there are appreciable density variations in the bulk.

We next checked that our model indeed includes the essential elements for the emergence of Min-oscillations. Because the MinD-switching between bulk and membrane is clearly indispensible, here we consider the relevance of diffusion on the membrane and of cooperative recruitment. All simulations in this section are carried out using geometry A in a $5\,\mu {\rm{m}}\times 2.5\,\mu {\rm{m}}\times 0.5\,\mu {\rm{m}}$ compartment geometry. We chose this geometry since we expect it to give rise to stable longitudinal pole-to-pole oscillations and hence can closely monitor any deviations from this default oscillation mode. All other parameters were kept fixed and we used parameter set A for this test.

Figure 5(A) shows three kymographs of the system along the major axis. Switching off membrane diffusion entirely (first panel) leads to a loss of stable oscillation patterns. Instead small striped-like MinD-patches emerge erratically on the membrane. The second panel in figure 5(A) shows the default diffusion coefficient of ${D}_{\mathrm{bound}}$ = 0.01 $\mu {{\rm{m}}}^{2}\,{{\rm{s}}}^{-1}$ where stable pole-to-pole oscillations emerge along the major axis. By further increasing the membrane-diffusion-coefficient the pole-to-pole oscillations still robustly emerge but increasingly smear out. This behavior is most clearly illustrated by looking at the time-averaged density profiles as shown in figure 5(C). The inset in figure 5(C) shows that without membrane diffusion, no oscillations emerge at all (square symbol), while a slow diffusivity as used in parameter set A seems optimal.

Zoom In Zoom Out Reset image size

In a similar fashion we analyzed the influence of the cooperative membrane recruitment of MinD. Figure 5(B) shows three kymographs for no cooperative recruitment (${\sigma }_{{dD}}=0$), the default value from parameter set A and a two-fold increase, respectively. Without cooperativity in this process no oscillations emerge at all. The second panel shows again the default oscillation mode while already a two-fold increase also leads to unstable behavior without any patterns emerging. This sensitivity with respect to the cooperative recruitment rate ${\sigma }_{{dD}}$ is also summarized in figure 5(D), which illustrates that only in a vary narrow range stable oscillations emerge. Although a complete parameter scan is out of the question at the current stage for reasons of computer time, we conclude that stable Min-oscillations emerge only for certain parameter values and that parameter set A performs very well in this respect.

As we have seen in the preceding section, the simulation compartment height has little influence on the emerging oscillation modes when the membrane only covers the bottom area of the simulation compartment (geometry A). However, when we extend the membrane to cover bottom and side walls of the compartment (geometry B), we find that the oscillation modes change with increasing wall height. First, we consider a compartment that only exhibits longitudinal pole-to-pole oscillations along the major axis (wall height of 0.5 μm and bottom area of $5\,\mu {\rm{m}}\times 2.5\,\mu {\rm{m}}$). In contrast to the flat membrane geometry A, here the MinD protein density decreases in the vicinity of the membrane edges between the side walls and the bottom area, as shown in figure 6. This is due to the decreased volume per membrane area ratio in the vicinity of the membrane corners, leading to a decreased density of bound MinD proteins in these regions. This effect is clearly visible by comparing the time-averaged density profile in figure 6 with the one shown in 3(C).

Zoom In Zoom Out Reset image size

The full oscillation pattern is shown in 7(A) for the lowest height value (black label) and resembles a pole-to-pole oscillation for geometry A, but now between the two opposing walls along the major axis. Now we increase the wall height in increments of 0.5 μm. The oscillation changes at a wall height of $3.5\,\mu {\rm{m}}$ to a new oscillation mode (green label, compare also supplemental movie S4). There the proteins start to oscillate between the middle of the bottom area and the upper edges of the walls. This newly identified edge oscillation with one polar zone at the bottom and one top polar zone at each of the four side walls persists until the wall height reaches 9.5 μm. Thereafter the oscillation mode changes to a striped edge oscillation along the walls (yellow label).

In figure 7(B) we show results for a simulation volume that gives rise to striped oscillations along the major axis (bottom area of the length of 9 μm and width of 4 μm). At a wall height of $z=0.5\,\mu {\rm{m}}$ the simulation gives rise to longitudinal striped oscillations along the major axis or transverse pole-to-pole oscillations along the minor axis. This is the same kind of bistability that we already observed in geometry A using a compartment of $10\,\mu {\rm{m}}\times 5\,\mu {\rm{m}}\times 0.5\,\mu {\rm{m}}$. After increasing the wall height to $z=1.5\,\mu {\rm{m}}$ the pole-to-pole oscillations along the minor axis disappear and only striped oscillations along the major axis are observed. At the wall height of $z=3.5\,\mu {\rm{m}}$ the oscillation mode changes to a large pole-to-pole oscillation along the minor axis and no edge or striped oscillations were observed. Figures 7(C) and (D) show the kymographs for the membrane-bound proteins corresponding to figures 7(A) and (B), respectively. These results demonstrate that in geometry B, both the mode selection and the detailed shape of the polar regions can be controlled by compartment height.

Zoom In Zoom Out Reset image size

For an edge oscillation the free MinD protein densities in the bulk of the compartment along one of the bottom area axis (y-axis) and along the walls (z-axis) are shown in figure 8(A). We see that spatial oscillations in the bulk only take place along the z-axis. On the y-axis kymograph (figure 8(A) top kymograph) the change of the density only takes place along the temporal axis. This corresponds to the change of total numbers of free and bound MinD proteins during the oscillation and no spatial oscillation takes place along the y-axis. On the z-axis kymograph (figure 8(A) bottom kymograph) we observe a spatial pole-to-pole oscillation. In figure 8(B) we show the densities of the free MinD proteins in the bulk of the compartment during the large pole-to-pole oscillation for $z=4\,\mu {\rm{m}}$. We see that the bulk proteins oscillate only along one axis (here the y-axis) and no spatial oscillations occur along the walls (z-axis).

Zoom In Zoom Out Reset image size

From our study of geometry B, we conclude that in non-flat membrane geometries the oscillations of the Min-proteins do not only take place along the membrane, but are rather defined by the canalized transfer of the proteins through the bulk. This shows the importance of three-dimensional simulations for non-flat membrane geometries in order to determine the self-organized oscillation modes.

We have also monitored the amount of membrane-bound MinD and MinE proteins and compared it to the total amount of proteins in the system. Interestingly, we observed that these two quantities have a linear dependence ${N}_{\mathrm{bound}}\propto {N}_{\mathrm{total}}$ as shown in figure 9(A). Surprisingly, this relation seems to hold for all geometries and dimensions, as indicated in figure 9(A) by the different symbols. Since in our simulations we keep density constant, the total amount of proteins scales linearly with the compartment volume (${N}_{\mathrm{total}}\propto V$). Overall, we conclude the following relation for the mean bound Min-protein density

$$\begin{eqnarray}&&{\bar{\rho }}_{\mathrm{bound}}=\displaystyle \frac{{N}_{\mathrm{bound}}}{A}\propto \displaystyle \frac{{N}_{\mathrm{total}}}{A}\propto \displaystyle \frac{V}{A},\end{eqnarray} \tag{ 7 }$$

where A denotes the total reactive membrane area. Thus the mean bound Min-protein density increases linearly with the volume-to-area ratio, as verified by figure 9(B). Strikingly, we again observe a dramatic difference between MinD and MinE. Although this scaling is the same for both, MinE is much closer to being completely bound (gray lines), indicating that once a stable oscillation emerges, MinE is almost completely depleted from the bulk.

Zoom In Zoom Out Reset image size

We next assessed the physiological relevance of these observations. In order to investigate how the volume-to-area ratio V/A changes during cell growth, we approximate the shape of an E. coli bacterium by a spherocylinder of length L and width W and evaluate V/A analytically. Since E. coli bacteria mainly grow in length and stay constant in width [59], we show the evolution of V/A as function of the cell length L for various fixed cell widths in figure 10(A). One clearly sees that V/A remains nearly constant as a function of the cell length L. Recalling that our previous observation stated that the mean membrane-bound Min-protein densities ${\bar{\rho }}_{\mathrm{bound}}$ scale linearly with V/A, this would suggest that living bacteria grow in a fashion that keeps both V/A and thus consequently ${\bar{\rho }}_{\mathrm{bound}}$ constant, which might be advantageous for the stability and robustness of the Min-oscillations and related processes, such as formation of the FtsZ-ring. A similar qualitative behavior is observed if we instead of a spherocylinder would assume an ellipsoidal shape as is shown in figure 10(B).

Zoom In Zoom Out Reset image size

Above we have seen that multistability is a recurrent phenomenon in the Min-system, both in geometries A and B. We next turn to a systematic investigation of the stochastic switching between two different oscillation modes. An oscillation mode transition of this type occurs frequently in a $8\,\mu {\rm{m}}\times 2\,\mu {\rm{m}}\times 0.5\,\mu {\rm{m}}$ compartment using geometry B. In this geometry the Min-system gives rise to both pole-to-pole and striped oscillations along the same (major) axis (kymograph in figure 11(A)). We measure the lifetimes (oscillation duration before a mode switch occurs) of the two modes during a 50 000 s long simulation trajectory. The resulting histograms of the lifetimes are shown in figure 12. We can approximate the mode switching as a Poisson process by assuming that the switching probability p(t) obeys $p(t)\propto \exp (-{kt})$. Fitting an exponential to the switching times histograms we obtain the switching rates of these processes. Here we neglect the measurements of short lifetimes below ${\tau }_{c}=100\,{\rm{s}}$ since our oscillation mode detection algorithm can miss a mode transition if its lifetime is shorter. The switching rate for the pole-to-pole oscillation is ${k}_{{\rm{p}}}=0.004\,22\,{{\rm{s}}}^{-1}$ and for the striped oscillation we find ${k}_{{\rm{s}}}=0.004\,06\,{{\rm{s}}}^{-1}$, thus the two modes seem to be equally frequent in this case.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We have also performed the same analysis with identical compartment geometry for the two other parameter sets (set B and C) as presented in table 1. Parameter set B gives rise to stable longitudinal pole-to-pole oscillations along the major axis (kymograph in figure 11(B)) in agreement with the results from [30], while parameter set C shows a qualitatively similar switching behavior as parameter set A (data not shown) with frequent mode transitions.

To analyze the effect of the Min-protein concentration on the oscillation mode switching, we have performed the same simulation with increased Min-protein particle densities using again parameter set A. The resulting fractions of the two oscillation modes during the 50 000 s simulation trajectory are shown in figure 13, and their corresponding switching rates k_p and k_s are given in table 2. We note that the fraction of striped oscillations increases monotonously with increasing particle density. In agreement with this, the rate k_s decreases with increasing particle density. In contrast, the rate k_p does not show a systematic change and seems to fluctuate strongly. The shift to the striped oscillations suggests that the Min-oscillations can be used not only to sense geometry, but also to sense concentrations.

Zoom In Zoom Out Reset image size

Wu et al have experimentally measured the fractions of different oscillation modes of Min-proteins in rectangular cell geometries of E. coli bacteria of various sizes with constant height [30]. In order to address these observations, we now turn to geometry C, which is a rectangular and fully membrane-covered compartment as sketched in figure 1(B). Throughout this section we keep the compartment height fixed at $z=1\,\mu {\rm{m}}$, while we vary the length (major axis) between 2 and $10\,\mu {\rm{m}}$ and the width (minor axis) of the compartment between 1 and $5\,\mu {\rm{m}}$ in steps of $1\,\mu {\rm{m}}$, respectively. To determine the oscillation mode fractions we calculate an ensemble average of the oscillation modes after 500 s simulation time. For each compartment size a sample of 20 independent simulations is used.

In the following analysis we only focus on the three main oscillation modes, which here are longitudinal pole-to-pole and striped oscillation (along the major axis) and transverse pole-to-pole oscillation (along the minor axis). Other oscillation modes as they have been reported in [30] were also observed using our framework, however, due to their low probability, a much higher sample size would be necessary to analyze their occurrence frequency with sufficient statistics. The results for the three oscillation modes are shown in figures 14(A)–(C). In general, each of the three oscillation modes dominates in one region of the phase diagram, but the transitions are fuzzy and therefore bistabilities occur. For compartments with length below $7\,\mu {\rm{m}}$ and width below $4\,\mu {\rm{m}}$, only pole-to-pole oscillations are observed. Most of those oscillations occur along the major axis of the system. The transverse pole-to-pole oscillations along the y-axis emerge most frequently in compartments with quadratic bottom area. Increasing the width further increases the fractions of transverse pole-to-pole oscillations at the expense of longitudinal pole-to-pole oscillations. Increasing the length for a fixed width shows a sharp transition from longitudinal pole-to-pole to longitudinal striped patterns at around 6 μm. For compartments with length larger than $7\,\mu {\rm{m}}$, the longitudinal pole-to-pole oscillations vanish almost entirely. In the region of both large long sides (around 7–10 μm in length) and large short sides (around 4–5 μm in width), the oscillation mode fractions are rather equally shared between longitudinal striped and transverse pole-to-pole oscillations, which is also in line with the bistability that we observed between these two patterns in the $10\,\mu {\rm{m}}\times 5\,\mu {\rm{m}}\times 0.5\,\mu {\rm{m}}$ compartment of geometry A as shown in section 3.1. Figure 14(D) summarizes these findings in a phase diagram that considers only the dominating mode (expect in the regions of clear bistability). Overall we find the determined oscillation mode fractions based on parameter set A and as shown in figure 14 to be in excellent agreement with the experimental results as reported by Wu et al [30].

Zoom In Zoom Out Reset image size