Boundary-Induced Pattern Formation from Temporal Oscillation: Spatial Map Analysis

Boundary-induced pattern formation from a spatially uniform state is investigated using one-dimensional reaction-diffusion equations. The temporal oscillation is successively transformed into a spatially periodic pattern, triggered by diffusion from the fixed boundary. We introduced a spatial map, whose temporal sequence, under selection criteria from multiple stationary solutions, can completely reproduce the emergent pattern, by replacing the time with space. The relationship of the pattern wavelength with the period of oscillation is also obtained. The generality of the pattern selection process and algorithm is discussed with possible relevance to biological morphogenesis.

Spatial patterns are ubiquitous in non-equilibrium systems, and the process of spontaneous pattern formation has been extensively studied as one of the main issues in non-linear dynamics. Topic areas in pattern formation include fluid, solid-state, optical, geophysical, chemical, and biological systems [1][2][3][4][5].
In particular, the reaction-diffusion system, which was pioneered by Turing [6], has been a focus of non-equilibrium studies. In his celebrated study, Turing showed how a spatial pattern or temporal rhythm is spontaneously formed from a spatially homogeneous and temporally stationary state, which is unstable against perturbations. Turing classified such pattern formation processes into six cases, one of which is known as the Turing pattern; stationary periodic pattern with the finite wavelength given by linear stability analysis.
While Turing originally proposed his theory as a model of morphogenesis, it has been applied to general pattern formation dynamics beyond developmental biology [7,8].
Spontaneous pattern formation from temporally dynamic states is also possible in a reaction-diffusion system. Destabilization of a spatially homogeneous and temporally oscillatory state by diffusive interaction often results in spatiotemporal dynamics such as waves, spirals, and turbulence [9]. Spatiotemporal dynamics by combination of Turing instability for pattern formation and Hopf bifurcation for temporal oscillation have been studied as Turing-Hopf bifurcation [10][11][12].
However, an important factor that has not been fully explored is the influence of boundary conditions on spatial patterns. In particular, pattern formation from temporally dynamic states may be crucially influenced by the introduction of a fixed boundary. The influence may be globally propagated to alter the pattern dynamics. (See also [13,14] for the relevance of boundary conditions in pattern dynamics.) In the present letter, we are concerned with such boundary-induced pattern formation, namely, how a fixed boundary condition destabilizes the stable, temporally periodic, spatially uniform state and transforms it into a temporally stationary, spatially periodic pattern.
The question of boundary-induced pattern formation is not only of theoretical interest but also of experimental interest in developmental biology. In the somitogenesis of vertebrate development, temporal oscillation in protein expression is fixed into a spatial pattern [15,16].
As development progresses, the intracellular concentration of the protein c-hairy1, which initially oscillates in the system, is fixed into a striped pattern along aligned cells. So far, such spatial pattern formation has been studied by introducing external inputs that move in space with time development (clock and wavefront) [17]. Diffusive interaction has been introduced by Meinhardt in addition to external input using a spatial gradient [18].
For these theoretical studies, the use of an external input is essential to destabilize the homogeneous state for pattern formation. Some recent experiments, however, suggested that an external input may not be essential, and intrinsic instability due to diffusion might result in somitogenesis [19,20]. If boundary-induced pattern formation is possible without imposing an external input throughout the space, it will provide a plausible mechanism for vertebrate somitogenesis.
Here, we demonstrate that temporal oscillation in a one-dimensional reaction-diffusion system is fixed into a stationary periodic pattern by introducing fixed boundary condition.
In contrast to the celebrated Turing pattern, the wavelength of the generated pattern cannot be obtained using linear stability analysis around the fixed point. Instead, to predict the selected pattern, we introduced a one-dimensional spatial map, whose attractor gives the one-dimensional pattern by replacing time with space. The generality of this oscillation fixation and the pattern selection mechanism will be discussed. Now, consider a one-dimensional reaction-diffusion system of two components X and Y .
We assume that the diffusion of X is much faster than that of Y , and the diffusion of the latter is neglected for simplicity (the formalism to be discussed is valid even without this approximation). Then, the equation is written as where f (X, Y ) and g(X, Y ) are the reaction functions for X and Y , D is the diffusion constant of X, and the attractor of the dynamical system without D is a limit cycle [21] .
We consider the case in which the spatially uniform, limit cycle state is stable against perturbations so that, under Neumann or periodic boundary conditions, this state is the attractor. However, under a fixed boundary condition, the variable X close to the boundary cannot oscillate, which may destabilize the oscillatory attractor. Indeed, we have found that a fixed periodic pattern is often generated in this case.
As a specific example, consider the following system, which describe the protein expression dynamics with two genes, where X inhibits the expression of Y and Y activates the expressions of both X and Y [22,23]. Without the diffusion term, this system has one unstable fixed point (X * , Y * ) = (1/2, 1/2) if β > 8, and the limit cycle is an attractor. Indeed, from initial conditions close to a spatially homogenous state, a uniform, limit-cycle state is reached if the boundary condition is Neumann or periodic.
In contrast, if a fixed boundary condition is adopted for at least one end, i.e., X(0) = X 0 , temporal oscillation is replaced by a fixed, spatially periodic pattern [24] (see Fig. 1A and 1B where ∂X ∂x | x=L = 0 is adopted at the other end). Note that a pattern of the same wavelength is organized independently of X 0 as long as X 0 is between [0,1](see Fig. 3). Here, oscillation ceased in the vicinity of x = 0, which works as a boundary for (X(x), Y (x)) for slightly larger x. The oscillation is successively fixed for larger x.
We analytically examined how the organized pattern is determined. First, we confirmed  Figure 1A). Then, we contrasted the marginal values k 1 and k 2 with the wavenumber of the emergent pattern by changing β. The observed wavelength agrees neither with k 1 nor with k 2 (see Supplemental Figure 1B). Even the dependence of the parameters on β does not agree. Indeed, this discrepancy is natural: pattern formation occurs from the homogeneous limit cycle, which is far from the unstable fixed solution.
Hence, the standard analysis for the Turing pattern does not work in this system. Now we need to find a procedure to determine the spatial pattern without utilizing linear stability analysis. We note that once the diffusion term is given, the stationary solution (X st (x), Y st (x)) can be obtained from the dynamical system of two variables X, Y , where standard nullcline analysis works. As shown in Fig. 2B, only the nullcline of X (not that of Y ) is horizontally shifted with the diffusion term. Accordingly, the fixed point Once the fixed point (X st (x), Y st (x)) is given, the gradient term is obtained. For a fixed pattern, the fixed point (X st (x), Y st (x)) and the gradient term are determined selfconsistently.
This self-consistent condition, given by ∂X(x, t)/∂t = 0, dY (x, t)/dt = 0, can be solved explicitly by adopting spatial discretization to compute the diffusion term. In general, this condition is given by where l is discretized space index with a as the discretized unit length, i.e., x = al and From these equations, we can derive the spatial map as follows: where g −1 X l is the inverse function of g(X, Y ) given that X is fixed to X l , which corresponds to the cross-points of the nullclines. With these equations, X l+1 is determined from X l , X l−1 , and Y l , while Y l+1 is determined as g −1 X l+1 (0) (for the use of a "spatial map" in the analysis of a spatial pattern, see also [25,26]). However, there can be multiple candidate solutions for g −1 X l+1 (0). In fact, not all of the above solutions are stable in terms of dynamical systems (1). Furthermore, this solution must be attracted from the uniform oscillatory state. These properties lead to the following selection principles.

Selection principles
1. Ignore candidates in which the resulting state (X l+1 , Y l+1 ) is not stable against small perturbations in X l+1 .
2. Ignore candidates in which the resulting value X l+2 is out of the basin of the original limit cycle.
3. Calculate X l+2 from all remaining candidates, and choose the solution that gives the smallest value of |X l+2 − 2X l+1 + X l |.
The first criterion is necessary for the stability of the novel fixed point. The second criterion is necessary for the pattern to emerge from the limit cycle. The third criterion implies the smallest net diffusion, which minimizes the distance from a uniform state. With the second and third criteria, the solution that is attracted from the uniform oscillatory state is selected.
An example of the pattern obtained from the spatial map with the above selection criteria is shown in Fig. 3, which agrees well with the numerically obtained pattern. The above procedure works independently of the parameter values, demonstrating its validity.
Note that the spatial periodic pattern is obtained as an attractor of the spatial map (4).
Hence, independently of the initial condition in the map, the same periodic pattern is selected as long as the initial condition in the spatial map belongs to the basin of the above attractor.
The initial condition in the spatial map corresponds to the value of fixed boundary. Hence, a pattern of the same wavelength is reached independently of the boundary value, while the phase of wave pattern is different, which is consistent with the numerical result (see Fig. 3).
In the state space, pattern selection using the spatial map can be described as follows.
When the diffusion term is small, the shift of the nullcline is small. The same stable fixed point, i.e., Y l on the same branch, continues to exist and is selected according to the above criteria. With iteration of the spatial map, the shift of the nullcline due to the diffusion term is accumulated and at some point, the fixed point of the same branch vanishes with bifurcation, or do not fulfill selection principles. The selected fixed point is replaced by that of the new branch, as the solution according to the above criteria. Repeating the above procedure, a periodic pattern is obtained.
This procedure suggests there can be a relationship between the period of oscillation T and the wavelength λ. First, by rescaling the spatial scale, the wavelength λ is proportional to √ D. (The period is scaled by τ if f and g are similarly scaled by 1/τ , but this is already set to unity in (1).) Now, following the vector field (f (X, Y ), g(X, Y )) in the state space ( Fig. 2A), (X(t), Y (t)) oscillates with the period T . On the other hand, when a stripe is formed, the "spatial orbit" (X(x), Y (x)) along space x orbits the state space similarly to the limit cycle orbit. Indeed, except for the vicinity of the Hopf bifurcation, we numerically confirmed that one stripe is generated in correspondence with one period of oscillation. With the parameter change, the (magnitude of) vector (f, g) in the state space changes, which changes the rate of change along the limit cycle as well as that of the "spatial orbit" in the same way. If these orbits do not have a small radius in the state space, both speeds are expected to change in proportion, implying that the period T and wavelength λ (or λ/ √ D) change in proportion. Of course, this is a rough estimate, but this proportionality approximately holds against changes in the parameter values β (and D), as shown in Fig.   4.
The fixation of the temporal to the spatially periodic pattern as well as the analysis using the spatial map can be generally applied to a one-dimensional reaction-diffusion equation that has a uniform limit-cycle attractor. As another illustration, we examined the so-called Brusselator [1]: In this case, a uniform periodic oscillation is replaced by a periodic spatial pattern by applying the fixed boundary condition X(0) = X 0 for a certain range of parameters, which agrees with the prediction using the corresponding spatial map (see Supplemental Figure   2).
The approximate proportionality between the period and wavelength is worse than the model (2) (see Fig. 4), possibly because the change in the flow (f (X, Y ), g(X, Y ) against the parameters is highly dependent on the state (X, Y ).
In the present letter, we have studied the formation of a periodic pattern from a uniform oscillatory state induced by the boundary conditions. The emergent pattern is predicted as an attractor of the spatial map under the selection principle. We have confirmed the generality of this pattern formation in reaction-diffusion systems. The periodic oscillation exists in a two-component reaction system with negative feedback. That is, in a system with an activator and an inhibitor, for a certain range of parameters, the spatially periodic pattern that exists as an attractor replaces the uniform oscillatory state, triggered by the fixed boundary condition, if the diffusion of the inhibitor is sufficiently large.
We also note that the pattern formation process as well as the nullcline analysis is still valid even for D Y = 0, although the spatial map explicitly includes the inverse function of Y , which would make the analysis more difficult.
Pattern formation generally depends on the boundary and initial conditions. In the spatial map, the former is represented by the initial condition in the map, while the latter is considered by the second and third criteria of the selection principle, which indicate that the initial condition is not far from the uniform oscillatory state. If the initial condition is far from uniformity, local inhomogeneity can grow. In the state-space representation, whether such growth occurs is determined according to the size of the spatial diffusion term, specifically whether the spatial diffusion term is large enough to go across the nullclines and to induce saddle-node bifurcation. In other words, the condition close to uniformity represents the absence of such inhomogeneity. On the other hand, for such initial conditions to allow for the local growth in inhomogeneity, the selection criteria are replaced so that the choice of the fixed-point solution is close to a given initial condition, leading to a requested switch to a different branch of fixed-points. Hence, the choice of the selection principle corresponding to the initial condition can predict the emergent pattern. Control of an emergent pattern is thus possible by manipulating the initial condition.
Here, one stripe is formed with one-to-one correspondence of one oscillation period in time. Near the Hopf bifurcation point, however, a complex pattern with one stripe is observed per few periods of oscillation. Indeed, the spatial map can have an attractor with such complex (or quasiperiodic) oscillation in the vicinity of the bifurcation and also by modifying the selection principle. Analysis of such complex patterns in terms of the spatial map (see also [26]) will be of interest in the future.
Experimental confirmation of the present pattern formation will be possible in reactiondiffusion systems. In particular, relevance of boundary condition to pattern selection should be of importance because the real experimental system is finite and often under a fixed boundary condition. By carefully examining the boundary effects, we can confirm the present pattern formation mechanism, where the period-wavelength relationship (Fig. 4) will be confirmed.
As mentioned in the introduction, spatial pattern formation based on temporal oscillation is often observed in biological morphogenesis [15,16] as well as in the numerical evolution of morphogenesis [27,28], where the relevance of cell-cell interaction has been recently discussed in addition to the external morphogen gradient [20,28]. Considering the simplicity in our mechanism, which only requires diffusion of the inhibitor and a fixed boundary, we expect that it could be adopted in biological development, which will be confirmed by examining cell-cell interaction and boundary effects.    X l X 0 =0.0 X 0 =1.5 X 0 =3.0 X 0 =0.0 (predicted) X 0 =1.5 (predicted) X 0 =3.0 (predicted) SF2 Stationary pattern (X st (l), Y st (l)) of the reaction-diffusion (5) and the predicted pattern of (X l , Y l ) obtained from spatial map (4). The predicted pattern (given by •) agrees well with the stationary pattern (lines). Three patterns with boundary values of X 0 = 0.0 (purple), X 0 = 1.5 (green), and X 0 = 3.0 (yellow) are considered.