Cross-diffusion induced patterns for a single-step enzymatic reaction

Several different enzymes display an apparent diffusion coefficient that increases with the concentration of their substrate. Moreover, their motion becomes directed in substrate gradients. Currently, there are several competing models for these transport dynamics. Here, we use mathematical modeling and numerical simulations to analyze whether the enzymatic reactions can generate a significant feedback from enzyme transport onto the substrate profile. We find that this feedback can generate spontaneous spatial patterns in the enzyme distribution, with just a single-step catalytic reaction. However, patterns are formed only for a subclass of transport models. For such models, nonspecific repulsive interactions between the enzyme and the substrate, or attractive interactions between the enzyme and the product, cause the enzyme to accumulate in regions of low substrate concentration. Reactions then amplify local substrate and product fluctuations, causing enzymes to further accumulate where substrate is low. Experimental analysis of this pattern formation process could discriminate between different transport models. While current mathematical models aim at understanding how enzymes move in fixed substrate gradients, the feedback of the enzymatic reaction on the substrate distribution has not been addressed yet. Here, the authors show that, in the presence of cross-diffusion due to repulsive nonspecific interactions, such feedback can lead to spontaneous spatial pattern formation

E xperiments performed during the last decade found that at least eight different enzymes display a higher diffusion coefficient when the concentration of the corresponding substrate in solution is increased [1][2][3][4][5][6][7] . These increases are in the range of 24-80% relative to the diffusion coefficients without substrate. Most of these experiments relied on fluorescence correlation spectroscopy (FCS) measurements. Although artifacts introduced by this technique have been pointed out 8,9 , recent findings using other techniques 10,11 validated the phenomenon, which is often referred to as "enhanced diffusion". The underlying mechanism is still under debate 12 . Some experiments suggest that catalysis plays a key role 3,5,13,14 , while others indicate that enhanced diffusion persists when the substrate is replaced by an inhibitor 1,4,15 . According to the latest experiments performed with the enzyme urease 11,13 , it appears that both the binding and the catalysis step of the reaction scheme contribute to enhanced diffusion.
A related question is how enzymes behave in the presence of substrate gradients. The answer to this question appears to be complex. Some experiments suggest that enzymes drift downstream gradients of substrates, performing "antichemotaxis" 5,13 . Others suggest that enzymes move upstream gradients of substrates, performing "chemotaxis" 6,7,15 . Antichemotaxis can be explained based on just the enhanced diffusion 5,16 (enzymes accumulate in regions with low substrate concentration where they have a lower diffusion coefficient), chemotaxis cannot be generated by enhanced diffusion alone. However, cross-diffusion is a possible cause for enzyme chemotaxis 6 . Cross-diffusion describes the response of the enzyme to forces generated by gradients of substrate. Mathematically, it corresponds to an offdiagonal element in the diffusion matrix describing the combined motion of the enzyme and the substrate. It has been suggested that cross-diffusion can be due to specific interactions (ligand binding) between the enzyme and the substrate 6,15,17 or due to nonspecific interactions (e.g., steric, electrostatic, and van der Waals) 18 . Specific interactions only lead to chemotaxis, while nonspecific interactions can cause the enzymes to move both upor downstream the substrate gradient, depending on whether the interactions are attractive or repulsive, respectively. The model including nonspecific interactions 18 can be considered as a mathematical generalization of other existing models 5,6,17 .
While the existing models and experiments study how enzymes move in preimposed substrate gradients, they do not consider the feedback from the enzymatic reaction onto the substrate distribution. Here, we analyze the effects of this feedback starting from the most general transport model 18 . We show that spatial patterns can emerge in initially homogeneous systems if nonspecific interactions contribute to the accumulation of the enzyme in regions where concentration of substrate is low. Enzymes accumulating in these regions further deplete the substrate, causing the substrate gradient to become steeper, hence further increasing the accumulation of the enzyme. We obtain a set of conditions for the parameter range in which patterns form. We see that patterns arising from initial homogeneous concentrations can emerge only for nonspecific interactions driving the enzyme away from the substrate, hence only for the model proposed by Agudo-Canalejo et al. 18 , but not for the models proposed by Zhao et al. 6 , Mohajerani et al. 15 , and Jee et al. 5 , suggesting that the analysis of spontaneous pattern formation experiments can be used to discriminate between the different proposed models. Our findings imply that patterns can arise for a single-step enzymatic reaction even in the absence of autocatalytic activity or allosteric regulation. This is surprising given that the formation of conventional Turing patterns 19 with such simple reaction schemes requires at least a reaction network of three states, where forward and backward reactions are catalyzed by two different enzymes 20 .

Results
Model scenario. Our starting model system is depicted in Fig. 1. We consider a single-step enzymatic reaction in a narrow reaction chamber of length L connected via a permeable membrane to a large substrate reservoir. We assume only substrate and product molecules can diffuse through the permeable membrane, while enzymes are confined to the reaction chamber (since enzymes are typically larger than their substrates and products). We indicate with γ s and γ p the permeation rates of substrate and product, respectively, i.e., the rate at which substrate and product molecules diffuse through the membrane. Such rates are obtained by dividing the substrate and product permeabilities by the reaction chamber thickness. The reservoir has a fixed concentration of substrate s R and no products.
For the reaction occurring in the bulk of this effectively 1D system we assume a Michaelis-Menten scheme (Fig. 1, inset). The substrate S binds to the enzyme E with rate k on , forming a complex, and it can unbind with rate k off . The catalytic step of the reaction has rate constant k cat and catalysis is irreversible. This leads to a turnover rate per enzyme k cat FðsÞ :¼ k cat sðK M þ sÞ À1 , where s is the substrate concentration and K M ¼ ðk off þ k cat Þk À1 on . Following Agudo-Canalejo et al. 18 we assume that the diffusion coefficient of the enzyme depends on whether the enzyme is free, D f , or it is in its complexed form, D c . This change in diffusion coefficient can be due to changes in either the hydrodynamic radius or the conformational fluctuations of the enzyme upon substrate binding 4,21,22 . Typically, enzymes become more tightly folded upon substrate binding, i.e., D c > D f , consistent with the Fig. 1 Illustration of the system considered in this work. A single-step enzymatic reaction takes place in a narrow reaction chamber, which is coupled to a reservoir through a permeable membrane. The membrane allows for the exchange of only substrates and products with permeation rates γ s , γ p , respectively, but no exchange for the enzymes. Inset: Reaction scheme. The substrate S is converted into a product P by an enzyme E. The reaction follows a Michaelis-Menten scheme where substrate can bind to the free enzyme with rate k on , forming a complex, and unbind with rate k off . The catalytic step of the reaction has a rate k cat . We assume that the enzyme has diffusion constant D f or D c depending on whether the enzyme is free or in its complexed form. Nonspecific pairwise interactions between enzyme and substrate ϕ fs , ϕ cs can depend on the enzyme form. How nonspecific interactions between the enzyme and the product affect the results is discussed further below in this article. experimentally observed trend. We also assume that short-range nonspecific interactions can either cause the enzyme to move towards or away from substrate, effectively generating a phoretic drift velocity that can also be interpreted as cross-diffusion 18 . We denote these interactions via pairwise potentials ϕ fs , ϕ cs depending on whether the enzyme is free or in the complexed form, respectively, as depicted in Fig. 1(inset).
Under the assumptions that (a) the enzyme is very dilute and (b) the system is locally in chemical equilibrium (i.e., the timescales of diffusion and cross-diffusion are slower than the chemical reactions), one can derive an effective transport equation for the enzymes 18 , within the quasi-1D reaction chamber ( Fig. 1) oriented along the x-axis. Here, e(x, t) denotes the local enzyme concentration (regardless of free or complexed) and s(x, t) the substrate concentration. The effective diffusion coefficient of the enzyme, is a function of substrate concentration, interpolating between the diffusion coefficient of the free enzyme and the complex, with F(s) as defined above. For s ≪ K M , D e (s)~D f , whereas D e (s)~D c for s ≫ K M . Note that the disassembly of enzyme oligomers into monomers can also contribute to enhanced diffusion 9,10,23 . Equation (1) would then describe the motion of the enzyme irrespective of its oligomeric state. The cross-diffusion term D s xd ðs; eÞ of Eq. (1) describes how enzymes respond to gradients of substrate due to the short-range nonspecific and hydrodynamic interactions, Here we define what we will refer to as 'cross-diffusion strength' C s c=f :¼ N A k B Tλ s 2 c=f η À1 , where η is the viscosity of the fluid, k B the Boltzmann constant, N A the Avogadro number, T the temperature and λ s c=f the Derjaguin length 18,24 . The Derjaguin length is a parameter capturing the effective short-range interaction between the complex/free enzyme and the substrate. It is typically a few angstroms 25,26 and can at most be as large as the Debye length (screening length), which in typical buffer conditions is about 1 nm 12 . It is expressed via the integral λ s 2 c=f ¼ R 1 0 dhhðe Àϕ cs=fs ðhÞðk B TÞ À1 À 1Þ. λ s 2 c=f is positive (negative) when the interaction is attractive (repulsive) 18 . The derivation of λ s 2 c=f is similar to that of the second Virial coefficient for a real gas 27 , but it also includes hydrodynamic effects and is computed by assuming that the size of the enzyme is much larger than the interaction length. The sign of λ s 2 c=f determines the sign of C s c=f and it has opposite sign as compared to D s xd ðs; eÞ. For attractive interactions D s xd ðs; eÞ<0, i.e., the enzyme drifts towards higher concentrations of substrate (chemotaxis). The enzyme performs antichemotaxis for repulsive interactions, for which D s xd ðs; eÞ>0. Note that the effect of nonspecific interactions can also be written as a phoretic drift 18 by swapping the ∂ x s in Eq. (1) with e in Eq. (3), with a drift velocity directly proportional to the substrate gradient v ph ðs; The enzyme dynamics as given by Eq. (1) is valid if short-range nonspecific interactions between the product and the enzyme are neglected. We first neglect such interactions, but will discuss their effects further below. This allows us to derive simple analytical conditions for the determination of the relevant parameter ranges for the formation of spontaneous patterns in the enzyme, substrate and product profiles. This approximation is valid for γ p ≫ γ s . In such a regime products very quickly permeate out of the reaction chamber, impeding the formation of large product gradients and causing the cross-diffusion induced by the product to be negligible.
Note that the values of γ p and γ s can vary by several order of magnitudes. For a reaction chamber that is 1-10 μm thick, the permeation rates are in the range of γ s,p ≈ 10 −9 -10 4 s −1 28 . In the regime where the enzyme is dilute and Eq. (1) is valid, the crossdiffusion that enzyme molecules would induce on substrate and product molecules can be neglected. Moreover the cross-diffusions related to the substrate-product interaction can also be neglected as product and substrate molecules are typically much smaller than the enzyme and they do not affect each other as much as they affect the enzyme motion. Hence, the reaction chamber of Fig. 1 is described by the coupled reaction-transport equations where p(x, t) denotes the product concentration. The dynamics of the system of Eq. (4) is determined by the interplay of the different physical processes described above, each associated with a different timescale. The diffusion processes have timescales τ diff: ¼ L 2 D À1 e;s;p . The cross-diffusion process has a timescale τ xd ¼ L 2 e h ðD s xd s R Þ À1 and it represents the time it takes an enzyme to explore a length L, driven by a constant gradient s R L −1 , where s R is the maximal substrate concentration allowed in the system. Equivalently, we can say that the nonspecific interactions cause the enzyme to drift with velocity v ph (s, ∂ x s), with associated timescale τ ph: ¼ Lv À1 ph ¼ τ xd . The reaction has a timescale τ react: ¼ k À1 cat and the permeations have timescales τ perm: ¼ γ À1 s;p . Having the reaction chamber coupled to a reservoir avoids product accumulation and substrate depletion, generating a nonzero homogeneous steady-state, with concentrations e h , s h , and p h for enzyme, substrate, and product, respectively. For simplicity, we express e h and p h as functions of s h , Since F(s h ) is a monotonic increasing function of s h , it is possible to write s h and p h in terms of e h , which can be directly tuned in experiments via the total enzyme concentration [see Supplementary Eqs. (4) and (5)  Instability driven by substrate induced cross-diffusion. The homogeneous solution as given by Eqs. (5), (6) is stable for any positive values of the parameters for a well-mixed system, i.e., a system with no diffusion and no cross-diffusion (see Supplementary Note 1.2). Figure 2 shows the results of two simulations of the full system, Eq. (4), with periodic boundary conditions and parameters as given in the Supplementary Table 1. In Fig. 2a we see that the homogeneous solution is unstable and patterns form, for a value of λ s 2 f ¼ λ s 2 c ¼ À1Å 2 . A steady-state is reached after about 10 3 s, corresponding to the enzyme diffusive timescale that for the chosen parameters is the slowest timescale τ diff : ¼ L 2 D À1 f ¼ 10 3 s. In Fig. 2b we see that for λ s 2 homogeneous solution is instead stable. Hence, depending on the parameters, the system will spontaneously form patterns.
To characterize the instability of the homogeneous steady-state solution, v h = (e h , s h , p h ), we linearize Eq. (4), v → v h + δv, and make the exponential ansatz δv = v 0 e σt e iqx , with q the spatial frequency of the linear perturbation and σ the perturbation growth rate. If σ > 0, the perturbation grows with a timescale as given by σ −1 and patterns form (see Supplementary Note 1.5 for the determination of σ −1 for the patterns of Fig. 2a). Note that since the system of Eq. (4) is isotropic, the results below will also be valid for systems of higher physical dimensions. A positive σ can be found (see Supplementary Note 1.3) provided that with D 0 e ðs h Þ ¼ ∂D e ðsÞ=∂sj s h and F 0 ðs h Þ ¼ ∂FðsÞ=∂sj s h . Using the functional form of D e [Eq. (2)], the term given in square brackets equals to ÀD f F 0 ðD e FÞ À1 , which is negative. Hence for the homogeneous solution to be unstable we must have D s xd ðs h ; e h Þ > 0, i.e., enzyme molecules drift downstream gradients of substrate, performing antichemotaxis. We can have a positive crossdiffusion if the nonspecific interactions between enzyme and substrate molecules are repulsive, i.e. for negative cross-diffusion strengths C s c=f < 0. However, how repulsive do these interactions need to be? By using the definitions of D e (s), F(s) and D s xd ðs; eÞ is possible to show [see Supplementary Eqs. (20)(21)(22)(23)(24)(25) in the Supplementary Note 1.3] that Equation (7) is fulfilled if with β = γ s s R . This result holds in the strong depletion regime (s h ≪ s R ), where the effect of the reaction on the substrate dominates over the outflow to the reservoir. In this regime, the expressions are simpler and it is easier to pinpoint the driving mechanism behind the instability observed in Fig. 2. We refer the reader to the Supplementary Note 1.3 for the full analysis. The argument of the square root in the inequality Eq. (8) is positive if Eq. (9) is fulfilled. In the Supplementary Note 1.4 we show that the instability is a Type II instability 29 , meaning that σ = 0 at q = 0, see Supplementary Fig. 2. This is natural as the total amount of enzymes in our system is conserved and homogeneously increasing or decreasing perturbations, i.e., perturbations at q = 0, would correspond to changes in the total enzyme amount. By looking at the inequalities Eq. (9), (8), we can see that negative cross-diffusion strengths C s c=f < 0 are needed to have instabilities. Given that C s c=f :¼ N A k B Tλ s 2 c=f η À1 , this means that repulsive interactions (λ s 2 c=f < 0) are necessary. Diffusion tends to homogenize the concentration profiles and contributes with a positive term to the left hand side of inequality Eq. (9). In the case where the enzyme-substrate interaction is attractive, i.e., C s f > 0, also the second term on the left hand side of inequality Eq. (9) is positive. Hence the interaction between the enzyme and the substrate would need to change sign upon substrate binding to have an instability (C s c < 0). This could happen if, upon substrate binding, the electrostatic surface charge distribution of the enzyme changes significantly, for instance due to a conformational transition. There can be cases in which the change in nonspecific interactions is less abrupt and both C s f < 0 and C s c < 0. Even in the case for which C s c ¼ C s f ¼ C s , i.e., there is no change in interaction upon substrate binding, it is possible to have an unstable homogeneous solution for where we rewrote relation Eq. (9) with the use of the Stokes-Einstein relation, k B T = 6πηD f R f , and the definition of C s . It is interesting to note how relations Eq. (9), (10) depend on the substrate concentration. One could ask, given certain nonspecific interactions between the enzyme and the substrate, at which substrate concentration s Ã h should we begin to observe instabilities? From relations Eq. (9)-(10), we find that with λ 2 s < 0. The value of s Ã h given by Eq. (11) is similar to the value of s h above which the phoresis generated by the nonspecific interactions v ph (s, ∂ x s) dominates over the drift induced by the enhanced diffusion v bi ≔ −∂ x D e (s) 18 . They differ only by a prefactor (D c − D f )/D f in the second term inside the square root of Eq. (11).
Having v ph > v bi corresponds to D s xd ðs; eÞ > eD 0 e ðsÞ. However, s Ã h is not only determined by this inequality, but also from the inequality Eq. (7), where reaction also plays a role. By considering biologically relevant ranges, such as R f~1 -10 nm 30,31 , K M~1 0 −2 -10 6 μM 32 , |λ s | being smaller than the Debye length |λ s | ≈ 10 −2 -10 Å, we This can happen only for enzymes for which K M = 0.1-1M, which is only a small fraction of enzymes 32 . For the vast majority of enzymes jC s jK M ð4D f Þ À1 ( 1 and s Ã h ) K M , meaning that we are in the saturated regime of the reaction, for which F(s) ≈ 1. In such a regime, we find that s Ã h % ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Fig. 3 we plot the phase diagram of the system of Eq. (4) where on the abscissa we have s h and on the ordinate we have λ 2 s . Positive feedback causing pattern formation. What is the physical mechanism underlying the instability given by the Eqs. (9), (8)? We illustrate the feedback mechanism generating the pattern in Fig. 4. We have already seen that to have instabilities the cross-diffusion D s xd ðs; eÞ must be positive, see Eq. (7). The enzymatic current induced by the cross-diffusion is given by J e xd ¼ ÀD s xd ðs; eÞ∂ x s, which for a negative slope of substrate concentration generates a positive current for the enzyme, i.e., the enzyme moves away from a high substrate concentration. The current generated by the enhanced diffusion J e D ¼ À∂ x ½D e ðsÞe also consists of a motion of the enzyme away from high substrate concentrations. In Fig. 2a we can see how in a regime where patterns form, more product is generated in locations where the substrate concentration is low and the enzyme concentration is high. Repulsive nonspecific interactions and enhanced diffusion cause the enzyme to accumulate in such regions (Fig. 4 from (a) to (b)). This accumulation then generates a higher reaction flux in these regions. Having a stronger reaction flux where substrate is already low, as compared with regions where substrate is abundant, causes substrate gradients to become steeper [ Fig. 4 from (b) to (c)]. A steeper substrate gradient, in turn, causes both J e xd and J e D to increase, hence generating a further accumulation of the enzyme in substrate depleted regions [ Fig. 4 from (c) to (d)]. This positive feedback between reaction and enzyme accumulation, leads to the formation of patterns. The feedback cycle halts when the substrate concentration is too low and the reaction is balanced by the influx of substrate from the reservoir. Then the substrate gradient stops getting steeper and the system approaches a steady-state.
The role of enhanced diffusion. What happens without enhanced diffusion, with only nonspecific interactions? The diffusion function D e (s) in Eq. (4) becomes a constant, i.e., D e = D c , and the interval as given by Eq. (8) for the unstable wave vector q is slightly affected. The Eq. (9) characterizing the instability is unaffected. This result suggests that enzyme patterns for a singlestep reaction form only if driven by repulsive nonspecific interactions between the substrate and the enzyme. Why cannot patterns be generated simply by enhanced diffusion? By considering the strong depletion regime (s h ≪ s R ) and D s xd ðs; eÞ ¼ 0, Eq. (7) becomes Both D e (s) and F(s) have a Michaelis-Menten dependence on s for all the models proposed so far for the enzyme motion. They differ only in the prefactors and a non-zero offset for D e (s). Equation (12) is never fulfilled for such D e (s) and F(s) and consequently any initial perturbation of the concentrations is smoothed out by diffusion (See Supplementary Note 2.1). However, Eq. (12) will apply to any model with a spatially dependent diffusion coefficient that is coupled to another diffusing and reacting species that is being depleted. In case the diffusing species is being generated, the sign of the above inequality is flipped. Systems of this type are used to study bacterial motion. Interestingly synthetic bacterial populations show stripe patterns as they grow on semi solid agar plates 33 . These bacteria produce a certain signaling molecule through which they sense their own Fig. 3 Phase diagram of pattern formation. We plot the instability curve of the full model of Eq. (4) (light blue line) and the approximated curve valid in the strong depletion regime, where the reservoir concentration s R is much larger than the substrate homogeneous concentration s h , as given by relation Eq. (10) (green dashed line). As a proxy to measure the presence of patterns, we use the ratio of the maximum of the steady-state profile of the enzyme concentration over the minimumẽ max =ẽ min . Below the curve, the system of Eq. (4) is unstable and patterns arise, the ratioẽ max =ẽ min can be as high as 10 9 . Above the curve, the homogeneous solution is stable and the ratioẽ max =ẽ min ¼ 1. The white cross corresponds to the simulation shown in Fig. 2a. Fig. 4 Positive feedback mechanism behind the pattern forming process. a We start with a substrate gradient and a homogeneous enzyme profile. b The enhanced diffusive current J e D and the cross-diffusive one J e xd cause the enzyme to accumulate where substrate is low, moving from left to right as indicated by the orange arrow. c The substrate is depleted at a rate k cat eF (s) and the enzyme accumulated in the region of low substrate cause the substrate gradient to get steeper, as illustrated by the purple arrows. d This leads to a further increase of J e D and J e xd causing further accumulation of the enzyme, as indicated by the thicker orange arrow. This process repeats itself determining the patterns. concentration and regulate their mobility. If we neglect bacterial growth, we have a system of equations similar to the equations describing the enzyme and the substrate motion in Eq. (4). The enzyme would correspond to the bacteria and the substrate to the signaling molecule. For this system FðsÞ ¼ const: < 0 and D e (s) is a Hill function with D 0 e ðsÞ < 0. Hence patterns can form for D 0 e ðsÞD e ðsÞ À1 < q 2 D s β À1 , where now β < 0. This shows that such synthetic population of bacteria can form patterns even in the absence of growth. The Eq. (12) specifies the minimal ingredients for pattern formation for such systems. It implies that, if F(s) > 0, patterns form whenever D e (s) is more sensitive than F(s) to perturbations in the substrate concentration. Having a more sensitive D e (s) than F(s) causes a more sensitive response in the enzyme motion than the depletion due to the reaction. Consider a local increase in substrate concentration and that both D e (s) and F (s) are monotonically increasing functions of s. Having a more sensitive D e (s) than F(s) implies a higher increase in the current J e D due to enhanced diffusivity away from the substrate, as compared to the depletion of substrate due to the reaction. Molecules E then migrate to regions with low substrate and if they do so with a high enough rate they can cause substrate gradients to get steeper, as in step (c) of Fig. 4. However, for the enzyme model as given by Eq. (4), D e (s) and F(s) alone cannot generate instabilities (see Supplementary Note 2.1). The accumulation of the enzyme in low substrate region at a high enough rate can be guaranteed only via the cross-diffusive term D s xd ðs; eÞ. This is because D e (s) and F(s) are Michaelis-Menten functions with the same K M and therefore are equally sensitive to substrate perturbations. In the Supplementary Note 2.3 we show that if the Michaelis-Menten constant associated to D e is much larger than the one associated to F(s), Eq. (12) can be fulfilled, see also Supplementary Fig. 3. In this scenario, there can be a range of s h concentrations for which D e (s h ) is in the linear regime and sensitive to substrate perturbations, while F(s h ) is in the saturated regime and therefore insensitive to perturbations in s. In the Supplementary Note 2.2, we also show that if D e (s) is a Hill function and F(s) a Michaelis-Menten reaction, we can have spontaneous pattern formation for a Hill coefficient > 1.
Inclusion of enzyme-product interaction. Up to now we focused on the effects that short-range nonspecific interactions between the substrate and the enzyme have on the enzyme motion. We neglected the product-enzyme interaction. This approximation was helpful to derive the Eqs. (9)- (11) and check the validity of pattern formation for biologically relevant parameter ranges. Now we ask: How is pattern formation affected by the shortrange interaction between product and enzyme? The enzyme dynamics becomes affected by gradients of products (see Supplementary Note 3) and instead of Eq. (1) we get where D p xd ðs; eÞ is the cross-diffusion due to nonspecific shortrange and hydrodynamics interactions between the product and the enzyme. Similarly to the substrate induced crossdiffusion, D p xd ðs; eÞ ¼ À½C c=f η À1 is the cross-diffusive strength corresponding to the product-enzyme interaction for the complexed and free form, respectively. The enzyme drifts upstream (downstream) the product gradients for an attractive (repulsive) interaction, for which D p xd ðs; eÞ < 0 (D p xd ðs; eÞ > 0). Note that D p xd ðs; eÞ depends on the substrate concentration, because the fraction of enzymes in the complexed form is given by F(s).
It is possible to derive a condition characterizing the instability of the homogeneous solution in the presence of product induced cross-diffusion, similar to Eq. (7) (see Supplementary Note 3): As compared to Eq. (7), in Eq. (14) there is an extra term that involves the cross-diffusion D p xd ðs; eÞ. The prefactor of such term goes to zero in the limit of large γ p as expected. The product induced cross-diffusion contribute to the instability with opposite sign as compared to D s xd ðs; eÞ. We have seen in Fig. 4 that patterns form if the interaction between substrate and enzyme is repulsive (transition from (a) to (b)) and D s xd ðs; eÞ > 0. Then the reaction has the role of steepening preformed substrate gradients (transition from (b) to (c)), i.e., more products are formed where substrate is already low. This generates a gradient in the products that has a opposite sign as the substrate gradient, see Fig. 2a. Hence we can expect that patterns are favored by an attractive nonspecific interaction between the product and the enzyme, for which D p xd ðs; eÞ < 0. This favors the accumulation of enzymes in region of low substrate, further steepening the substrate gradient. In Fig. 5a we present the results of simulations where λ p f ¼ λ p c ¼ λ p . We vary Fig. 5 Pattern formation process as the short-range nonspecific interaction between product and enzyme is varied. The concentration profiles of substrate s(x, t), enzyme e(x, t), and product p(x, t) are plotted as a function of time and spatial coordinate. By varying the square of the Derjaguin length λ 2 p , we affect the nonspecific interaction between the enzyme and the product. a λ 2 p ¼ 1Å 2 : the interaction is attractive and the patterns are more pronounced as compared to the case where no interaction is present, which is shown in Fig. 2a. b λ 2 p ¼ λ 2 s ¼ À1Å 2 : the interaction is repulsive and equal in strength to the repulsive interaction between substrate and enzyme, which strength is given by λ 2 s ; we can see how the patterns are still present but less pronounced. c λ 2 p ¼ À10Å 2 : the interaction is so repulsive that the homogeneous solution is stable. ARTICLE COMMUNICATIONS PHYSICS | https://doi.org/10.1038/s42005-020-00427-w the value of λ p and fix the other parameters to be the same as the ones considered for Fig. 2a, as given in Supplementary Table 1. For λ 2 p ¼ 1Å 2 (Fig. 5a), i.e., for an attractive interaction between product and enzyme (D p xd ðs; eÞ < 0), we see that the steady-state pattern is more pronounced as compared to the case in which D p xd ¼ 0 (Fig. 2a). For λ 2 p ¼ λ 2 s ¼ À1Å 2 (Fig. 5b), we see how the pattern is still present, however, less pronounced. Note that we still get a pattern because we considered different permeation rates γ s = 1s −1 and γ p = 10s −1 . Finally for λ 2 p ¼ À10Å 2 (Fig. 5c) we see that the repulsion between the product and enzyme is so strong that the pattern disappears. In this case D p xd ðs; eÞ > 0 and the crossdiffusive term induced by the product prevails over the crossdiffusive term induced by the substrate. Therefore Eq. (14) is not fulfilled. Note that Eq. (14) can also be fulfilled if both crossdiffusions are negative, i.e., both substrate and product attract the enzyme. In this case, patterns form if the product attraction is stronger than the substrate one. As a result, enzyme would still antichemotax from the substrate and the intuitive mechanism leading to patterns illustrated in Fig. 4 still holds. Note that one could think of the intuitive mechanism leading to patterns in terms of the feedback between product induced cross-diffusion and product formation. Given an initial product gradient and a uniform enzyme profile, the enzyme would chemotax to the region of high product. Having more enzymes in this region would cause the formation of even more product. This would cause the steepening of the product gradient and the further accumulation of enzymes. This feedback mechanism is the 'time reversal' mechanism of the one illustrated in Fig. 4.

Discussion
We have seen that spontaneous patterns arising from homogeneous concentration profiles can form for a single-step catalytic reaction if cross-diffusive effects are present. Patterns form given sufficiently strong repulsive nonspecific interactions between the enzyme and the substrate, and/or sufficiently strong nonspecific attractive interactions between the enzyme and the product as indicated by Eq. (14). These interactions cause the enzyme to move away from regions of high substrate concentrations and to accumulate in regions of low substrate, performing antichemotaxis. The accumulated enzymes then deplete the substrate, steepening substrate, and product gradients. Steeper gradients further drive the accumulation of enzymes as illustrated in Fig. 4. This positive feedback cycle between enzyme accumulation and reaction is what generates the patterns.
The enzyme accumulation is driven by the antichemotaxis with respect to substrate gradients due to nonspecific interactions. The enzyme chemotaxis considered in some of the models 6,15,17 has a stabilizing effect. Enzymes accumulate in regions of high substrate concentrations. Then reactions flatten substrate gradients, breaking the feedback that leads to patterns. Hence for such systems spontaneous patterns cannot form. Even in the absence of short-range nonspecific interactions, antichemotaxis exists due to enhanced diffusivity 5,13 . Nevertheless, we have seen that enhanced diffusivity alone cannot generate spontaneous patterns for a simple enzymatic reaction. Hence, among all models proposed so far for the enzyme motion, only the model given by Eq. (1), first proposed by Agudo-Canalejo et al. 18 , can lead to spontaneous pattern formation. The analysis of the patterns generated can shed light onto the microscopic interactions between enzyme and substrate. In fact, from Fig. 3, we can see how the ratio of the maximum over the minimum of the steadystate enzyme profile is affected by the Derjaguin length jλ 2 s j and it does not depend on the substrate concentration s h .
The patterns observed in Fig. 2 are not generated via the common short-range activation and long-range inhibition mechanism 34 , as neither the enzyme nor the substrate have autocatalytic activity. It is also not a motility induced phase separation (MIPS) mechanism 35 , which relies on the slowing down of active particles in regions of high particle concentrations. Here a positive feedback mechanism between particle accumulation and reaction leads to the pattern formation. Note that for the system of Eqs. (4), patterns can form even if D e~Ds , because inequality (9) does not depend on D s , whereas for classical Turing patterns large differences in the diffusion coefficients of the different species are required 19 . Moreover, it is surprising to see that patterns can form for a single-step enzymatic reaction with no autocatalytic activity nor allosteric regulation. In fact, for a system where species have a constant diffusion coefficient and enzymatic reactions follow a simple Michaelis-Menten scheme, patterns form for a minimal network of three states, where forward and backward reactions are catalyzed by two different enzymes, respectively 20 . In our system patterns form for a single-step catalytic reaction because of cross-diffusion.
Our findings are consistent with recent studies analyzing the effects of cross-diffusion in pattern formation 36 . Moreover it has been shown that phase separation, formation of static or selfpropelled aggregates can be observed in mixtures of crossdiffusive species interacting via a fast diffusing chemical that can be produced or consumed 37 . The scenario considered here corresponds to the case of a single cross-diffusive species, i.e., the enzyme, that is able to consume the fast diffusing chemical, i.e., the substrate. Here we considered the full nonlinear forms of enhanced diffusion and cross-diffusion for the enzyme motion, whereas these other studies considered constant diffusion and constant cross-diffusion 36,37 . The nonlinear model permitted us to address the question why enhanced diffusion alone is not able to generate patterns for a simple enzymatic reaction. We found that enhanced diffusion would need to be more sensitive to perturbations in substrate concentrations than the reaction, see Eq. (12). Although this is not the case for enzymes, we believe that Eq. (12) can characterize the pattern formation of species presenting different enhanced diffusion functions, such as synthetic bacteria 33 . Moreover, the cross-diffusion terms characterizing the dynamics of E in Eq. (13) can correspond to chemotactic responses of organisms to gradients of S and P, which can be chemoattractants or chemorepellants 38,39 . Hence Eq. (14) can be relevant in the characterization of spontaneous pattern formation of species performing chemotaxis. As similarly shown in Fig. 4, patterns can form if a chemotacting organism is repelled from a chemorepellent and such chemorepellent is also depleted from the organism. Moreover, patterns can also form if the organism is attracted to a chemoattractant and the chemoattractant is produced by the organism. The latter scenario corresponds to the 'time reversed' mechanism of the one shown in Fig. 4 and is at the base of aggregation phenomena in chemotacting Amoebae 40 .

Methods
Simulation parameters. The numerical simulations of the system of Eqs. (4) have been carried out by using COMSOL Multiphysics v5.3. We used the parameters listed in Supplementary Table 1 for the simulations shown in Fig. 2. We used λ s 2 f ¼ λ s 2 c ¼ À1 Å 2 for the unstable homogeneous solution (Fig. 2a) and λ s 2 f ¼ λ s 2 c ¼ 1 Å 2 for the stable homogeneous solution (Fig. 2b). The initial homogeneous concentrations e h , s h , and p h were perturbed by adding white Gaussian noise with variance e h ⋅ 10 −2 , s h ⋅ 10 −2 , and p h ⋅ 10 −2 , respectively. We considered periodic boundary conditions and we used the "Time dependent" solver of COMSOL with a relative tolerance of 10 −6 and an element size of 0.01 μm, i.e., 10 4 lattice points.
For the COMSOL specific parameters, we used a relative tolerance of 10 −6 and a maximum element size of 0.02 μm, except when two peaks were observed as a final result of the simulation. In these cases, we repeated the simulations with a finer grid (element size of 0.01 μm) and again we observed a single peak for the concentrations as the final result of the simulations. All the other parameters were the same as of Supplementary Table 1 and the boundary conditions were periodic. For the simulations shown in Fig. 5 we modified the dynamics of the enzyme, following Eq. (13). We then considered the same simulation parameters as shown in the Supplementary Table 1, the same COMSOL specific parameters used for the simulations depicted in Fig. 2 and we considered the λ p values shown in Fig. 5.

Data availability
The data that support the findings of this study are available from the corresponding author upon reasonable request.

Code availability
The code that supports the findings of this study is available from the corresponding author upon reasonable request.