Edinburgh Research Explorer Kinetic roughening in active interfaces

. The essential features of many interfaces driven out of equilibrium are described by the same equation—the Kardar-Parisi-Zhang (KPZ) equation. How do living interfaces, such as the cell membrane, ﬁt into this picture? In an endeavour to answer such a question, we proposed in [F. Cagnetta, M. R. Evans, D. Marenduzzo, PRL 120, 258001 (2018)] an idealised model for the membrane of a moving cell. Here we discuss how the addition of simple ingredients inspired by the dynamics of the membrane of moving cells a ﬀ ects common kinetic roughening theories such as the KPZ and Edwards-Wilkinson equations.


Introduction
The theory of kinetic roughening describes the properties of surfaces which might appear smooth or rough depending on the scale of observation [1,2].Although usually applied to inanimate surfaces, such as solid-liquid interfaces, the theory can also describe biological, animate interfaces and this is the subject of the present paper.The Eden model, for instance, which represents one of the earliest attempt towards a probabilistic formulation of cluster growth, was introduced as an oversimplified portrayal of an expanding bacterial colony [3].
However, while some biological interfaces are simply understood in terms of scaling concepts, the features of some others have proven far more challenging to fathom.For instance, the issue of whether the low-frequency fluctuations in the shape of red blood cells are thermal or not has been resolved only recently [4], after a forty year debate [5].It was originally believed, in fact, that the observed 1/k 4 spectrum of fluctuations could be described by the energy equipartition principle [5].However, Prost and Bruinsma [6] showed that active fluctuations, generated by ATP-consuming processes, contribute to the dynamics of biological interfaces as well as thermal fluctuations.This initiated the field of active interfaces [7][8][9]: whithin such framework, the characteristics of fluctuations spectra (such as the 1/k 4 behaviour in flickering) can be attributed to various membrane activities (such as that of ion channels).It was indeed shown, in [4], that membrane response and fluctuations violated the fluctuation-dissipation relation, indicating the non-equilibrium nature of fluctuations.All this indicates that an understanding of active process is crucial for any quantitative description of membrane dynamics.
In the present work we are concerned with active processes that are related to cell locomotion.The perspective we take [10,11], as in usual kinetic roughening theories, is focused on the scaling of the interface width w with the * e-mail: F.Cagnetta@ed.ac.uk system size [1].We present a brief review of such theories (Edwards-Wilkinson and Kardar-Parisi-Zhang equations) in Section 2. In Section 3 we introduce our modelling approach, inspired by the leading edge of eukaryotic cells crawling on two-dimensional substrates: the protruding force which sets the membrane in motion is exerted by the actin cytoskeleton below the membrane, but directed by membrane proteins which collect signals from the exterior of the cell.The first model we discuss consists indeed of a moving interface whose growth is effected by a number of diffusing particles, which we refer to as activators.The model is a linear theory for active interfaces which can be thought of as an extension of the Edwards-Wilkinson equation.In Section 4 we show that lattice simulations of growth with diffusing activators presents different scaling behaviour, which we ascribe to a non-linear effect analogous to that considered in the Kardar-Parisi-Zhang equation.Finally, in Section 5, we consider in detail the coupling between the interface shape and the activators distribution generated by active growth.

Theory of kinetic roughening
Consider an interface described by a time-dependent height function h(x, t) over some substrate x.The average of h(x, t) over x yields the mean distance from the substrate, while the variance gives the squared width w 2 (L, t) about the mean profile, function of the system linear size L. The width is related to the time-dependent structure factor (square modulus of the height Fourier modes h k (t)) by In a system of linear size L, the smallest wavenumber allowed is 2π/L: hence, the small-k scaling of the structure factor translates into the scaling of the width with L, a fact often referred to through the statement that the interface roughness depends on the scale of observation [1].In fact, EPJ Web of Conferences 230, 00001 (2020) https://doi.org/10.1051/epjconf/202023000001FisMat 2019 not only the width depends on system size, but also the time it takes for the width to develop, starting from a flat interface with w = 0.Both these facts are summarised in the Family-Vicsek scaling hypothesis [12], with roughness exponent α, dynamic exponent z and scaling function f (x) which is constant at large x and behaves as x α/z for small x.
The simplest model satisfying the Family-Vicsek scaling is the Edwards-Wilkinson (EW) equation [13], with η a gaussian, white, space-time noise with zero mean and unit variance.ν is the surface tension of the interface, D h quantifies the intensity of fluctuations.Upon identifying D h with the temperature k B T , the right hand side of (3) represents the competition between the smoothing effect of surface tension and random height fluctuations generated by thermal noise.By solving Eq. ( 3) for the modes of h(x, t) in one spatial dimension, with flat initial condition h(x, t = 0) = 0, a simple calculation (see e.g.[2]) shows that Eq. ( 2) is obeyed with α = 1/2 and z = 2. Additionally, z = 2 holds for all dimensions, while α = 0 in d ≥ 2, implying a smooth interface.When the interface is driven out of equilibrium, (3) is augmented with additional terms.On the one hand, a uniform driving force, in the form of a constant term λ added to the right-hand side of Eq. (3), can be removed with a galilean shift of the height h → h − λt, thus it would not alter the EW scaling.On the other hand, geometric considerations [14] imply that driving forces are generically directed along the local normal to the interface, so that the λ-term aquires a projection factor 1 + (∇h) 2 (cf.Fig. 2).Expanding the square root for ∇h small yields the celebrated Kardar-Parisi-Zhang (KPZ) equation [14] characterised by scaling exponents α = 1/2 and z = 3/2 in one spatial dimension.

Growth by diffusing activators: modified EW equation
On the front of a moving cell, or leading edge, the driving force is localised around the positions of specific membrane proteins [6] which catalyse membrane growth [8].
As they are the source of active growth, we call these proteins activators and denote their position with X i (t), i = 1, . . ., N. We begin by assuming that i) the plasma membrane can be described by an height function h(x, t); ii) the height fluctuations are described by (3) in the absence of activators; iii) the activators' motion within the interface is purely diffusive.We then write the following field equation for our active interface (and activators), with ξ i an independent Gaussian white noise with zero mean and variance 2D a for each i.We stress that the diffusion coefficient of the activators D a is not related to the coefficient of the interface noise η, as the activator noise regards fluctuations of the positions X i within the interface, while the interface noise refers to fluctuations of the d-dimensional interface in the d + 1 dimensional space.
Eq. (5a) can be solved via Fourier transform with and integrating from a flat interface at t = 0, we get Upon averaging w.r.t the noise distribution, the η k term vanishes, while the average of the complex exponential coincides with the characteristic function of the variable s 0 du ξ(u).As the latter variable is Gaussian with average zero and variance 2Ds, the characteristic function is e −2D a sk 2 .Thus, after inverting the Fourier transform, the average height profile reads at variance with the vanishing average height of the EW equation.Eq. ( 8) represents the convolution of the EW response function (1−e −νk 2 t )/νk 2 , with ν replaced by (ν−D a ), with the probability density function p a (x, t) of the position X(t) of the activator.Contributions due to different activators add up linearly.Eq. (7) gives us access to the structure factor too, which reads, The second contribution in the square brackets coincides with the characteristic function of the Gaussian variable s 1 s 2 du ξ(u), whose variance equals 2D|s 1 − s 2 |.By plugging in also the η noise correlations in the k-space, we get The generalisation to multiple activators is simple, at least for the modes with k 0 (the only modes relevant for the width, according to Eq. ( 1)).Specifically, it suffices to multiply the second term in Eq. ( 10) by the number of activators N, so that the overall density ρ 0 = N/L d appears.
As time is always multiplied by k 2 in Eq. ( 8), we conclude that z = 2.To compute α, we consider the t → ∞ limit, where all the exponential factors vanish and we are left with The small k behaviour of Eq. ( 11) is dominated by the 1/k 4 term, so that, in the L → ∞ limit, w 2 (L) ∼ L 4−d , i.e. α = (4 − d)/2 for d < 4, 0 otherwise.Thus, in one spatial dimension the roughness exponent (3/2) is much higher than that of the standard EW and KPZ classes (1/2).

Growth by diffusing activators: modified KPZ equation
We now study a stochastic lattice model with diffusing activators which should display similar scaling to the field equations discussed in the previous section.There are, however, two crucial differences.First, the difference in height between neigbouring lattice sites is fixed to 1 (solidon-solid condition) the maximum possible width is O(L), thus the maximum allowed value of α is 1.Second, when a continuum description is derived from the lattice model, non linear KPZ like terms may be generated.We shall discuss these issues below, but first we define the lattice model.
The model, introduced in [10], consists of a discrete interface, made of L unit-slope segments, and a collection of N activators.Both the interface and the activators live on the one-dimensional lattice, with periodic boundary conditions enforcing the ring topology.Being made of unitary slopes, the interface can be described with a set of height variables {h i } over the lattice points i = 1, . . ., L which obey the solid-on-solid condition |h i+1 − h i | = 1.Each activator is represented by a discrete random walk hopping between neighbouring lattice sites.The activator dynamics is thus specified by the hopping rate q, which also coincides with twice the diffusion coefficient 2D a .
The interface, in turn, evolves according to a singlestep dynamics.Due to the solid-on-solid condition, each site can be a peak (∧), a trough (∨) or a slope ( or ).Troughs can grow and become peaks at rate p + whereas peaks become troughs at rate p − , so that the solid-on-solid condition is preserved at all times.In order to account for the growth-stimulating action of the activators, we take the interface rates p ± to depend on the local number of activators n i , such that According to Eq. ( 12), each activator increases the growth rate of the interface by λ.As we assume no exclusion interaction among the activators, n i = 0, . . ., N.
In the lattice model, the width can be directly measured as the variance of the height variables h i over the lattice.We then track the time-dependent width of active interfaces of various sizes in Monte Carlo simulations of the update rule Eq. ( 12).The results are plotted in Fig. 1 Figure 1.Family-Vicsek scaling of the width w(L, t) for the lattice model of active interface with diffusing activators, for q = p = 1 and λ = 1, system size L as in the key.Averages are performed over 100 independent realisations of the stochastic dynamics.The best collapse of the curves is achieved with z = 3/2, as in the KPZ class, and α = 1, the maximal roughness of solid-on-solid models.The black dashed line is a guide to the eye for the ∼ t α/z law.
according to the Family-Vicsek scaling.Specifically, we found the best collapse to be achieved for α = 1 and The scaling is different from that of the linear model discussed in the previous section.First the roughness exponent α takes its maximal value 1, which is less than the value 3/2 for Eq.(5a), though still higher than the EW-and KPZ-class exponent 1/2.Second, the value 3/2 for the dynamic exponent suggests that a KPZ like nonlinearity is present.This can be traced back to the interface transitions occurring only at troughs and peaks in the solid-on-solid model, which implies a factor 1 − (∇h) 2 multiplying the driving force in Eq. (5a): This term has the typical form (albeit with a difference in sign) of the normal projection factor of the KPZ equation, illustrated in Fig. 2, panel A. In fact, the added non-linear term can be shown to be relevant in the renormalisation group sense [15], thus it is expected to change the scaling properties of the model.

Slope-coupling due to normal growth
In this section we introduce a further modification of Eq. (5a), which includes additional forces on the activators coming from geometric considerations [16].An explanation for the origin of such terms is sketched in Fig. 2. First, as in the KPZ equation [14], local growth forces are exerted along the normal to the interface (see Fig. 2, panel A).The infinitesimal displacement λδt due to a force λ over time δt should then be increased by a factor [1 + (∇h) 2 ] 1/2 1 + (∇h) 2 /2 (Fig. 2A).The same reasoning applies to a driving force which depends on the local density of activators along the interface, λρ(x, t), where ρ(x, t) = N i=1 δ(x − X i (t)).According to the same argument, the matter which constitutes the interface is displaced horizontally by δx = λδt(∇h)[1+(∇h) 2 ] −1/2 (Fig. 2, panel B), generating an effective coupling of the activators positions with the interface slope.We now write down equations for the density of activators ρ(x, t) and an interface h(x, t) driven by normal force λρ(x, t).By keeping only terms which are at most quadratic in the fields, the resulting equations read

A B
with ρ 0 = [0,L] d d d x ρ(x, t)/L d the global density of activators.We can represent the effective coupling of the activator positions with the interface slope via a minor modification of the lattice model.Instead of having the activators hop left or right on the lattice at the same rate q, we let the rate depend on the difference in height between arrival and departure site, i.e.
The parameter γ = q + −q − controls the rate of slope advection, while q = (q + + q − )/2 determines the mobility of the activators.The model defined by the rates of Eq. ( 12) and Eq. ( 15) can be shown to be described by field equations analogous to Eq. ( 14), for the special choice of parameters λ = 2γ [17].Here we present the scaling of the interface width and the variance of the density of activators measured from Monte Carlo simulations of the lattice model.

Width scaling
Let us begin with the width.The initial power-law increase, computed when starting from a flat initial condition, proceeds as in the case described in Section 4, Fig. 1.Let us recall that the activators of the simulations shown in Fig. 1 diffuse freely with no slope-coupling (γ = 0).By turning on the slope-coupling, after some time which decreases for increasing γ, the width decreases and begins oscillating, until it reaches saturation.This property is manifest in Fig. 3 The Family-Vicsek scaling (Eq.( 2)) is still obeyed for γ 0, with a different set of exponents.As shown in Fig. 3, bottom panel, the roughness exponent α goes back to the EW and KPZ value 1/2, typical of an interface with no correlations among the slopes in steady-state.However, the coupling of the activators position with the slopes results in an emergent ballistic behaviour, highlighted by the dynamic exponent z = 1.This is a peculiar property of the active interface model introduced in [10] and can be explained by solving the inviscid limit of the field equations ( 14) [17].

Density scaling
Figure 4. Scaling of the density variance for a passive (top) and an active (bottom) interface.Parameters are γ = 0.5, ρ 0 = 1, while λ and L are given in the key.In the passive case, the density variance obeys a scaling hypothesis analogous to the Family-Vicsek one, with dynamic exponent z = 2 and "roughness" α = 1/3.For the density, the power-law initial growth in time and the power-law dependence on the system size of the steady-state value are evidences of coarsening in the system.In the active case, the scaling hypothesis is still obeyed albeit with trivial exponents α = z = 0.
As the variance of the height profile provides a picture of the interface roughening dynamics, the variance of the activator density measures the activators coarsening.In the lattice model, where an occupation number n i (t) specifies the instantaneous number of activators at the i-th site, the density variance is defined as Under appropriate circumstances, δρ obeys a scaling hypothesis analogous to that of Family-Vicsek for the width, Eq. (2), i.e. δρ(L, t) = L χ g(t/L z ).As in the Family-Vicsek scaling, the dynamic exponent describes the dependence of the density relaxation time on the size of the system.The exponent χ, instead, by representing the scaling of the steady-state density variance with the systems size, indicates the extent of clustering in the system.If, for instance, macroscopic clustering takes place -i.e., a finite number of sites hosts a finite fraction of the available particlesthen a few n i 's scale as N = ρ 0 L d while all the others are close to zero, so that δρ 2 will scale as L d and χ = d/2.By contrast, for an homogeneous distribution of activators, δρ 2 does not depend on the system size, i.e. χ = 0.The scaling hypothesis of the density variance is obeyed in problems of non-interacting particles sliding down the slopes of a fluctuating interface [18][19][20][21].In the λ = 0 limit of our lattice model, corresponding to passive particles sliding on a EW interface, we find z = 2 and χ = 1/3 (cf.Fig. 4, top panel).When activity is turned on (bottom panel of Fig. 4), the density variance saturates at a finite value which is independent of the system size, indicating a homogeneous distribution of activators at large scales, i.e. χ = 0.The saturation time does not depend on the system size either, i.e. z = 0. Nevertheless, the initial power-law growth of the density variance, especially visible at smaller values of λ, indicates that an initially homogeneous distribution of activators coarsens in time as in the passive sliders case.However, the size of the aggregates remains finite rather than growing with the system size-a phenomenon interpreted as microphase separation in [10].

Conclusions
In this paper we have reviewed how interface growth equations and the resultant scaling are modified in the presence of diffusive particles which activate the growth.We first studied analytically a modified EW equation with the addition of activators and found novel scaling.However, this scaling is not observed in a simple simulation model due to the presence of KPZ-like non-linearities.We then presented numerical results for the dynamics of activators coupled to the interface shape, which generalises the problem of passive scalar advection by fluctuating interfaces.Our results indicate that a wealth of new dynamical behaviours are possible in the kinetics of active interfaces: scaling concepts might then be crucial in sorting them into different classes, as it is done for generic critical phenomena [22].

Figure 2 .
Figure 2. Pictorial representation of the effects of normal growth on the active interface.Both pictures show the interface profile at two infinitely close instants, with the cartesian and local normal-tangent reference frames displayed in the bottom left corner, respectively in black and gray.Panel A: The displacement δh due to a normal force λ over time δt reads δh = λδt[1 + (∇h) 2 ] 1/2 .Panel B: The activator (red disk) is displaced by δx = λδt(∇h)[1 + (∇h) 2 ] −1/2 due to the normal force acting on the interface.

Figure 3 .
Figure 3. Width of the active interface with slope-coupling γ.Top: comparison of the roughening profiles of interfaces with various γ, at λ = 1, system size L = 213 and density of activators ρ 0 = 1.The width of an EW interface is also shown for comparison as a teal solid line.The black dashed line is a guide to the eye representing the roughening law discussed in section 4, w ∼ t 2/3 .Bottom: Family-Vicsek scaling of the active interface with γ = 0.5, λ = 1, density ρ 0 = 1 and L as in the key.Here the oscillations are clearly visible.The best collapse is achieved for α = 1/2 and z = 1.The exponents z = 1 and α = 1/2 are typical of the whole λ, γ > 0 region of the parameter space.
, top panel, where the time-dependent widths of interfaces with different values of γ ∈ [0, 1] are compared.The width of an EW interface (γ = λ = 0) is also shown for comparison (solid line in the figure).