Generic transport mechanisms for molecular traffic in cellular protrusions

TL;DR

This paper introduces a lattice gas model to analyze molecular motor transport in cellular protrusions, revealing mechanisms for motor localization and task sharing that optimize energy use and are governed by strong correlations.

Contribution

It presents a novel analytical approach to understanding how active and diffusive transport mechanisms interact in cellular protrusions, highlighting the role of correlations in motor dynamics.

Findings

Identification of a mechanism for motor tip localization.

Demonstration of task sharing between active and diffusive transport.

Analytical calculation of reduced active currents due to correlations.

Abstract

Transport of molecular motors along protein filaments in a half-closed geometry is a common feature of biologically relevant processes in cellular protrusions. Using a lattice gas model we study how the interplay between active and diffusive transport and mass conservation leads to localised domain walls and tip localisation of the motors. We identify a mechanism for task sharing between the active motors (maintaining a gradient) and the diffusive motion (transport to the tip), which ensures that energy consumption is low and motor exchange mostly happens at the tip. These features are attributed to strong nearest-neighbour correlations that lead to a strong reduction of active currents, which we calculate analytically using an exact moment-identity, and might prove useful for the understanding of correlations and active transport also in more elaborate systems.

Figures

Equations

$$f_i={\rho }_i^{\mathrm{T}}+\omega \sum _{j=0}^i({\rho }_j^{\mathrm{T}}-{\rho }_j^{\mathrm{S}}).$$

$$\omega \sum _{j=0}^{L-1}\left({\rho }_j^{\mathrm{T}}-{\rho }_j^{\mathrm{S}}\right)=-\beta \left({\rho }_L^{\mathrm{T}}-⟨n_L^{\mathrm{T}}{n}_L^{\mathrm{S}}⟩\right),$$

$${\mathrm{Cov}}_i:=⟨n_i^{\mathrm{T}}{n}_{i+1}^{\mathrm{T}}⟩-{\rho }_i^{\mathrm{T}}{\rho }_{i+1}^{\mathrm{T}}=f_i-{\rho }_i^{\mathrm{T}}{\rho }_{i+1}^{\mathrm{T}}$$

$$\mathrm{Cov}(x)={\rho }^{\mathrm{T}}(x)\left[1-{\rho }^{\mathrm{T}}(x+a)\right]+\Omega {\int }_0^x\mathrm{d}y({\rho }^{\mathrm{T}}-{\rho }^{\mathrm{S}}).$$

$${\rho }^{\mathrm{T}}(x)\approx \frac{\mathrm{erf}\left((x-z)/W(z)\right)+\mathrm{erf}\left(z/W(z)\right)}{\mathrm{erf}\left((1-z)/W(z)\right)+\mathrm{erf}\left(z/W(z)\right)},$$

$$\frac{\mathrm{d}\text{Prob}\{n_i^{\mathrm{T}}\}}{\mathrm{d}t}$$

$$\frac{\mathrm{d}\text{Prob}\{n_i^{\mathrm{S}}\}}{\mathrm{d}t}$$

$$+\omega \left(\text{Prob}\{n_i^{\mathrm{T}},\overline{n_i^{\mathrm{S}}}\}-\text{Prob}\{n_i^{\mathrm{S}},\overline{n_i^{\mathrm{T}}}\}\right),$$

$${\partial }_t{\rho }_i^{\mathrm{T}}={\rho }_{i-1}^{\mathrm{T}}-f_{i-1}-{\rho }_i^{\mathrm{T}}+f_i+\omega \left({\rho }_i^{\mathrm{S}}-{\rho }_i^{\mathrm{T}}\right),$$

$${\partial }_t{\rho }_i^{\mathrm{S}}=ϵ\left({\rho }_{i+1}^{\mathrm{S}}+{\rho }_{i-1}^{\mathrm{S}}-2{\rho }_i^{\mathrm{S}}\right)+\omega \left({\rho }_i^{\mathrm{T}}-{\rho }_i^{\mathrm{S}}\right),$$

$${\partial }_t({\rho }_i^{\mathrm{T}}+{\rho }_i^{\mathrm{S}})={\rho }_{i-1}^{\mathrm{T}}-f_{i-1}-{\rho }_i^{\mathrm{T}}+f_i+ϵ\left({\rho }_{i+1}^{\mathrm{S}}+{\rho }_{i-1}^{\mathrm{S}}-2{\rho }_i^{\mathrm{S}}\right)$$

$${\partial }_t{\rho }_0^{\mathrm{T}}=-{\rho }_0^{\mathrm{T}}+f_0+\omega \left({\rho }_0^{\mathrm{S}}-{\rho }_0^{\mathrm{T}}\right)$$

$${\partial }_t{\rho }_0^{\mathrm{S}}=\alpha (1-{\rho }_0^{\mathrm{S}})-ϵ{\rho }_0^{\mathrm{S}}+ϵ\left({\rho }_1^{\mathrm{S}}-{\rho }_0^{\mathrm{S}}\right)+\omega \left({\rho }_0^{\mathrm{T}}-{\rho }_0^{\mathrm{S}}\right)$$

$${\partial }_t{\rho }_L^{\mathrm{T}}={\rho }_{L-1}^{\mathrm{T}}-f_{L-1}-\beta \left({\rho }_L^{\mathrm{T}}-⟨n_L^{\mathrm{T}}{n}_L^{\mathrm{S}}⟩\right)$$

$${\partial }_t{\rho }_L^{\mathrm{S}}=ϵ\left({\rho }_{L-1}^{\mathrm{S}}-{\rho }_L^{\mathrm{S}}\right)+\beta \left({\rho }_L^{\mathrm{T}}-⟨n_L^{\mathrm{T}}{n}_L^{\mathrm{S}}⟩\right)$$

$${\partial }_t{\rho }_L^{\mathrm{T}}={\rho }_{L-1}^{\mathrm{T}}-f_{L-1}+\omega \left({\rho }_L^{\mathrm{S}}-{\rho }_L^{\mathrm{T}}\right)$$

$${\partial }_t{\rho }_L^{\mathrm{S}}=ϵ\left({\rho }_{L-1}^{\mathrm{S}}-{\rho }_L^{\mathrm{S}}\right)+\omega \left({\rho }_L^{\mathrm{T}}-{\rho }_L^{\mathrm{S}}\right)$$

$${\partial }_t({\rho }_L^{\mathrm{T}}+{\rho }_L^{\mathrm{S}})={\rho }_{L-1}^{\mathrm{T}}-f_{L-1}+ϵ\left({\rho }_{L-1}^{\mathrm{S}}-{\rho }_L^{\mathrm{S}}\right)=J_{L-1}^{\mathrm{T}}+J_{L-1}^{\mathrm{S}},$$

$${\partial }_t({\rho }_0^{\mathrm{T}}+{\rho }_0^{\mathrm{S}})=-{\rho }_0^{\mathrm{T}}+f_0+\alpha (1-{\rho }_0^{\mathrm{S}})-ϵ{\rho }_0^{\mathrm{S}}+ϵ\left({\rho }_1^{\mathrm{S}}-{\rho }_0^{\mathrm{S}}\right)=-J_0^{\mathrm{T}}+\alpha (1-{\rho }_0^{\mathrm{S}})-ϵ{\rho }_0^{\mathrm{S}}-J_0^{\mathrm{S}}$$

$$J_i^{\mathrm{T}}={\rho }_i^{\mathrm{T}}-f_i$$

$$J_i^{\mathrm{S}}=ϵ\left({\rho }_i^{\mathrm{S}}-{\rho }_{i+1}^{\mathrm{S}}\right)$$

$${\partial }_t({\rho }_i^{\mathrm{T}}+{\rho }_i^{\mathrm{S}})=J_{i-1}^{\mathrm{T}}+J_{i-1}^{\mathrm{S}}-J_i^{\mathrm{T}}-J_i^{\mathrm{S}}.$$

$$J_{-1}^{\mathrm{S}}=\alpha (1-{\rho }_0^{\mathrm{S}})-ϵ{\rho }_0^{\mathrm{S}}$$

$${\partial }_t({\rho }_0^{\mathrm{T}}+{\rho }_0^{\mathrm{S}})=-J_0^{\mathrm{T}}+J_{-1}^{\mathrm{S}}-J_0^{\mathrm{S}}$$

$$0={\partial }_t\sum _{i=0}^L({\rho }_i^{\mathrm{T}}+{\rho }_i^{\mathrm{S}})=J_{-1}^{\mathrm{S}}$$

$${\rho }_0^{\mathrm{S}}=\frac{\alpha }{\alpha +ϵ}:=\sigma ,$$

$${\rho }_i^{\mu }\to {\rho }^{\mu }\left(x=x_i\right)\mathrm{where}{x}_i=ia.$$

$$J^{\mathrm{T}}(x)$$

$$J^{\mathrm{S}}(x)$$

$$0=-a{\partial }_x(J^{\mathrm{T}}+J^{\mathrm{S}})+\frac{1}{2}{a}^2{\partial }_x^2(J^{\mathrm{T}}+J^{\mathrm{S}})+\mathcal{O}(a^3)=$$

$$=-a{\partial }_x({\rho }^{\mathrm{T}}(x)-f(x))+\mathcal{O}(a^2)$$

$${\rho }^{\mathrm{T}}(x)\equiv f(x)$$

$$\text{Prob}\{n_i^{\mathrm{T}}\}={\rho }_i^{\mathrm{T}}=f_i=⟨n_i^{\mathrm{T}}{n}_{i+1}^{\mathrm{T}}⟩=\text{Prob}\{n_i^{\mathrm{T}},n_{i+1}^{\mathrm{T}}\}=\text{Prob}\{n_{i+1}^{\mathrm{T}}|n_i^{\mathrm{T}}\}\text{Prob}\{n_i^{\mathrm{T}}\}$$

$${\rho }_i^{\mathrm{T}}=\text{Prob}\{n_i^{\mathrm{T}}\}=0$$

Full text

Linear protrusions of cells, as for instance filopodia or stereocilia, perform multiple tasks in living organisms, ranging from cell migration and signal transduction to wound healing. They contain actin filaments cross-linked into bundles by actin-binding proteins Revenu et al. (2004); Mattila and Lappalainen (2008); Les Erickson et al. (2003); Nambiar et al. (2010), and molecular motors of the myosin family which interact with actin filaments and walk on them in a persistent, uni-directional fashion towards the tip of the protrusion Les Erickson et al. (2003); Nambiar et al. (2010); Salles et al. (2009); Hartman and Spudich (2012); Kerber and Cheney (2011); Kambara et al. (2006); Bird et al. (2014). These motors play an important role in the biological function of protrusions Revenu et al. (2004); Nambiar et al. (2010); Salles et al. (2009). In particular, they are known to localise to the tips of filopodia and stereocilia, and are (jointly) responsible for length control Mattila and Lappalainen (2008); Les Erickson et al. (2003); Nambiar et al. (2010); Salles et al. (2009); Hartman and Spudich (2012); Kerber and Cheney (2011); Kambara et al. (2006); Bird et al. (2014); Rzadzinska et al. (2004); Manor et al. (2011); Belyantseva et al. (2003).

Motivated by these observations, various models have been investigated. Some are detailed mathematical models addressing specific biological issues. These include the role of motor transport in shaping the concentration profile of G-actin at the base of protrusions Zhuravlev et al. (2012), the localisation of different proteins along stereocilia Naoz et al. (2008), and the effect of myosin X on filopodial growth Wolff et al. (2014). Complementary, simplified conceptual models have been studied asking how the interplay between active and diffusive transport in open systems may lead to non-equilibrium phase transitions and ensuing steady states with interesting correlations and nontrivial density profiles Krug (1991); Lipowsky et al. (2001); Parmeggiani et al. (2003, 2004); Evans et al. (2011). The latter are based on the totally asymmetric simple exclusion process (TASEP) MacDonald et al. (1968); Spitzer (1970), a lattice-gas model that, despite its simple structure, has become a paradigm for non-equilibrium dynamics Blythe and Evans (2007); Chou et al. (2011).

Here we present and analyse a conceptual model capturing two basic properties of the motion of persistent, plus-end directed motors inside narrow, elongated cellular protrusions. First, there is an interplay between two genuinely different types of dynamics: directed (active) transport with steric hindrance along actin filaments, and diffusive motion in the cytoplasm. These are coupled by particle exchange between the filament and the cytoplasm. Second, the half-closed geometry of cellular protrusions is special: At one end, the protrusions are connected to the cell body and thus to a reservoir, whereas everywhere else protein diffusion is confined by the cell membrane, so that mass conservation and resource limitation play an important role there. The combination of mass-conservation (closure) on the one hand and the interplay of equilibrium (diffusion) and non-equilibrium (active transport) processes on the other hand is intrinsically interesting to study as closure in a system entails no-flux boundary conditions that oppose currents from active transport. Here, we want to examine the interplay of these mechanisms with the help of an abstract model that is motivated biologically but has a level of description that makes it possible to understand all the processes accounted for. We identify generic mechanisms based on correlations and non-equilibrium physics that could be of importance for biological systems as cellular protrusions but, inevitably, predictions for biological systems are qualitative.

Specifically, we consider a two-lane lattice-gas model in a half-closed geometry [Fig. 1] in steady state, similar to Ref. Müller et al. (2005). One lane represents the actin filament and the second lane the cytoplasm. While there is a rich literature on the non-equilibrium dynamics of two-lane systems Klumpp and Lipowsky (2003); Popkov and Schütz (2003); Pronina and Kolomeisky (2004); Schmittmann et al. (2005); Müller et al. (2005); Pronina and Kolomeisky (2006); Reichenbach et al. (2007); Pronina and Kolomeisky (2007); Jiang et al. (2007); Reichenbach et al. (2008); Wang et al. (2008); Evans et al. (2009); Schiffmann et al. (2010); Evans et al. (2011); Melbinger et al. (2011); Saha and Mukherji (2013); Gupta and Dhiman (2014); Johann et al. (2014); Pinkoviezky and Gov (2014); Curatolo et al. (2016), very few of these studies address how the physics of non-equilibrium steady-states is affected by a half-closed geometry Müller et al. (2005); Pinkoviezky and Gov (2014). We are interested in the limit where the actin filament and the cytoplasm are coupled weakly by attachment and detachment processes, while Ref. Müller et al. (2005) focuses on the strong coupling limit. Our analyses show that, due to the closure of the system at the tip, there is only one type of density profile, namely, domain walls (DW) separating a high- and a low-density region. The limit where the width of the high-density region becomes microscopically small (of the order of a few lattice sites) corresponds to tip localisation. Furthermore, correlations in such systems have not received much attention, and we want to investigate nearest-neighbour (NN) correlations on the filament. In a biological context, this is related to efficient transport on actin filaments and to the significance of steric hindrance of motors there. We find that those correlations reach high values close to the DW, and the transport current along the filament is strongly reduced compared to its mean-field (MF) prediction. This suggests an important role for the cytoplasm, namely to transport the proteins to the tip. Conversely, active transport effectively sets up and maintains a gradient of motor proteins.

Our model consists of two coupled sublattices [Fig.  1], namely a TASEP and a SSEP (symmetric simple exclusion process) lane of $L+1$ sites $\,{\in}\,\{0,1,...,L\}$. The dynamics on lane $1$ (TASEP) are governed by a rate $\nu$ at which particles jump forward one site towards the right (tip) provided that the site in front of them is empty (exclusion). This corresponds to the directional motion of the motors on the filament that is oriented towards the tip. By convention, we measure time in units of $\nu$ (i.e. set $\nu\,{=}\,1$) and length in units of the system size (i.e. the lattice spacing is $a\,{=}\,\frac{1}{L}$ and the total length $La\,{=}\,1$). In lane $2$ (SSEP) particles jump non-directionally between neighbouring sites, at rate $\epsilon$ again respecting exclusion. This represents the diffusive motion in the cytoplasm that is taken to be effectively one-dimensional in the thin cylinder-like protrusions (the steady-state behaviour of a system with several lanes for diffusion arranged on a cylinder around the TASEP can be reduced to the steady-state behaviour of this model Sup ). Particles enter or exit the system only at site [math] of lane $2$ (base) but not at site 0 of lane 1. At rate $\alpha$ a particle is injected provided the site is empty and at rate $\epsilon$ a particle leaves the system. This reflects the exchange of motors between the protrusion base and the cell body. In the bulk both lanes are coupled via a rate $\omega$ at which particles jump from site $i\,{\in}\,\{0,1,...,L-1\}$ of lane $1$ to site $i$ of lane $2$ or vice versa (each respecting exclusion). Since the biochemistry at the tip is only poorly understood, at site $L$ we consider two extreme cases: In model A, particles jump from site $L$ of lane 1 to site $L$ of lane $2$ at rate $\beta\,{\gg}\,\omega$ (respecting exclusion) but not in the opposite direction. This describes a scenario where motors at the tip detach mainly due to the lack of a filament subunit in front of them. In model B, the exchange rates between the lanes at sites $L$ are the same as in the bulk. The comparison of both models stresses that seemingly minor changes in a non-equilibrium system may have a strong influence on the dynamics Reese et al. (2014), and highlights the relevance of the biochemistry at the filament tip for the motor density profile. In the following, when considering the continuum limit $a\rightarrow 0$, we focus on the mesoscopic limit for $\omega$ Parmeggiani et al. (2003, 2004), i.e. we keep $\Omega\,{=}\,\omega/a$ fixed, thus ensuring the number of jumps between lanes is of the same order as that of entry or exit events (persistent motors) and implementing weak coupling between the diffusive and the directed motion. For simplicity, we take the attachment and detachment rates to be equal. However, the qualitative results do not change for different attachment and detachment rates $\omega_{A}\neq\omega_{D}$ as long as both are still taken in the mesoscopic limit (not shown here).

A single TASEP exhibits three phases, namely a maximal current (MC), a low- (LD) and a high-density (HD) phase Krug (1991). Moreover, on the phase boundary between the LD and HD phase the steady-state profile is given by a DW that performs a random walk. Due to the closure at the tip, we do not find a MC phase in our system Sup , but instead observe localised DWs Lipowsky et al. (2001); Parmeggiani et al. (2003). That is, the generic steady-state TASEP profile $\rho^{\mathrm{T}}(x)$ for both models is given by a localised DW separating low density at the base from high density at the tip [Fig. 2]. For generic parameters, the filament current is thus comparatively small, and restricted to a small part of the system. This might be beneficial from a biological point of view, since every motor step on the filament consumes ATP. The position of the DW depends on the model parameters, and is shifted towards the tip (base) for smaller (higher) values of $\sigma\,{:=}\,\frac{\alpha}{\alpha+\epsilon}$ Sup . This ratio can be thought of as the motor density in the cell body, and thus the value for the cytoplasmic density at the protrusion base. The width of the DW decreases with increasing $L$ and with increasing distance from the tip Sup . Generically, when describing the DW as a random walker (see below), we observe that its motion is mainly confined to a small part of the system.

Hence we assume for the moment that the DW is fully localised and adopt a step-function ansatz for $\rho^{\mathrm{T}}(x)$ that holds to lowest order in $a$: $\rho^{\mathrm{T}}(x)\,{=}\,0$ for $x\,{\in}\,{[0,z}[$, and $\rho^{\mathrm{T}}(x)\,{=}\,1$ for $x\,{\in}\,]z,1]$. First, we determine the position $z$ of the DW to find out whether tip localisation (i.e.  $z\,{\approx}\,1$) occurs for certain parameter ranges or not. For this purpose we derive a mass-balance equation relating the total average occupancies of the TASEP ($\mathrm{T}$) and the SSEP ($\mathrm{S}$). We denote by $n^{\mu}_{i}$ the occupation number on site $i\,{\in}\,\{0,\ldots,L\}$ of lane $\mu\,{\in}\,\{\mathrm{T},\mathrm{S}\}$, i.e. we write $n^{\mu}_{i}\,{=}\,0$ if site $i$ of lane $\mu$ is empty and $n^{\mu}_{i}\,{=}\,1$ if it is occupied. Since $n^{\mu}_{i}\,{\in}\,\{0,1\}$, we have $\langle n^{\mu}_{i}\rangle\,{=}\,\mathrm{Prob}\{n^{\mu}_{i}\,{=}\,1\}$; ensemble averages are denoted by $\langle\ \rangle$. The steady-state condition in the bulk corresponds to a flux balance Sup : $(\rho^{\mathrm{T}}_{i-1}\,{-}\,f_{i-1})\,{-}\,(\rho^{\mathrm{T}}_{i}\,{-}\,f_{i})\,{=}\,\omega(\rho^{\mathrm{T}}_{i}\,{-}\,\rho^{\mathrm{S}}_{i})$ with the correlator $f_{i}\,{=}\,\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\rangle$ and the average occupancy $\rho^{\mu}_{i}\,{=}\,\langle n^{\mu}_{i}\rangle$. The difference of the TASEP currents, from site $i\,{-}\,1$ to site $i$ and from site $i$ to site $i\,{+}\,1$, (left-hand side) must equal the current between sites $i$ of the TASEP and SSEP, $\omega(\langle n^{\mathrm{T}}_{i}(1-n^{\mathrm{S}}_{i})\rangle\,{-}\,\langle n^{\mathrm{S}}_{i}(1-n^{\mathrm{T}}_{i})\rangle)$, (right-hand side). It follows that $f_{i}\,{=}\,f_{0}\,{+}\,\rho^{\mathrm{T}}_{i}\,{-}\,\rho^{\mathrm{T}}_{0}\,{+}\,\omega\sum_{j=1}^{i}(\rho^{\mathrm{T}}_{j}\,{-}\,\rho^{\mathrm{S}}_{j})$. At the base, there is no direct flux from the cell body into the filament; so the boundary condition is $f_{0}\,{=}\,\rho^{\mathrm{T}}_{0}+\omega(\rho^{\mathrm{T}}_{0}-\rho^{\mathrm{S}}_{0})$ and we find the following exact moment-identity

[TABLE]

The systems are both closed at the tip but, due to the different attachment and detachment behaviour, the boundary conditions read $f_{L-1}\,{=}\,\rho^{\mathrm{T}}_{L-1}\,{-}\,\beta\langle n^{\mathrm{T}}_{L}(1{-}n^{\mathrm{S}}_{L})\rangle$, and $f_{L-1}\,{=}\,\rho^{\mathrm{T}}_{L-1}\,{+}\,\omega(\rho^{\mathrm{S}}_{L}\,{-}\,\rho^{\mathrm{T}}_{L})$ for models A and B, respectively. Combining these with Eq. (1), the following exact mass-balance equations can be derived for model A

[TABLE]

and similarly for model B: $\omega\sum_{j=0}^{L}(\rho^{\mathrm{T}}_{j}{-}\rho^{\mathrm{S}}_{j})\,{=}\,0$. These equations relate the average occupancy on the two lanes in such a way that the total influx into the TASEP lane equals the total outflux from it Klumpp et al. (2007). Interestingly, these equations reveal that a global detailed balance holds for the total exchange between the two lanes, rather than local detailed balance for any pair of sites (cf. adsorption equilibrium in Ref. Müller et al. (2005)). The moment-identity and the mass-balance equations are useful in two ways: (i) By using the DW ansatz, we are able to find an analytic formula for the DW position. (ii) The moment-identity enables us to express covariances with respect to NNs on the TASEP lane in terms of densities.

To address the first issue, we determine the average density $\rho^{\mathrm{S}}(x)$ on the SSEP lane corresponding to the fully-localised-DW ansatz by solving the bulk equation for the SSEP. For that, we implement the continuum limit and, for model A, assume $\langle n^{\mathrm{T}}_{L}n^{\mathrm{S}}_{L}\rangle\,{\approx}\,\rho^{\mathrm{T}}_{L}\rho^{\mathrm{S}}_{L}$ Sup . Note, that we only need this MF approximation for the tip densities of model A. This is due to the otherwise symmetric attachment and detachment rates and the diffusive motion in the cytoplasm, for both of which the correlations drop out in the dynamical equation. With the resulting equation for $\rho^{\mathrm{S}}(x)$ we can then conclude $z\,{=}\,1\,{-}\,l\cosh^{-1}\left[\sigma\cosh(1/l)\right]$ and $z\,{=}\,1\,{-}\,l\sinh^{-1}\left[\sigma\sinh(1/l)\right]$ for models A and B, respectively, where $l\,{:=}\,\sqrt{{D}/{\omega}}$ with $D\,{:=}\,\epsilon a^{2}$ being the diffusion constant in the cytoplasm Sup ; $l$ can be understood as the typical dwell length for motors in the cytoplasm before attaching to the filament. Comparing these expressions with our stochastic simulation results Gillespie (1976) we find excellent agreement [Fig. 3].

The phase diagrams for the two models are qualitatively different. For model A, one can switch between the DW phase ($z\,{\ll}\,1$) and the tip-localisation phase ($z\,{\approx}\,1$) by only slightly increasing the dwell length $l$. In contrast, for model B, tip localisation is attained only as $\sigma\,{\rightarrow}\,0$, even for large $l$ Sup . This can be understood from the symmetry between the attachment and detachment processes in model B, which should also be reflected in a symmetry in occupancy between the filament and the cytoplasm. For large $l$, the cytoplasm becomes well-mixed with constant density $\sigma$. So, the average density on the filament is $\sigma$ as well, which is realised for a DW position $z\,{=}\,1\,{-}\,\sigma$. For model A, fast diffusion in the cytoplasm leads to rapid diffusion of the motors away from the tip. Therefore, the exit rate at the tip is high and motors quickly leave the tip of the filament, thus enabling tip localisation but reducing large jamming. For small but finite $\sigma$, there is tip localisation also if the typical dwell length is smaller than the system size, $l\,{\leq}\,1$. That is, even if the motors have a relatively small $D$, tip localisation can occur if the tip has an enhanced detachment rate (model A). For model B, even by increasing $l$ well beyond the system size, tip localisation occurs only for very low motor density $\sigma$ at the base.

These features may have interesting biological implications: A higher detachment rate at the actin filament’s end would promote tip localisation, but simultaneously avoid jamming. As a result, motor exchange between filament and cytoplasm occurs primarily near the tip (where the DW is located) and the motors at the tip can be continuously replenished by new ones delivered through the cytoplasm. Furthermore, energy consumption (ATP hydrolysis) is low, as the filament current is kept small and mainly restricted to the tip area. Transport to the tip is facilitated mainly by diffusion in the cytoplasm, which does not consume chemical energy. In summary, energy from ATP hydrolysis could efficiently be used to localise motors to the filament tip, while material transport is facilitated by diffusion in the cytoplasm Sup .

In this regard, the model is also interesting from a theoretical point of view, as it allows for the calculation of the filament current $J^{\mathrm{T}}_{i}\,{=}\,\rho^{\text{T}}_{i}\,{-}\,f_{i}$ that depends on NN correlations or, equivalently, the NN covariances

[TABLE]

on the TASEP lane. To go beyond MF, we use Eq. (1), which relates $f_{i}$ to the average densities, and find $\mathrm{Cov}_{i}\,{=}\,\rho^{\mathrm{T}}_{i}(1\,{-}\,\rho^{\mathrm{T}}_{i+1})\,{+}\,\omega\sum_{j=0}^{i}(\rho^{\mathrm{T}}_{j}\,{-}\,\rho^{\mathrm{S}}_{j}),$ which in the continuum limit translates to

[TABLE]

The value of the first term in Eq. (4) depends sensitively on the width and shape of the DW (the density profile $\rho^{\mathrm{T}}(x)$ is increasing with $x$, so that $\rho^{\mathrm{T}}(x)(1\,{-}\,\rho^{\mathrm{T}}(x{+}a))$ is maximal if both $\rho^{\mathrm{T}}(x)$ and $\rho^{\mathrm{T}}(x{+}a)$ are close to $0.5$). Therefore, one needs to refine the fully-localised-DW ansatz, by taking into account the stochastic dynamics of the DW position. Following Refs. Derrida et al. (1995); Kolomeisky et al. (1998) we consider the DW as a random walker with (site-dependent) hopping rates depending on currents and densities in the low- and high-density regions. For small $a$ we find

[TABLE]

where $W(z)\,{=}\,\sqrt{2\sigma al\sinh\left(z/l\right)}$  Sup . This fluctuation-corrected DW profile, as well as the covariance obtained from it, agree very well with our simulation data [Fig. 2]. If one uses a refined MF method instead, accounting for second-order spatial derivatives, the DW width and the strength of the correlations are both markedly underestimated Sup .

In general, covariances are non-zero only close to the DW, but there they can reach quite high values of around $0.2$. We have $\mathrm{Cov}_{i}\,{=}\,\rho^{\mathrm{T}}_{i}(1\,{-}\,\rho^{\mathrm{T}}_{i+1})\,{-}\,J_{i}^{\mathrm{T}}$, where the first term corresponding to the MF current dominates. Therefore, the actual current $J^{\mathrm{T}}_{i}\,{=}\,\omega\sum_{j=0}^{i}(\rho^{\mathrm{S}}_{j}-\rho^{\mathrm{T}}_{j})$ can be orders of magnitude smaller than the MF current, and scales with the particle exchange rate $\omega$ Sup . That is, density correlations dominate in such a way that they suppress the TASEP current substantially. This demonstrates that, even if basic properties of the single TASEP are captured well by MF theory, correlation effects in TASEP-based systems should be studied more closely and might lead to unanticipated features Popkov et al. (2003); Tsekouras and Kolomeisky (2008); Reithmann et al. (2016).

The covariances are non-negative everywhere, $\mathrm{Cov}_{i}\,{\geq}\,0$, so the conditional probabilities obey the inequalities $\mathrm{Prob}\{n^{\mathrm{T}}_{i}{=}1|n^{\mathrm{T}}_{i+1}{=}1\}\,{\geq}\,\mathrm{Prob}\{n^{\mathrm{T}}_{i}{=}1\}$ and $\mathrm{Prob}\{n^{\mathrm{T}}_{i+1}{=}1|n^{\mathrm{T}}_{i}{=}1\}\,{\geq}\,\mathrm{Prob}\{n^{\mathrm{T}}_{i+1}{=}1\}$ at any site, implying that particles typically form clusters. As a result, from the DW region onwards (where the covariances are highest), the mean time a particle spends at a certain site is increased considerably compared to freely moving motors. This is due to the effective jump rate on the filament, which is decreased by excluded-volume effects. Thus, in case of tip localisation, particles spend more time near the tip than in the main part of the filament. This prolongation of the residence time is further enhanced by the exclusion in the cytoplasm that prevents motors from detaching if the cytoplasmic tip density is high. This is especially important for model A, where the cytoplasmic tip density takes a value of 1, which is much higher than in the bulk Sup . Biologically, the extended residence time at the tip might facilitate the tasks of the motors or their cargo at the tip.

Our results all essentially rely on the exact moment-identity and the exact mass-balance equations. The derivation of both depends on the TASEP dynamics and the coupling between the two lanes, but not at all on the dynamics of the second (here SSEP) lane. Hence, those equations do not change if this dynamics on the second lane is modified, and we expect them to be useful for other model systems in which a TASEP lane is coupled to another lattice via attachment and detachment kinetics. To our best knowledge, the moment-identity has not been mentioned before. We believe that it might open doors in the understanding of correlations and the prediction of active currents also in more elaborate models. Furthermore, both equations could easily be generalised to the case of different attachment and detachment rates.

Our models could also be varied in other interesting ways. One could account for the three-dimensional geometry or for polymerisation and depolymerisation of the filament and the accompanying changes in length. Nevertheless, we expect that some of the phenomena seen here should be robust against such modifications. Tip localisation, which is mostly based on mass conservation at the tip, should still be present. Second, the TASEP current might still be suppressed, and the roles of the TASEP and the diffusive lane in being responsible for tip localisation and motor transport, respectively, should remain untouched. This seems to be supported by more elaborate models in a related context Zhuravlev et al. (2012); Wolff et al. (2014).

I Calculation of the density profiles

To set the stage, let us denote by $n^{\mu}_{i}$ the occupation number on site $i\in\{0,\ldots,L\}$ of lane $\mu\in\{\mathrm{T},\mathrm{S}\}$, i.e. we will write $n^{\mu}_{i}=0$ if site $i$ of lane $\mu$ is empty and $n^{\mu}_{i}=1$ if it is occupied ($\mathrm{T}$: TASEP, $\mathrm{S}$: SSEP). The common bulk master equations of the Markov processes, corresponding to model A and B as shown in Fig. 1, are given by

[TABLE]

where $\text{Prob}\{n^{\mu}_{i}\}$ denotes the probability that $n^{\mu}_{i}=1$ and $\text{Prob}\{n^{\mu}_{i},\overline{n^{\nu}_{j}}\}$ the one that $n^{\mu}_{i}=1$ and $n^{\nu}_{j}=0$. The term $\text{Prob}\{n^{\mathrm{T}}_{i-1},\overline{n^{\mathrm{T}}_{i}}\}-\text{Prob}\{n^{\mathrm{T}}_{i},\overline{n^{\mathrm{T}}_{i+1}}\}$ is due to the jump process on the TASEP lane that respects the exclusion property and occurs at bare rate $\nu=1$. The terms proportional to $\omega$ describe the exchange between the lanes, again respecting the exclusion. And finally, the term proportional to $\epsilon$ describes the diffusion on the SSEP lane. Note that we assume exclusion not only on the filament but also in the cytoplasm. This is based on the idea that, due to the finite size of particles, there should be a maximal number inside any finite volume element, introducing a carrying capacity $N_{\mathrm{max}}$. If we assume that the maximal effective attachment and detachment rate stay the same, i.e. if we assume that attachment happens at rate $\omega\left(n^{\mathrm{S}}\left(1-n^{\mathrm{T}}\right)\right)$ and detachment at rate $\omega\left(n^{\mathrm{T}}\left(N_{\mathrm{max}}-n^{\mathrm{S}}\right)\right)$ where $n^{\mathrm{S}}\in\{0,1,\dots,N_{\mathrm{max}}\}$, the case $N_{\mathrm{max}}$ finite but arbitrary can be reduced to $N_{\mathrm{max}}=1$ by redefinition of the parameters, and in the following we will focus only on the case $N_{\mathrm{max}}=1$. We will, however, come back to this case again in the last paragraph of the Supplemental Material when we discuss the case of several lanes for diffusion.

The bulk master equations can be rewritten in terms of averages over the occupation numbers by using that $\langle n^{\mu}_{i}\rangle=\mathrm{Prob}\left(n^{\mu}_{i}=1\right)=\mathrm{Prob}\{n^{\mu}_{i}\}$ and $\langle n^{\mu}_{i}\left(1-n^{\nu}_{i}\right)\rangle=\mathrm{Prob}\{n^{\mu}_{i},\overline{n^{\nu}_{i}}\}$ hold, as each site can be either empty or occupied by one particle/motor:

[TABLE]

where $f_{i}=\left\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\right\rangle$ is the nearest-neighbour correlator for the TASEP. Summing both equations we find

[TABLE]

for the time evolution of the combined density $\rho^{\mathrm{T}}_{i}+\rho^{\mathrm{S}}_{i}$ in the bulk. At the left boundary (base) we find

[TABLE]

for both models. For the TASEP lane there is only outflux from site [math] to site $1$ or exchange with site [math] of the SSEP lane. On the SSEP lane, there is influx at rate $\alpha$ (respecting the exclusion at site 0 of the SSEP), outflux with the diffusion rate $\epsilon$, diffusion between site [math] and site $1$ and exchange with site 0 of the TASEP. The behaviour at the right boundary (tip) differs between the models and is given by

[TABLE]

for model A, and by

[TABLE]

for model B. For model A, at site $L$ of the TASEP lane there is influx from the neighbouring site $L-1$ and outflux to site $L$ of the SSEP lane at bare rate $\beta$ respecting the exclusion. For site $L$ of the SSEP lane there is diffusion between sites $L-1$ and $L$ of the SSEP lane and influx from site $L$ of the TASEP lane. For model B, we have the same behaviour except that the asymmetric exchange between sites $L$ of the TASEP and SSEP lane is replaced by the symmetric exchange $\omega\left(\rho^{\mathrm{T}}_{L}-\rho^{\mathrm{S}}_{L}\right)$. Since the exchange terms drop out when considering the time derivative of $\rho^{\mathrm{T}}_{L}+\rho^{\mathrm{S}}_{L}$, at the tip of both models it holds

[TABLE]

and at the base

[TABLE]

where we introduced the local currents

[TABLE]

and

[TABLE]

for $i=0,\dots,L-1$ on the filament and in the cytoplasm, respectively. With that Eq. (S3) translates to

[TABLE]

Note that we can introduce currents $J^{\mathrm{T/S}}_{-1}$ and $J^{\mathrm{T/S}}_{L}$ as well. However, as the system is closed at the tip $J^{\mathrm{T/S}}_{L}\equiv 0$, and since there is no direct influx from the left into the TASEP $J^{\mathrm{T}}_{-1}=0$. In order to be consistent with the structure of Eq. (S11) where the currents appear as $J_{i-1}-J_{i}$, we define

[TABLE]

so that

[TABLE]

holds.

As mentioned above, we are interested in the steady-state behaviour of the system and therefore set $\partial_{t}\rho^{T/S}_{i}\equiv 0$ for all $i$. In particular, we have

[TABLE]

where we used Eqs. (S7), (S8), (S11). As a result,

[TABLE]

and the parameter $\sigma$ as used for the phase diagrams in Fig. 3 can be identified as the motor density in the cell body.

To proceed, let us now use the continuum limit where $a=\frac{1}{L}\rightarrow 0$. That is, we will replace the discrete lattice by the continuous space $[0,1]$ and the average occupation numbers are replaced by a continuous density

[TABLE]

With that, we can substitute $\rho^{\mu}(x_{i})\pm a\partial_{x}\rho^{\mu}(x_{i})+\frac{1}{2}a^{2}\partial_{x}^{2}\rho^{\mu}(x_{i})+\mathcal{O}(a^{3})$ for $\rho^{\mu}_{i\pm 1}$ and the currents are

[TABLE]

where $f(x)$ is the continuous version of $f_{i}$. In the steady-state Eq. (S11) translates to

[TABLE]

to first order in $a$ and thus, we have $\rho^{\mathrm{T}}(x)-f(x)=const$ to lowest order in $a$. Using Eq. (S7) or (S8) we conclude that $const=0$ or that the combined current $J^{\mathrm{T}}_{i}+J^{\mathrm{S}}_{i}$ must be constant and zero everywhere. This corresponds to the fact that the system is closed at the tip so that in steady state on average there is no influx into the system from the base. Therefore,

[TABLE]

or, equivalently,

[TABLE]

holds to lowest order in $a$ where we defined $\text{Prob}\{n^{\mathrm{T}}_{i},n^{\mathrm{T}}_{i+1}\}=\text{Prob}\left(n^{\mathrm{T}}_{i}=1,n^{\mathrm{T}}_{i+1}=1\right)$ and $\text{Prob}\{n^{\mathrm{T}}_{i+1}|n^{\mathrm{T}}_{i}\}=\text{Prob}\left(n^{\mathrm{T}}_{i+1}=1|n^{\mathrm{T}}_{i}=1\right)$ analogously to before. This implies that either

[TABLE]

or

[TABLE]

The latter implies that, whenever site $i$ is occupied, site $i+1$ is occupied with probability 1 as well. Hence, if there is some site $j$ occupied at some time, any site $i>j$ is occupied at that time as well. It follows that to lowest order in $a$ the density suddenly jumps from zero (S14) to maximal density

[TABLE]

This means that to lowest order the steady-state TASEP profile is given by a step function separating a region of [math] density on the left (towards the base) from a region of density $1$ on the right (towards the tip) and, in particular, there is no maximal current phase as it occurs for TASEP alone. Instead, the generic TASEP profile is given by a domain wall and we make the following fully-localised-DW ansatz:

[TABLE]

where $z\in[0,1]$ is the position of the domain wall (step) that we will determine later.

We will now use this ansatz to determine the steady-state density profile for the SSEP lane depending on $z$. For this purpose, let us go back to Eq. (S2) that is given by

[TABLE]

in the continuum limit in the steady state. This equation can be solved in the two regions $x<z$ and $x>z$ as

[TABLE]

The constants $A_{1}$, $A_{2}$, $B_{1}$ and $B_{2}$ can be determined from the boundary conditions: Eq. (S5) yields $\beta\left(\rho^{\mathrm{T}}_{L}-\left\langle n^{\mathrm{T}}_{L}n^{\mathrm{S}}_{L}\right\rangle\right)=0$ to lowest order in $a$ for model A, and assuming that $n^{\mathrm{T}}_{L}$ and $n^{\mathrm{S}}_{L}$ are uncorrelated (mean-field assumption that is supported by the simulations) we conclude that

[TABLE]

for model A, unless $\rho^{\mathrm{T}}_{L}=0$ and so $z=1$. For model B, we can use Eq. (S6) and find

[TABLE]

again unless $z=1$. This latter case needs to be treated separately. This can be done for instance by regarding the ”domain wall“ as a boundary layer with slope $(\rho^{\mathrm{S}}(1)_{\mathrm{right}}-\rho^{\mathrm{S}}(1)_{\mathrm{left}})/a$ where $\rho^{\mathrm{S}}(1)_{\mathrm{right}}$ is the density at the very last site $L$ and $\rho^{\mathrm{S}}(1)_{\mathrm{left}}$ is the limit $x\rightarrow 1$ of the low-density phase. Basically, this case can be understood as the limit where the domain wall is shifted to the right out of the system. As long as $\epsilon\in\mathcal{O}(a^{0})$, the generic density profile is given by a domain wall separating density [math] on the left from density $1$ on the right, that is $\rho^{\mathrm{T}}\in\{0,1\}$ holds to lowest order in $a$. In the limiting cases $\rho^{\mathrm{T}}\equiv 1$ the density reaches the value of $1$ in only a few lattice sites from the base, whereas in the case $\rho^{\mathrm{T}}\equiv 0$ the density has a very small spike (boundary layer) only at the tip. The details of this calculation are, however, out of the scope of this letter and we will carry on with the treatment of the case $z\neq 1$.

For both models we have Eq. (S13) for the boundary condition at the base. Requiring that both the diffusive steady-state density profile and the derivative thereof is continuous (it can be seen from Eq. (S2) that $\left(\rho^{\mathrm{S}}_{i+1}-\rho^{\mathrm{S}}_{i}\right)-\left(\rho^{\mathrm{S}}_{i}-\rho^{\mathrm{S}}_{i-1}\right)=\mathcal{O}(a)$ holds), the following expression for the SSEP profiles can be derived:

[TABLE]

for model A and

[TABLE]

for model B where

[TABLE]

with the diffusion constant $D=\epsilon a^{2}$ and

[TABLE]

To illustrate the differences between the two models, Fig. S1 shows typical density profiles on the SSEP lane for both cases. As one can see, the density profiles differ significantly between the two models. For model A, the density generically increases considerably from the protrusion base with density $\sigma$ towards the tip with density 1 (unless $z=1$). In contrast, the density for model B at the tip reaches a value that is only slightly larger than the value at the base. The fact that for model A the cytoplasmic density at the tip reaches a value very close to 1 is also the reason why for model A, the exclusion in the cytoplasm has a much higher influence on the residence time at the tip compared to model B where the cytoplasm is occupied quite homogeneously and where exclusion at the tip is not more important than in the bulk.

Using Eqs. (S17) and (S18) we now have expressions for the steady-state density profiles that only depend on one quantity, namely $z$, the position of the domain wall. Using the mass-balance equations [Eq. (2) and the corresponding one for model B], this enables us to find an analytic formula for $z$.

II Position of the domain wall

In order to determine the position of the domain wall, we go back to Eq. (S1). In the steady state this reduces to

[TABLE]

Hence, we can write

[TABLE]

Using Eq. (S4) we can rewrite $f_{0}=\rho^{\mathrm{T}}+\omega\left(\rho^{\mathrm{T}}_{0}-\rho^{\mathrm{S}}_{0}\right)$ and, thus,

[TABLE]

This is the moment-identity [Eq. (1)]. It is an important result for us as it gives an exact relation between the nearest-neighbour correlator $f_{i}=\left\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\right\rangle$ and the average densities. We will use it later to predict the covariances and the TASEP current in the system. This will allow us to make a sharp distinction between the actual TASEP current $J^{\mathrm{T}}_{i}$ and the mean-field current $J^{\mathrm{T}}_{\mathrm{MF},i}=\rho^{\mathrm{T}}_{i}\left(1-\rho^{\mathrm{T}}_{i+1}\right)$.

From the moment-identity, Eq. (S19), we can also easily deduce the mass-balance equations [Eq. (2) and the corresponding one for model B] by using the boundary conditions at the tip [Eq. (S5) for model A and (S6) for model B] and the moment-identity for $i=L-1$:

[TABLE]

In order to determine $z$ we will now insert the domain wall ansatz, Eq. (S16), for the TASEP together with the respective SSEP profile into the respective mass-balance equation. For this purpose, we need to write the mass-balance equations in the continuum limit: Using that $\int_{0}^{1}\mathrm{d}x\ \rho^{\mathrm{T}}(x)=1-z$ and that $\sum_{i=0}^{L}\rho^{\mu}_{i}\approx\int_{0}^{L}\mathrm{d}i\ \rho^{\mu}(ia)=\frac{1}{a}\int_{0}^{1}\mathrm{d}x\ \rho^{\mu}(x)$ we find to lowest order in $a$:

[TABLE]

where for model A the mean-field assumption $\left\langle n^{\mathrm{T}}_{L}n^{\mathrm{S}}_{L}\right\rangle\approx\rho^{\mathrm{T}}_{L}\rho^{\mathrm{S}}_{L}$ was used again, together with Eq. (S5). Inserting the explicit densities (S17) or (S18) into these expressions, integrating and cancelling common factors results in

[TABLE]

Those directly lead to the expressions for the domain wall given in the main text:

[TABLE]

We see that $z$ is shifted towards the tip (base) for smaller (higher) values of $\sigma:=\frac{\alpha}{\alpha+\epsilon}$. It is interesting to note that for model A the limit $z\nearrow 1$ is reached for finite values of the parameters $\sigma$, $\omega$ and $\epsilon$ when $\sigma\cosh\left(\frac{1}{l}\right)=1$ or $l=1/\cosh^{-1}(1/\sigma)$, whereas for model B, $z\nearrow 1$ is only possible in the limit $\sigma\rightarrow 0$. This is also reflected in the fact that $\sigma\cosh\left(\frac{1}{l}\right)=\cosh\left(\frac{1-z}{l}\right)$ lacks a real solution for $z$ if $\sigma\cosh\left(\frac{1}{l}\right)<1$. That is, for model A there are parameter regimes where there is no domain wall and where to lowest order the density is zero everywhere. For model B, this case only occurs in the limit where the cytoplasmic density at the base is zero. Even for large $l$ we have $z=1-l\sinh^{-1}\left(\sigma\sinh\left(\frac{1}{l}\right)\right)\approx 1-l\sinh^{-1}\left(\frac{\sigma}{l}\right)\approx 1-\sigma$ and thus, $z<1$ unless $\sigma=0$.

III Refined mean-field TASEP density profile

Certainly, the TASEP profile is not just given by a plain step function but has a smooth form that separates the two density regions via an intermediate region where the density increases strongly but whose width is finite. This is due to fluctuations in the stochastic system that soften this transition. Our first approach to resolve this finite width relies on an idea similar to the method used in Saha and Mukherji (2013), namely to split the system into a part around the domain wall and the parts further away from it. By this, only parts of the terms contribute, respectively, and assuming a steep domain wall for the TASEP one can neglect the terms stemming from the exchange with the diffusive lane in this narrow region: Reconsidering equation (S1) we realize that in the steady state it can be approximated as

[TABLE]

where we used the continuum limit and a mean-field assumption for the nearest-neighbours on the TASEP: $f_{i}=\left\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\right\rangle\approx\rho^{\mathrm{T}}_{i}\rho^{\mathrm{T}}_{i+1}$. If we want to investigate the vicinity of the domain wall $x\approx z$ this can be done by considering $\tilde{x}:=(x-z)/a$ that is very large away from the domain wall. Using this coordinate system the above equation looks like follows:

[TABLE]

Dropping the term of order $\mathcal{O}(a)$, we end up with

[TABLE]

where $const$ can be estimated from the boundary conditions for $\tilde{x}\rightarrow\pm\infty$: we have $\rho^{\mathrm{T}}(x=0)=1-\rho^{\mathrm{T}}(x=1)=0$ and $\left.\partial_{x}\rho^{\mathrm{T}}(x)\right|_{x=0/1}\approx 0$, and hence, using $\tilde{x}\rightarrow\pm\infty$ we find that $const\approx 0$. As a result, after integrating $\rho^{\mathrm{T}}(\tilde{x})\left(\rho^{\mathrm{T}}(\tilde{x})-1\right)+\frac{1}{2}\partial_{\tilde{x}}\rho^{\mathrm{T}}(\tilde{x})=0$ with respect to $\tilde{x}$, using that $\rho^{\mathrm{T}}(z)=\frac{1}{2}$, and rescaling we find

[TABLE]

for the refined mean-field TASEP density profile. Contrary to the step profile, Eq. (S16), it exhibits a finite width that scales like $a$ or inversely proportionally to $L$. But as one can infer from Fig. 2 it still underestimates the actual width due to fluctuations that are ignored by the mean-field assumption for the nearest-neighbours on the TASEP.

IV Domain wall theory

This is why we also pursue another approach to refine our prediction for the domain wall, namely by going beyond mean-field theory and treating the domain wall as a random walker that moves in a non-uniform potential with reflecting boundaries. This technique has been introduced in Ref. Derrida et al. (1995) for the TASEP and, more generally, for shock waves in Ref. Kolomeisky et al. (1998). It has been applied to TASEP-LK in Refs. Parmeggiani et al. (2003, 2004) using non-uniform hopping rates. We follow their ideas and describe the domain wall by site-dependent hopping rates

[TABLE]

to the left and right from site $i$, respectively. Here, $J_{\mathrm{LD},i}$ ($J_{\mathrm{HD},i}$) denotes the steady-state TASEP current at site $i$ under the assumption that site $i$ is in the low-density (high-density) phase, i.e. on the left (right) of the domain wall. $\rho^{\mathrm{T}}_{\mathrm{LD},i}$ and $\rho^{\mathrm{T}}_{\mathrm{HD},i}$ are the respective TASEP densities in the low- and high-density phase. The heuristic reason why the hopping rates have the above form, Eq. (S24), is that the low-density current $J_{\mathrm{LD},i}$ should just correspond to the current arriving at the left of the domain wall and causing the domain wall to move one step to the left. Thereby, at site $i$ the density increases by $\rho^{\mathrm{T}}_{\mathrm{HD},i}-\rho^{\mathrm{T}}_{\mathrm{LD},i}$ at rate $w_{l,i}$. Similarly, the high-density current $J_{\mathrm{HD},i}$ should correspond to the current leaving the region at the right of the domain wall, and causing the domain wall to move one step to the left. This happens at rate $w_{r,i}$ and the density decreases by $\rho^{\mathrm{T}}_{\mathrm{HD},i}-\rho^{\mathrm{T}}_{\mathrm{LD},i}$.

From the detailed balance condition $p_{s,i}w_{r,i}=p_{s,i+1}w_{l,i+1}$ it is easy to see that the stationary distribution of the shock position $p_{s,i}$ satisfies

[TABLE]

We will show next that $w_{l,i_{z}}=w_{r,i_{z}}$ holds for $i_{z}$ being the discrete domain wall position ($i_{z}=z/a$): From Eq. (S1) as well as from the definition of the local TASEP current, Eq. (S9), we know that $J^{\mathrm{T}}_{i}=J^{\mathrm{T}}_{i-1}+\omega\left(\rho^{\mathrm{S}}_{i}-\rho^{\mathrm{T}}_{i}\right)$. Hence, we have

[TABLE]

As mentioned after introducing $J_{\mathrm{LD/HD},i}$ in Eq. (S24), $J_{\mathrm{LD/HD},i}$ denotes the TASEP current at site $i$ assuming that this site is in the low-density (high-density) region. That is, we can calculate $J_{\mathrm{LD},i}$ ($J_{\mathrm{HD},i}$) by using the above iterations starting from $i=0$ ($i=L-1$) and assuming that between site [math] (site $L-1$) and site $i$ the densities are given in the low-density (high-density) phase:

[TABLE]

We want to compare those two currents right at the domain wall $i=i_{z}$. For this purpose, we use that $J_{\mathrm{LD},j}=J^{\mathrm{T}}_{j}$ for $j<i_{z}$ and $J_{\mathrm{HD},j}=J^{\mathrm{T}}_{j}$ for $j>i_{z}$ since we chose $z$ in such a way that the densities to the left (right) of the fixed domain wall $i_{z}$ are in the low-density (high-density) phase. As a result,

[TABLE]

for both models: $J^{\mathrm{T}}_{L-1}=\beta\left(\rho^{\mathrm{T}}_{L}-\left\langle n^{\mathrm{T}}_{L}n^{\mathrm{S}}_{L}\right\rangle\right)$ for model A and $J^{\mathrm{T}}_{L-1}=\omega\left(\rho^{\mathrm{T}}_{L}-\rho^{\mathrm{S}}_{L}\right)$ for model B and, hence, the above holds due to the respective mass-balance equation. It follows that

[TABLE]

and, thus, $\ln\left(\frac{w_{l,i_{z}}}{w_{r,i_{z}}}\right)=0$. Similarly, we have $\ln\left(\frac{w_{l,j}}{w_{r,j}}\right)<0$ for $j<i_{z}$ and $\ln\left(\frac{w_{l,j}}{w_{r,j}}\right)>0$ for $j>i_{z}$ such that the domain wall preferably walks towards the (fixed) domain wall position $i_{z}$ than away from it. Therefore, the exponential in Eq. (S25) has a maximum at $i=i_{z}$ and taking the continuous version of this equation

[TABLE]

we can use the method of steepest descent to approximate

[TABLE]

where $w_{l}(x)$ and $w_{r}(x)$ is the continuous version of $w_{l,i}$ and $w_{r,i}$, respectively. Since $\frac{w_{l}(x)}{w_{r}(x)}=\frac{J_{\mathrm{LD}}(x)}{J_{\mathrm{HD}}(x)}$, we conclude

[TABLE]

where we used that $J_{\mathrm{LD}}(z)=J_{\mathrm{HD}}(z)$ [Eq. (S28)] and defined the width $W(z)$ by

[TABLE]

To find a formula for $\rho^{\mathrm{T}}(x)$ we realize that

[TABLE]

since as along as the shock position is left (right) of $x$ there is a high-density (low-density) region at $x$ and since we can approximate $\rho^{\mathrm{T}}_{\mathrm{HD}}(x)\approx 1$ and $\rho^{\mathrm{T}}_{\mathrm{LD}}(x)\approx 0$ to lowest order. Using Eq. (S30) we can determine $\rho^{\mathrm{T}}(x)$ as

[TABLE]

where $\mathrm{erf}$ is the error function. Here, we used that we need to normalize $p_{s}(x)$ such that $\int_{0}^{1}\mathrm{d}\tilde{x}\ p_{s}(\tilde{x})=1$.

In order to use this formula we still need to determine the width $W(z)$ more concretely. For this purpose, we must make some assumption on how we treat the diffusive (SSEP) lane when considering the low- and high-density currents. We are not aware of previous attempts that apply domain wall theory in case where a TASEP is coupled to another lattice that is occupied stochastically as well. Then the attachment and detachment rates for the TASEP not only depend on the TASEP occupancy but also on the occupancy on the other lattice. The difficulty with this situation is that it is a priori not clear if one should assume that the occupancy on the coupled lattice is fixed, that is, does not change very much if the domain wall is shifted towards the left or right, or if one should look at the momentary steady-state density corresponding to a certain position of the domain wall. For simplicity, we chose the first ansatz that seems to be working well. That means, we use Eqs. (S17) and (S18) for the SSEP density profile with the calculated fixed position of the domain wall $z$ given by Eqs. (S20) and (S21) and we assume that this SSEP density does not change considerably when the domain wall is shifted shortly to some other position $\neq z$, i.e. $\rho^{\mathrm{S}}_{\mathrm{LD}}(x)=\rho^{\mathrm{S}}_{\mathrm{HD}}(x)=\rho^{\mathrm{S}}(x)$.

With that in mind, using Eqs. (S26) and (S27), we find that

[TABLE]

and obtain

[TABLE]

where we used that $\rho^{\mathrm{T}}_{\mathrm{LD}}(x)\approx 0$ and $\rho^{\mathrm{T}}_{\mathrm{HD}}(x)\approx 1$. As a result, $J_{\mathrm{LD}}^{\prime}(x)=\Omega\rho^{\mathrm{S}}(x)$ and $J_{\mathrm{HD}}^{\prime}(x)=\Omega\left(\rho^{\mathrm{S}}(x)-1\right)$ so that

[TABLE]

where the last equality results from integrating $\rho^{\mathrm{S}}$ and then using the defining equation for $z$. Interestingly, the result in this form agrees for both models.

We observe that the width $W(z)$ increases with increasing $z$ or decreasing distance from the tip. Furthermore, in the limit $a\rightarrow 0$ it holds that $l\rightarrow 0$ (if we keep $\epsilon$ constant) so that we can approximate $\mathrm{arccosh}\left(\sigma\cosh\left(\frac{1}{l}\right)\right)\approx\frac{1}{l}+\ln\left(\sigma\right)$ and, similarly, $\mathrm{arcsinh}\left(\sigma\sinh\left(\frac{1}{l}\right)\right)\approx\frac{1}{l}+\ln\left(\sigma\right)$. Thus, $z\approx-l\ln\left(\sigma\right)$ for both models in the limit $a\rightarrow 0$ [Eqs. (S22), (S23)]. With that we find

[TABLE]

in the limit $a\rightarrow 0$ and so the width of the domain wall decreases with decreasing $a$ or increasing number of sites $L$.

Note that in order to calculate the position of the domain wall $z$ by the mass-balance equation, we can use the step function instead of the refined profiles. This is due to the fact that there the errors more or less cancel since only the sum $\sum_{i=0}^{L}\rho^{\mathrm{T}}_{i}$ over the densities on the left and on the right of the domain wall enters. However, this approximation gets worse the closer the calculated $z$ is to 0 or 1 since then, the errors are not symmetric anymore. As a result, also the approximation, Eq. (S31), deteriorates. Certainly, if there is no solution for $z$, e.g. for some parameters in model A, we can not use this ansatz either. Furthermore, if $\epsilon$ is too large, we expect our approximation that $\rho^{\mathrm{S}}(x)$ does not change when the domain wall is shifted to deteriorate as well.

V Covariances and Currents

Another point we want to emphasize is that in the system the nearest-neighbour correlations of the TASEP significantly modify the TASEP current in the sense that the mean-field prediction for it overestimates the current by orders of magnitude. The reason lies in the covariances that reach very high values of around $0.2$ in the region of the domain wall (see for example Fig. 2). As we have shown above in Eq. (S19) the nearest-neighbour correlator $f_{i}=\left\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\right\rangle$ is given by $f_{i}=\rho^{\mathrm{T}}_{i}+\omega\sum_{j=0}^{i}\left(\rho^{\mathrm{T}}_{j}-\rho^{\mathrm{S}}_{j}\right)$. Therefore, the covariances are given by

[TABLE]

Due to the first term, the mean-field current, which is large ($\approx 0.2$) in the region of the domain wall, also the covariances are non-zero and large there. However, since the mean-field current gives the main contribution to the covariances, this implies that the actual current

[TABLE]

is much smaller than predicted from a mean-field theory. To illustrate this further, in Fig. S2 we show a comparison of the mean-field current on the left as obtained from simulations using the average densities, and of the actual current on the right as obtained directly from the simulation (filled red (dark grey) circles), from the average densities and nearest-neighbour correlations in the simulation (filled orange (grey) circles) and from Eq. (S32) (dotted green (light grey) curve). Certainly, this discrepancy scales with $\Omega$ (and the other parameters) and the case shown here is extreme as $\Omega=0.001$ was chosen. However, the fact that both currents differ by orders of magnitude is robust and also occurs for $\Omega$ being of the order of the other jump rates [Fig. S3]. This shows that a mean-field description fails to capture essential properties of our model system.

Furthermore, from the fact that the covariances are non-zero everywhere $\mathrm{Cov}_{i}\geq 0$ and that we can rewrite

[TABLE]

we can conclude that $\mathrm{Prob}\{n^{\mathrm{T}}_{i}|n^{\mathrm{T}}_{i+1}\}\geq\mathrm{Prob}\{n^{\mathrm{T}}_{i}\}$. And from using a similar argument and only interchanging the roles of $i$ and $i+1$ we find $\mathrm{Prob}\{n^{\mathrm{T}}_{i+1}|n^{\mathrm{T}}_{i}\}\geq\mathrm{Prob}\{n^{\mathrm{T}}_{i+1}\}$. This tells us that, in particular from the domain wall region onwards where covariances are high, the motors preferentially cluster and hinder each other effectively so that the mean time a particle spends at a certain site is increased considerably with respect to freely moving motors. Eventually, this leads to the substantial difference between the actual current and the mean-field current.

Potentially, those results are also important from a biological point of view. In case of an enhanced detachment rate at the filament’s tip (model A), tip localisation is facilitated [Fig. 3] and large traffic jams can be avoided. In this case there is only high filament density in a small region around the tip so that motors spend more time in the tip region than in the other part of the filament, a feature that might be favoured biologically as then the motors or their cargo might have more time to perform the necessary tasks at the tip. Furthermore, transport to the tip might be strongly promoted by diffusion in the cytoplasm (SSEP) whose currents in both the direction of the tip and of the base significantly exceed the filament (TASEP) current [Fig. S3]. Certainly, in case of our original model where exclusion in the cytoplasm occurs at occupancy 1, the motors cannot really bypass each other, so that the motors do not circulate very often. This, however, should be greatly enhanced in the generalised model with carrying capacity $N_{\mathrm{max}}>1$ (see next chapter) or in a real biological system where the motors can overtake each other in the cytoplasm. The fact that the filament current is strongly suppressed by excluded volume effects might also be beneficial from a biological point of view as every motor step on the filament consumes ATP contrary to diffusion in the cytoplasm. Thus, transport of motors by the cytoplasm rather than by the filament might be advantageous energetically. Moreover, both the attachment and detachment current $J_{\mathrm{on}}$ and $J_{\mathrm{off}}$ are mainly restricted to the tip area [Fig. S3] so that motor exchange between the filament and the cytoplasm occurs primarily near the tip where the cargo is used. Taken together, tip localisation and the suppression of the filament current might be beneficial from a biological point of view as then energy consumption is low and motors could be efficiently transported to the tip by the cytoplasm. Near the tip they attach to the filament, have an enhanced residence time on the filament due to steric hindrance between the motors and then detach at the tip back into the cytoplasm.

VI Cylindrical geometry with several lanes for diffusion

Finally, we want to deal with a generalisation of our model where instead of one lane for diffusion we have several lanes for diffusion arranged on a cylinder around the filament. For an illustration of the case with $N_{\mathrm{diff}}=4$ lanes for diffusion please refer to Fig. S4. The three-dimensional view is shown on the left-hand side, the profile is illustrated on the right-hand side. As before the dynamics on the filament (red) are given by a TASEP lane with jump rate $\nu=1$ but now the dynamics in the cytoplasm is modelled by several lanes for diffusion (blue), each with diffusive rate $\epsilon$ and respecting the exclusion. The lanes for diffusion are arranged in a cylinder-like fashion around the TASEP lane, and each can interact with the TASEP lane by attachment/detachment processes at rate $\omega$ respecting the exclusion. Apart from the diffusion along the cylinder axis there is also lateral diffusion between neighbouring lanes of diffusion. This diffusion happens at rate $\epsilon_{\textbf{lat}}$ and again respects the exclusion. At the base (not shown) there is influx at rate $\alpha$ into every lane for diffusion and outflux at the diffusive rate $\epsilon$.

In the following, we will show that we can reduce the steady-state behaviour of this more elaborate model to the steady-state behaviour of our model by scaling the parameters for diffusion and attachment/detachment by the number of lanes for diffusion $N_{\mathrm{diff}}$: $\omega\rightarrow\omega N_{\mathrm{diff}}$, $\epsilon\rightarrow\epsilon N_{\mathrm{diff}}$. Here, we only consider the case where we have the same attachment and detachment kinetics at the tip than in the bulk (model B), but an analogous generalisation and argument can be done for model A or different rates for attachment and detachment. Let us denote by $n^{\mathrm{T}}_{i}$ the occupancy at site $i$ of the TASEP (as before) and by $n^{\mathrm{S},m}_{i}$ the occupancy at site $i$ of the $m=1,\ldots,N_{\mathrm{diff}}$-th lane for diffusion. We start by $m=1$ at an arbitrary lane for diffusion and then consecutively number the lanes for diffusion in clockwise order. Since the lanes are arranged on a cylinder we have periodic boundary conditions and identify $n^{\mathrm{S},N_{\mathrm{diff}}+1}_{i}\equiv n^{\mathrm{S},1}_{i}$. With this convention we have the following bulk master equations, written straight away in terms of the averages over occupancies $\rho^{\mu,m}_{i}=\left\langle n^{\mu,m}_{i}\right\rangle$ and $f_{i}=\left\langle n^{\mathrm{T}}_{i}n^{\mathrm{T}}_{i+1}\right\rangle$:

[TABLE]

At the base we find

[TABLE]

and at the tip

[TABLE]

Due to the exclusion property on both the TASEP as well as on the lanes for diffusion the state space of the system is finite. Furthermore, it is an irreducible continuous-time Markov process so there exists a unique steady-state. This is important since then the cylindrical symmetry of the system must be reflected in this steady-state. This is why we can assume that in the steady-state $\left\langle n^{\mathrm{S},m}_{i}\right\rangle=\left\langle n^{\mathrm{S},n}_{i}\right\rangle$ and $\langle n^{\mathrm{S},m}_{i}n^{\mathrm{T}}_{i}\rangle=\langle n^{\mathrm{S},n}_{i}n^{\mathrm{T}}_{i}\rangle$ holds for all $m,n$. As a result, the terms proportional to $\epsilon_{\mathrm{lat}}$ drop out in Eqs. (S33)-(S35) and we can define the total occupancy at site $i$ in the cytoplasm

[TABLE]

and the corresponding average $\rho^{\mathrm{S}}_{i}$ as the sum of the (average) occupancies at sites $i$ of all the lanes for diffusion taken together. Using this quantity we can rewrite Eqs. (S33)-(S35) in steady-state as follows:

[TABLE]

for the bulk master equation,

[TABLE]

for the base and

[TABLE]

for the tip. From these equations it becomes apparent that the cylindrical system with $N_{\mathrm{diff}}$ lanes for diffusion is equivalent to the above mentioned generalisation of our model where exclusion on the filament does not happen at occupancy of 1 but at a maximal occupancy or carrying capacity of $N_{\mathrm{max}}=N_{\mathrm{diff}}$ that is reflected in the term $N_{\mathrm{diff}}-n^{\mathrm{S}}_{i}$. To show that the steady-state behaviour of those two generalisations of our model are cast by the steady-state behaviour of our model, we next introduce the quantity

[TABLE]

and the corresponding average occupancy $\tilde{\rho}^{\mathrm{S}}_{i}$ at site $i$ of one lane for diffusion. In terms of $\tilde{n}^{\mathrm{S}}_{i}$ we find

[TABLE]

for the bulk and

[TABLE]

for the base and the tip where we used that $\langle\tilde{n}^{\mathrm{S}}_{i}(1-n^{\mathrm{T}}_{i})\rangle-\langle n^{\mathrm{T}}_{i}(1-\tilde{n}^{\mathrm{S}}_{i})\rangle=\tilde{\rho}^{\mathrm{S}}_{i}-\rho^{\mathrm{T}}_{i}$.

In summary, if we replace $\omega\rightarrow\omega N_{\mathrm{diff}}$, $\epsilon\rightarrow\epsilon N_{\mathrm{diff}}$ and $\alpha\rightarrow\alpha N_{\mathrm{diff}}$ in our model we can deduce from the steady-state behaviour of our model the steady-state behaviour of these generalised models. Certainly, to stay within the scope of our considerations, this implies that we need $\omega N_{\mathrm{diff}}$ to be much smaller than 1 (order $a$) and $\epsilon N_{\mathrm{diff}}$ to be of the order of 1 implying that the ”total attachment/detachment rate“ $\Omega N_{\mathrm{diff}}$ and the ”total diffusion rate“ $\epsilon N_{\mathrm{diff}}$ should be of the order of the hopping constant $\nu=1$.