I. Introduction

The mammalian endoplasmic reticulum (ER) is the single largest intracellular structure (see sketch in Fig. 1a). The organelle is made up of membranous sheets interconnected with the nuclear envelope and branching out into a planar network of tubules extending throughout the cell periphery (Fig. 1b) [1]. The ER dynamics on a second scale include the cytoskeleton-assisted tubular network restructuring [2, 3] and an interconversion between two distinct forms, including a narrower form covered by membrane curvature-promoting proteins [4]. The ER morphology and its dynamics presumably enable and facilitate its functions: the ER is responsible for the production, maturation and quality-controlled folding of secretory and membrane proteins, which constitute approximately a third of the cell’s proteome [5]. The organelle’s membranes also harbour the lipid biosynthesis machinery, while its lumen stores calcium. The contiguous nature of the ER is believed to ensure an efficient delivery of all these components across the cell periphery.

FIG. 1:

Timely transport of the content within the ER is therefore integral to the function of the cell. The geometry and dimensions of several cell types with extensive ER-containing projections (e.g. neurons and strocytes) pose a kinetic challenge for material distribution with physiological timing. These considerations predict the need for an active luminal transport to ensure timely material homogenisation across the vast ER. Empirically, the active nature of the ER luminal transport is supported by the sensitivity of GFP mobility (measured by FRAP) to ATP depletion [6, 7]. Thus past measurements indicate that the transport of proteins in ER is not consistent with Brownian motion [6, 8]. Using single particle tracking of ER luminal markers, an active solute transport in the planar ER network was observed directly [7]. A photoactivation chase technique also measured a superdiffusive behaviour of luminal material spread through the ER network [9].

However, the mechanism for generating ER luminal flows remains unclear. Understanding the mode of material exchange across the organelle is crucial for rationalising the ER shaping defects-related neuronal pathologies [10], identifying factors controlling ER transport and informing the development of ER transport modulation approaches with health benefits. Based on ER marker velocity fluctuations measured in single particle tracking and the detection of transient narrowing points in the tubules by improved super-resolution and electron microscopy (EM) [7], it has been postulated that these active flows may result from the stochastic contractility (pinching and un-pinching) of ER tubules at specific locations along their lengths (see sketch in Fig. 1e); other plausible mechanisms for flow generation were also considered. However, measurements for testing this pinching hypothesis are currently inaccessible, due to limitations in space-time resolution of live cell microscopy; the live cell-compatible super-resolution techniques achieve resolution of ∼80 nm at the relevant speed, while the tubular radius is estimated in the range 30−60 nm [9, 11, 12]. Improvement or resolution currently can only be achieved by trading off speed. To circumvent these experimental difficulties, in the current study we use mathematical modelling to quantitatively analyse the relevant scenarios of actively contractility-driven flows. Specifically, we explore how various sets of spatio-temporal parameters of ER contractility may produce flows facilitating solute transport in quantitative agreement with experimental measurements.

Specifically, assuming that the flow is driven by tubule contractions (shown schematically in Fig. 1e), we construct a physical and mathematical model of the ER network and solve it for the flows inside the network (description in Methods §IV). We use our model to carry out numerical simulations to study the motion of Brownian particles carried by these flows and show that the tubule pinching hypothesis is not supported by the results of our model (§II A), a result independent of the network geometry (§II B). The failure of active pinching to drive strong flows can be rationalised theoretically by deriving a rigorous upper-bound on the rate of transport induced by a single pinch, (§II C), including possible coordination mechanisms (§II D). Only by increasing both the length of pinches and their rate can beyond admissible values can we produce particle speeds in agreement with measurements (§II E). We then explore in §II F two different hypotheses as possible explanations for the active ER flows, first the pinching of the junctions between tubules (Fig. 1c) and then the pinching of the two types of ER sheets, perinuclear (Fig. 1f) and peripheral (Fig. 1d). We investigate the conditions under which both would be consistent with the experimental measurements of active transport in Ref. [7].

The question of how the endoplasmic reticulum homogenises its content across cell expanses is an open problem which extends to the fundamentals of cell biology. Even though it is highly debated in the field, it is extremely difficult to study due to the enormous technical limits of directly observing intraorganellar nanofluidics. Thus, our study represents a meaningful effort to break through this impasse by conducting a meticulous analysis of nanofluidic scenarios. The findings we present constitute the best currently possible endeavor to shed light on this challenging problem. Our results do suggest that the biological origin of solute transport in ER networks remains open and call for extensive empirical exploration of the alternative mechanisms for flow generation.

II. Results

A. Tubule contractility-driven ER luminal motion yields inadequate transport kinetics

To assess the kinetics of particle motion in the lumen of tubular structures in response to their contractility, we generated an in-silico simulation model of the process. The model incorporates local calculations for the low-Reynolds number hydrodynamics of a contracting tubule, assuming in the first instance the no-slip boundary conditions at the tubule walls (i.e. Poiseuille hydrodynamics), into a global analysis using Kirchhoff’s Laws and standard graph theoretical results of the flows throughout the network geometry (see Methods §IV for details).

We initially implemented the model’s numerical simulations of particle transport in a reconstructed ER network of a COS-7 cell [7] (which we label C0, see §IV A) with tubules locally contracting stochastically according to the spatiotemporal parameters suggested by microscopy measurements [7]. We will therein refer to these contractions, illustrated in Fig. 1e, as ‘pinching’, with relevant parameters that include the duration and frequency of pinch events, and the length that the pinch sites occupy along the tubule (for details of the pinching kinematics, see Methods §IV B). An estimate reveals that these pinches are indeed afforded by biologically realistic forces, of the order of 30 pN (see Methods §IV G).

The fluid flows in the edges of the network (model in §IV C computed as detailed in §IV D) reveal a rapid direction alternation of luminal currents (on average with a frequency of approximately 50 Hz), as reflected in the changes of the axial velocity sign (Fig. 2a), with an average flow speed of 1.3 μm/s (see also Supplementary Video S1). Further, the resulting instantaneous speed distribution of Brownian particles advected by these flows (methodology in §IV E) is considerably shifted relative to the experimental counterpart [7] towards lower values (Fig. 2b). Similarly, the distribution of average edge traversal speeds (defined precisely in §IV F) from our simulations (solid blue line in Fig. 2c) is lower than the experimentally measured speed distribution by an order of magnitude (Fig. 2c, inset).

FIG. 2:

Moreover, the results for all measures of transport under pure diffusion, in the absence of pinching-induced flows, were virtually indistinguishable from those where the pinchingdriven flows are included (see Fig. 2b-c). Within the framework of transport theory, this conclusion can be rationalised by estimating the relevant Péclet numbers that measure the relative importance of advection and diffusion. Using the average value Ū ∼ 1.3 μm/s of the mean flow speeds over time and edges as a velocity scale, we may estimate a mean Péclet number as Pe ∼ Ū R/D. Using R = 30 nm and the measured diffusivity D ≈ 0.6 μm²s⁻¹ [7], this leads to the estimate Pe ∼ 0.07. Flows affect transport for Péclet numbers above order-one values [13], and therefore the pinching-driven flows have a negligible influence on the transport inside the ER network. In order for fluid motion to have a noticeable effect, one would thus need flows to be either stronger, or sustained in one direction for longer. The chaotic flows produced by the pinching events with stochastic parameters suggested by resolution-limited microscopy [7] appear too weak to generate fast, non-diffusive edge traversals.

Importantly, our conclusions remain unchanged upon relaxing the no-slip boundary conditions, an important point to consider since the tubule membranes themselves could be moving in response to the nanoscale fluid flow. This can be modelled by introducing a finite slip boundary condition to couple the membrane-bound lipids with luminal flows (see Methods §IV G for details). Our simulations show that slip boundary conditions have virtually no effect on the average edge traversal speeds (Supplementary Fig. 12a). Physically, while a non-zero wall slip does change the shape of the flow profiles in the tubules (Supplementary Fig. 12b), it does not modify the mean flow speeds which directly affect global particle transport. These results justify our use of Poiseuille hydrodynamics (no-slip boundary conditions) throughout this work.

B. Tubule pinching-induced transport is network geometry independent and fails to facilitate luminal homogenisation

To estimate the network geometry’s contribution to transport simulation outcomes, we compared results across four different ER structures (C1-C4, Fig. S8). The distributions of the average edge traversal speeds appeared insensitive to ER structure variations for both pinching-induced and exclusively diffusional transport. This is reflected in the small deviations from the mean of the data points averaged across the different structures (Fig. 2de). Further, pinching-induced flows inside a regular honeycomb network (Fig. 2f, inset) with a typical ER edge length (1 μm) appears within a reasonable variance compared to the natural networks. Therefore, this mathematical idealisation of the ER network geometry can be used for exploring the consequences of the network ultrastructure contractility on transport kinetics in a standardised manner.

To test the effectiveness of the particle velocities for facilitating luminal material exchange across the ER, we tracked homogenisation kinetics by measuring intermixing of particles of two distinct colours equally seeded in each half of a honeycomb network at t = 0 (Fig. 3; see also Supplementary Video S2 for the flows inside such a network driven by pinches with the parameters from [7]). The measure of homogenisation is given by the variance Var(ϕ(t)) ≡ Var(n_b(t) − n_r(t)) over twenty regions of the network (Fig. 3a, horizontal lines) of the difference between the numbers n_b (blue) and n_r (red) of particles of each colour in region; note Var = 0 represents perfect homogeneity. The measures of mixing over time for pure diffusion and active transport with parameters [7] from Fig. 2 show a nearly complete overlap (Fig. 3e; see also Supplementary Videos S3 and S4), distinct from faster mixing under stronger flows in a network driven by pinches whose lengths are increased to their maximum possible value (i.e. the length of the tubule) and which are in addition sped up by a factor of 10 (see Supplementary Video S5; see also §II E for a discussion of the effects of artificially strong pinch parameters on average edge traversal speeds). This indicates that even if the experimental particle speeds were overestimated, the presumed pinching parameters would be inadequate to facilitate ER luminal material exchange.

FIG. 3:

C. Theoretical analysis of advection due to a single pinch explains weak pinchinginduced transport

The slow luminal transport driven by the pinching-induced flows is intrinsically linked to the volume of fluid expelled by each pinch during a contraction. The fundamental reason underlying the weak pinching-induced transport is that individual pinches are very weak generator of flow; even in the best possible configuration, the volumes of fluid pushed by each pinch are too small to impact luminal transport. Specifically, in Methods (§IV J 1), we mathematically show that, given a pinch of length 2L, the maximum distance Δz_max a suspended particle may be advected instantaneously by an individual pinch is

Using the experimentally-measured average value of the pinch length 2L = 0.14 μm [7], a typical pinch can then propel a particle by a maximum distance of Δz_max ≈ 0.19 μm.

This may, equivalently, be framed in terms of velocities. A transported particle experiences an average velocity during the contraction of at most V_max = 8L/3T, where T is the duration of a contraction or a relaxation. Using the pinch length as above and the average values of 2T = 0.213 s [7] this leads to the estimate V_max ≈ 3.9 μm/s, which is an order of magnitude smaller than the measured typical edge traversal speed of ∼ 45 μm/s, consistent with the order of magnitude difference between the measurements and the predictions of the computational model.

Note that this theoretical argument relies solely on the magnitude of fluid volumes expelled out of (or driven into) the pinching regions. These volumes are unaffected by the boundary conditions at the tubule walls because the presence of a non-zero slip length only changes the shape of the flow profile and not the volume fluxes (see Supplementary Fig. 12); any potential slip is thus is inconsequential for averaged measures of particle transport.

D. Marginal addition of coordinated contractility to luminal transport

The theoretical upper bound in the previous section for the maximum luminal transport produceable by an individual pinch is realisable only in the hypothetical situation where the flow generated by the pinch is all directed to one end of the tubule i.e. when the other end is effectively blocked. However, content transport produced during the contraction of a single tubule would then be reversed when the tubule relaxes back to its original state, with a single pinch site expected to exhibit only reciprocal (i.e. time-reversible) motions. Any advection contributing to edge traversals must thus be dominated by non-reciprocal motions of multiple pinches resulting in net displacements of solute particles.

The simplest system capable of producing non-reciprocal motions consists of two pinches, and the optimal sequence of motions to maximise the resulting advective particle displacement is illustrated in Fig. 5. We show in Methods (§IV J 2) that this is indeed the optimal two-pinch coordination, which results in a time-averaged displacement equal to the upper bound derived in Eq. (1). Since this optimal sequence of motions involves one pinch site starting a pinch halfway through the pinching of the other site, it is reasonable to estimate its duration as 3T, and therefore an average particle speed of 8L/9T ≈ 1.3 μm/s. The low particle speed achievable by the optimal coordination between two pinches suggests that a very high level of coordination among multiple tens of pinches per tubule would be required to reproduce the measured edge traversals.

FIG. 4:

FIG. 5:

E. A combination of high frequency and pinch length is required to replicate experimental particle speeds

Since the magnitude of ER contractility suggested by microscopy [7] (pinch lengths and frequency) do not explain the measured speeds, we set out to explore different sets of parameters that may generate particle velocities closer to the experimental measurements. The currently achievable imaging spatiotemporal resolution limits the detection of tubular diameter contractility by microscopy. Therefore, it is reasonable to postulate that the relevant parameters may have been underestimated. We simulated ER transport varying individually or in combination the values of pinch duration 2T, waiting time T_wait between successive pinches on a tubule, and pinch length 2L.

We first decreased both the original [7] values of T and T_wait by the same factor of 1/α, where α ≥ 1, and simulated particle transport in the honeycomb network (Fig. 4a). In effect, this simply ‘fast-forwards’ the flows in the original active network by a factor of α, and Brownian particles of the original diffusivity are released into this sped-up flow. Instantaneous and edge traversal speeds exhibited corresponding increases when we increased the value of α (Fig. 4a). An extreme value of α = 100 produces an average edge traversal speed distribution that peaks at around 8 μm/s (Fig. 4a). Similar results are observed in the C0 network from a COS-7 cell (Fig. 4b). The longer tails of these distributions (compared to those from the honeycomb network) result from the variation in edge lengths in the real network, with shorter edges, across which edge traversals are correspondingly fast, contributing to the tails.

These results suggest that, in order to produce average edge traversal speeds of the same order as the experimental values, we would need an active network sped up by a unrealistic factors considerably greater than 100, probably on the order of α ∼ 𝒪(10³) or higher and corresponding to multiple thousands of pinches occurring on average per second on each tubule. Similarly, it took an extreme increase in pinch site size to span the entire length of an average tubule, only to yield an average speed of ∼ 10 μm/s (Fig. 4d).

Next, we attempted to obtain a better fit to experimental data by combining changes in both the pinches’ time and geometry parameters. In Fig. 4e-f we thus sped up these maximally long pinches by a factor of α = 5 and 10 respectively, which yielded speeds averaging around 30 μm/s and above 60 μm/s respectively. Notably, the tail of the speeds distribution appeared longer than that seen in the experiments.

F. Luminal transport kinetics derived from contractile ER tubular junctions and sheets

As shown above, establishing effective transport in a tubular constrictions-driven model based on realistic ER fluid dynamics required a set of questionable assumptions, compelling us to explore alternatives. Thus, we set out to investigate how ER luminal transport would be impacted by the contractility of its structural components with volumes larger than tubules: (i) the tubular junctions (Fig. 1c), (ii) the perinuclear ER sheets (Fig. 1f) and (iii) the peripheral sheets (Fig. 1d).

First, we run numerical simulations of transport driven by contracting junctions on an ER network (C1 from Fig. S8d; sketch of junctions in Fig. 1c). Since junctions contracting at the tubular pinch temporal parameters measured experimentally yielded inadequate transport, we consider contractions/relaxations with duration 2T exponentially distributed with a mean of α⁻¹ times the original value in Ref. [7] and the waiting time T_wait between subsequent contractions/relaxations exponentially distributed with a mean β⁻¹ times the original value (i.e. values of α > 1 and β > 1 reflecting a higher frequency and faster pinches). Naturally, the results depend on the choice of the volume ΔV expelled by a junction during each contraction. Using fluorescence microscopy images, we may estimate the volume in junctions (see Methods, § IV I 1). We assumed that the junction volumes are drawn from a normal distribution with the same mean and standard deviation as our data set i.e. the distribution N(0.0045 μm³, 0.0021 μm³). Volumes greater than the maximum value in our experimental estimate (0.0081 μm³) or less than the minimum value (0.0020 μm³) are rejected.

We show the results in Fig. 6a. The thick solid line illustrates a set of values of (α, β) which produces distributions of average edge traversal speeds with means similar to the experimental values in Ref. [7], and the dashed lines one standard deviation away. We are interested in values of (α, β) close to unity because this corresponds to junctions which pinch with similar pinch durations and frequencies as the experimentally observed tubule pinches and are therefore biologically plausible; (α, β) ≈ (2.5,2) is the closest pair on this line to unity. However, only for α ≈ 1 we obtained average edge traversal speed distributions of reasonable shapes (i.e. approximately Gaussian), which would require a large β = 5 to match experimental results quantitatively.

FIG. 6:

Further, the assumption that the entire junction volume is expelled during a contraction may not be realistic, as it would require extreme bending/extension of tubule walls and the ability of the entire junction to empty and fill out during each contraction. Reducing the volume ΔV expelled in each contraction to two-thirds of the estimated distribution of the junction’s volume (results shown in Fig. 6b) or to half (Fig. 6c) causes the average edge traversal speeds to drop considerably. The lower the proportion of the junction volume expelled during a contraction, the faster the pinches are required to be (i.e. large value of α) as well as the larger the frequencies of the pinch events (large β) in order to produce reasonably high average edge traversal speeds.

Next, we considered the particle transport generated by two types of ER sheets contracting and relaxing over a duration 2T (see Methods for details). The perinuclear sheets are shaped as contiguous layers of flat cisternae with a luminal volume larger than the tubules branching from these structures (see sketch in Fig. 1f). Accordingly, their contraction with 2T = 5 s and V_sheet = 10 μm³ yielded a mean average traversal speed of 35 μm/s and consistent with the single particle tracking experiments (see Fig. 7b). However, the speeds distributions tailed towards higher values (Fig. 7a-b), something that has not be observed experimentally so far; it may be that high velocities cannot be recovered in experiments due to the constraints on linkage distance imposed to ensure trajectory fidelity in particle tracking.

FIG. 7:

Further, our simulations reveal that the flow decays sharply with distance away from the sheet where it originated as it branches out into tubules. This is illustrated in Fig. 7c-d and we further quantify the spatial gradient of the average edge traversal speeds across Cartesian 2D coordinates in Fig. 7e, revealing the stark contrast between the homogeneous profile for contracting junctions (red dotted line) vs sheet-driven transport (blue solid line). The short range of influence afforded by contracting perinuclear ER sheets thus argues against its ability to sustain mixing flows and fast particle transport on the distal tubular network.

Instead, we sought to explore whether the the peripheral sheets (i.e. the smaller flat intertubular ER regions, see sketch in Fig. 1d) may overcome the range limit. The peripheral sheets have volume significantly larger than tubules and junctions, which we estimated at 0.12μm³ ± 0.04 μm³ (see Fig. S14 and Methods §IV I 4 for details). This suggests that the peripheral sheets could produce average edge traversal speeds compatible with experimental measurements. Assuming that a sheet typically occupies the area enclosed by a ‘triangle’ of tubules, we may incorporate the transport inside a sheet-driven network using our model for junctions-contraction, but with junction volumes set to the measured sheet volumes (see §IV I 4 for details) and with either (i) each node expelling one third of the volume expelled by a contracting sheet (since peripheral sheet are in contact with three nodes on average), or (ii) a contracting node expelling the entire volume expelled by a contracting sheet, but with only one-third of the junctions actively contracting at any one time. Further assuming that during each contraction an peripheral sheet expels half its total volume so that in simulation (i) node k expels a volume V_k/6 in each contraction and in simulation (ii) each active node k expels a volume V_k/3, nodes which contract at rates α⁻¹ = 2.5 and α⁻¹ = 5 times slower than the tubule pinches in Ref. [7] in simulations (i) and (ii) respectively give average edge traversal speeds in the correct range of 40 μm/s (Fig. 7f-g). The contraction of peripheral sheets may thus be offered as a plausible mechanism for fast luminal transport, provided they are able to contract with sufficiently high rates.

III. Discussion

The motion of solutes in cellular compartments is now understood to be facilitated by active components. This is evident from direct motion measurements and the dependence of motion speed on the availability of ATP-contained energy [6–8, 14]. The origin of these active driving forces is however challenging to identify. The cytoplasm’s currents are often believed to originate from the motion of large complexes such as ribosomes and large vesicle cargo shuttled by the cytoskeleton-mediated motorised transport [14]. In the case of enhanced transport in the lumen of the ER tubules, the contractility of tubules has been suggested as the flow generating mechanism, and indeed such tubule deformations have been observed in microscopy [7]. However, establishing the direct empirical link between tubule contraction and active flows, or experimentally testing other hypotheses for the driving mechanism behind ER solute transport, remain currently unattainable. In this study, we thus propose a physical modelling approach, which provides a platform to explore the nanofluidics behaviour of biological systems such as the ER network. The outcomes of our simulations for a contractile ER argue against the plausibility of local pinch-driven flow; pinches with frequency and size on the order of those estimated by microscopy yield significantly lower speeds than single particle tracking measurements, as well as no enhancement of mixing beyond that from passive diffusion. The deficit stems from the fact that the displaced fluid volume upon local contraction is too small to generate sufficient particle transport.

Given the uncertainty of the empirical measurements for ER tubule deformation, due to the limits in the spatio-temporal resolution of organelle structure imaging, the pinch parameters may have been significantly underestimated. This sanctioned exploring a wider range of spatio-temporal parameters in our simulations, which revealed that a combination of a higher frequency with a much larger pinch length may provide higher particle speeds that are comparable to the single particle tracking measurements. Further, our modelling results suggest the possibility of transport by luminal width contractions of larger volume ER subdomains (which are contiguously interconnected with the network). In that respect, the contraction of peripheral (i.e. inter-tubular) sheets, in particular yielded speed values in a plausible range provided they are able to contrast fast enough. In contrast, the contraction of a large-volume perinuclear sheets led to fast transport but with a limited spatial range and with flows that would not impact transport beyond a few microns into the peripheral network. Similarly, the alternative scenario of tubular junctions’ contractility appears unlikely as particle speeds similar to experiments could only appears for junctions pinching at up to 8 times faster than experimentally observed contracting tubules.

The modelling approach in this study, although focused on the ER network, provides a step forward towards understanding intra-organellar fluid dynamics. Our simulation results have been used to rule out several scenarios, which seemed physically intuitive but are nevertheless unable to explain the observed enhanced luminal transport. Furthermore, our results generated a set of potentially testable predictions that can be used to validate or refute each of the envisaged transport mechanisms. For example, as fluid expulsion from peripheral sheets appears to be in broad in agreement with single particle tracking measurements data under even conservative assumptions, future measurements may explore whether an active luminal motion is more pronounced and faster in proximity to peripheral sheets. Moreover, a set of improved spatio-temporal resolution measurements of particle tracking and structural contractions will be needed to complete the physical picture of ER luminal transport.

While our study allowed us to test different plausible scenarios, questions remain open as to the force-generation mechanisms responsible for the observed pinching dynamics. We have estimated in Sec. IV G that forces on the order of 30 pN would be required to periodically bend the tubules as seen empirically. The force exerted by a single molecular motor has been estimated to be on the order of 6 pN [15], so the required forces for pinching may be provided by several motors working together [16]. Recently discovered hydro-osmotic instabilities could also contribute to shape fluctuations and pinching, although they are predicted to have much longer wavelengths than the typical size of a pinching region [17]. In contrast, the mechanisms involving topological remodelling of the ER, such as the welldocumented process of ring closure [18], occur over timescales of minutes and therefore cannot account for the millisecond-scale transport measured here.

Although we have used in silico fluidics modeling for our conclusions, particularly in identifying new physically permitted mechanisms of luminal propulsion, our simulations should be regarded as theoretical predictions demonstrating physical plausibility rather than mere speculation. We have explored theoretically the consequences of structural fluctuations which presumably do take place. Fluctuations in width of the flat ER areas, in particular, are expected since the structures appear dynamic in live light microscopy and narrowing/extension points are observable in electron microscopy [19, 20]. Further, variability in the areas of flat ER that are observable in electron microscopy can be explained by capturing the structures in different states of fluctuations. It should be expected that the mobile elastic structure such as ER sheets would not exhibit rigidity required to prohibit fluctuations.

In conclusion, it is worth emphasising that our in-silico fluid dynamical modelling reveals that for structural fluctuation-based mechanisms to facilitate luminal motion require assumptions currently not supported by empirical data. This warrants explorations of alternative mechanisms for ER luminal transport, for example anomalous diffusion driven by the fluctuations of macromolecular complexes [14] or by osmotic forces, as previously suggested [7].

IV. Methods

Here we describe the mathematical and physical model for the fluid mechanics and transport driven by active contractions in the ER. This requires the introduction of a network model for the geometry of the ER (§IV A) and individual pinches (§IV B) and a framework for the hydrodynamics of pinching tubules (§IV C). Our solution method for the flows inside our network is then described in §IV D. From simulations of Brownian particles advected by these flows (§IV E) quantitative measures of particle transport (§IV F) are extracted for comparison with experiment. Modifications of our model to explore alternative flow generation mechanisms are discussed in §IV I. Finally in §IV J we present in detail the derivations of the theoretical results discussed in §II C (transport upper bound by a single pinch) and §II D (coordination of pinches).

A. Network modelling

We represent the ‘skeleton’ of a two-dimensional ER network as a planar graph with each node assigned a position x ∈ R². Given an edge of the network labelled (i, j) and of length |x_i − x_j| = l, we model the lumen of the corresponding tubule to occupy a cylinder of radius R whose axis lies along the edge and has length l. This assumption avoids the intrinsic difficulty in defining a precise boundary between a tubule and a tubular junction, as well as leading to a simplified model of the intra-nodal dynamics of a solute particle (see below); since the tubules are long compared to the size of the junctions, the impact of their small overlap can be safely neglected.

We mathematically reconstruct a model ER network, which we refer to as C0, using the skeleton image of the COS7 ER network given in Fig. 4 of the Supplementary Material from Ref. [7]; this network is reproduced in Fig. 8a. We use the multi-point tool in ImageJ [21] to place numbered points at the positions of the nodes on the source skeleton image, and then obtain a list of the indices and position coordinates of each node. The edges are then manually tabulated as a list of pairs of nodes; this gives us all the information required to construct a mathematical graph as shown in Fig. 8b, with the original network superimposed on the mathematical model in Fig. 8c. Note that in what follows we work with the largest connected component of the source skeleton image in order to study transport in a fully connected network. We also use the same procedure to extract the graph structures from microscopy images, of four smaller ER networks which we label C1-C4 (original network with mathematical graph superimposed in Fig. 8d). In Supplementary Fig. 10, we show the distributions of the edge (tubule) lengths in each of the C0-C4 networks, as well as the mean edge lengths and the mean degrees (the degree of a node is the number of edges connected to it). The mean edge lengths are around 1 μm and mean degrees are approximately 3. In order to compare the biological network to an idealised ER system, we also consider below a honeycomb network, i.e. one where every node (apart from those at the boundaries) has a degree of 3, with all edge lengths exactly 1 μm.

FIG. 8:

FIG. 9:

FIG. 10:

B. Pinch modelling

To describe the pinches (Fig. 1e), we consider a model where the kinematics of each pinch is fully prescribed in time. We therefore assume that the active biochemical forces responsible for the deformation of the tubules balance with the elastic resistance of the tubules and with the dissipative forces in the fluid in such a way that the pinches occur as described.

The geometrical model for a pinch is illustrated in Fig. 9. Each tubule is assumed to have a pinch site at its midpoint (i.e. with L₁ = L₂); we make this simplifying assumption in our simulations, after verifying that more general pinch locations (L₁ /= L₂) have virtually no effect on particle transport. The pinching events occur at the pinch sites stochastically and independently of each other. Each pinch is defined geometrically by three parameters from a random distribution (see below): (i) the duration of a pinch 2T, (ii) the time T_wait between the end of a pinch and the beginning of a new one on the same site, and (iii) the length of a pinch 2L (see Fig. 9). We assume for simplicity that all pinches are axisymmetric so, using the notation of Fig. 9, each tubule remains a cylinder of time-varying radius denoted by r = a(z, t). Pinches are assumed to have reflectional symmetry about the plane in the crosssection through the centre of the pinch (i.e. Fig. 9 each pinch is characterised by the same length L on either side of it). The total length of a tubule is denoted by l = L₁ + 2L + L₂, so the portion before the pinch has length L₁, and that after the pinch is of length L₂.

We model the geometrical profile of each pinch of length 2L as following a linear radius change (see Fig. 9). Within a pinch located at z = z₀, the radius of the cylinder is therefore given by

The time-varying function b(t) ≤ R is therefore the minimum pinch radius in the centre of the tubule. As the simplest modelling choice, we assume in our simulations that b changes in time sinusoidally and thus will use b(t) = (R + b₀)/2 + (R − b₀) cos(πt/T)/2, where t is time after the pinch begins and b₀ is the value of b halfway in time through the pinch. Choosing a smooth time variation for the function b(t) will ensure the continuity of fluxes in time (see Eq. (10) later). We have verified that changing the pinch shape to a smoother geometry has essentially no impact on the fluid volume expelled/taken in, and therefore no significant effect on the flows/transport.

In our simulations the stochastic pinch parameters are drawn from the distributions measured in Ref. [7]. The pinch duration 2T is therefore drawn from an exponential distribution with rate parameter λ = ln 10/0.167 s⁻¹ (so the mean value is s), the time between pinches T_wait from from an exponential distribution with rate parameter λ_wait = ln 10/0.851 s⁻¹ (mean value s) and the pinch length 2L from a uniform distribution with mean μ = 0.12 μm and variance σ² = (0.06²/ ln 10) μm² ≈ (0.040 μm)². Throughout we set the tubule radius R = 30 nm and b₀ = 0.01R (recall from Fig. 9).

C. Hydrodynamic modelling

1. Hydrodynamics of network

We assume that the fluid occupying the ER network is Newtonian. Flow inside the network occurs at low Reynolds number, which can be justified in the following fashion. From experiments in Ref. [7], we know that velocity scale relevant to ER flows is of the order U ∼ 10⁻⁵ ms⁻¹. With a typical tubule radius R ≈ 10⁻⁸ m and a kinematic viscosity of at least that of water ν ≈ 10⁻⁶ m²s⁻¹, we obtain a Reynolds number of the order Re = UR/ν ∼ 10⁻⁷, so the flow is indeed Stokesian and inertial effects in the fluid can be safely neglected.

At a typical instance in time, a tubule undergoing contraction (or relaxation) causes a net volume flux to exit (or enter) the corresponding tubule. In the context of our graph theoretical model, we therefore model each pinch site as a ‘pinch node’ generating a net hydrodynamic source/sink whenever the tubule contracts/relaxes.

However, since it is not guaranteed that the total volume created by all pinching events always adds up to zero, we need a mechanism for the corresponding net volume to exit, or enter, the network. We achieve this through a number of ‘exit nodes’ that allow mass to be globally conserved. We locate the exit nodes at the exterior of the network and choose their number randomly. From a hydrodynamic standpoint, we impose the pressure condition p = p₀ at each exit node to model their connection to a large fluid reservoir.

We numerically tested the robustness of our results to the details of the exit nodes by repeating simulations with different configurations. The exact choices of exit nodes turn out to not affect the transport results shown below provided there are sufficiently many of them to avoid channeling the entire network’s worth of pinch-induced flow into a few tubules towards the exterior, thereby producing artificially strong flows.

2. Hydrodynamic model for a pinch

a. Flow rate

The velocity field in a straight cylindrical tubule at low Reynolds number is the classical parabolic Poiseuille flow [22]. When integrated through a cross section of the tubule, this flow yields the Hagen-Poiseuille law relating the pressure change Δp across a length l of a tubule to the a net volume flux Q in the positive axial direction,

where μ is the dynamic viscosity of the Newtonian fluid. Note that when Q > 0 our notation leads to Δp < 0, meaning that the pressure decreases across the length of the channel.

The result in Eq. (3) is valid for a straight (i.e. not pinching) tubule, and we need to generalise it to the case of a pinching tubule. Consider first a more general axisymmetric pipe whose radius a(z, t) varies with axial position z and time. We may use the longwavelength (lubrication) solution to Stokes’ equations for the streamwise velocity u(z, r, t) and flux Q(z, t) inside such a pipe into which a flux Q₁(t) enters at z = 0 [23],

Where

The equality in Eq. (5) can be derived using an intuitive mass conservation argument, independently of the inspired ansatz in Eq. (4). Conservation of mass inside a small section [z, z + δz] of the cylinder requires Q(z) − Q(z + δz) = 2πaa?δz + 𝒪 (δz²). Considering the limit δz → 0 and integrating the resulting expression for ∂Q/∂z from 0 to z yields Eq. (5).

In order for the no-slip boundary conditions at r = a(z, t) to be satisfied, and also to satisfy the incompressibility condition ∇ · u = 0, the radial component of the velocity, v(z, r, t), is necessarily non-zero and given by [23]

Our geometrical model for each pinch in §IV B yields straightforwardly a piecewise linear expression for a(z, t) in different regions of the tubule, with the time dependence entering only through the value of the pinch radius b(t). Denoting by Q_i (1 ≤ i ≤ 4) the fluxes in the four regions of the pinched tubules shown in Fig. 9, we may substitute the linear shape functions into Eq. (5) and obtain

Note that the final expression may be rearranged as

which may be interpreted as the instantaneous volume source/sink during a contraction/relaxation at a pinch site.

b. Pressure drop

We next need to compute the pressure drop in the pinches. We integrate the z-component of the Stokes equation

along 0 ≤ z ≤ L₁ + 2L + L₂, and use the solution for u, to obtain

Here, the second term on the right-hand side of Eq. (12a) has vanished because ∂_zu ∝∂_z(Q/a²) but Q and a are approximately constant at the entrance and exits of the tubule when the pinch site is sufficiently far from the ends of the tubule so that the flow is fullydeveloped there.

Using the expression for Q₂, an integration yields

and, by symmetry,

Subtracting these two results (and noting the integrals cancel out by symmetry) we obtain the modified Hagen-Poiseuille expression for a pinching tubule as

Note that this relationship is linear in each of Q₁ and Q₄, and is to be solved alongside Eq. (10) to relate the flow rates and pressure drops to the change in size of the pinches. Importantly, the classical Hagen-Poiseuille law is recovered as b → R since Eq. (15) becomes in that limit

which agrees with Eq. (3) when taking l = L₁ + 2L + L₂.

D. Solving the hydrodynamic network model

The incorporation of pinches as ‘dummy nodes’ into a graph theoretical framework of §IV A along with Eqs. (10) and (15) for the necessary pinch-related quantities allow us to reduce the problem of determining the time-dependent flows in an active pinching network into the simpler problem of solving at each instant for the instantaneous fluxes inside a ‘passive’ network with newly added nodes, appropriate sources/sinks, and modified pressure drops. Note that since the flows at these sub-cellular scales are inertialess (i.e. Stokes flows), we are able to effectively decouple time from our problem and solve the problem in the quasi-steady limit.

For each edge (i.e. tubule) (i, j) in the network, we define Q_ij to be the flow rate from node i to node j, with the sign convention such that flow is from i to j if Q_ij > 0. For mathematical convenience, we define Q_ij = 0 in all cases where (i, j) is not an edge in the graph. The goal is to solve for the values of the Q_ij’s corresponding to each edge.

After the incorporation of the dummy nodes, we denote by N the number of nodes and E the number of edges. We label the nodes such that {1, …, M } denotes the M exit nodes. Let q_i be the source or sink carried by the i^th node (so that q_i = 0 if i is a normal node, q_i is as specified by the RHS of Eq. (10) if i is a pinch node, and q_i is a quantity to be determined if i is an exit node.). Our M + E independent variables are therefore {q_i|i = 1, …, M } i.e. the sources/sinks carried by the exit nodes, and the Q_ij’s corresponding to each edge. To obtain their values, we employ the viscous hydraulic analogues of Kirchoff’s Laws.

1. Kirchhoff’s First Law (K1)

The first equation is that mass is conserved at each junction i.e. for each node i we have

which gives us therefore N equations. Note that these equations together imply global conservation of mass,

2. Kirchhoff’s Second Law (K2)

The second equation is a consistency of pressure, namely that the pressure change around any cycle (i.e. closed loop) of the network is zero. Therefore, in a given cycle C = {v₁, v₂, …, v_n, v_n+1 = v₁}, if denotes the pressure change from node v_i to node v_i+1 we necessarily have

Note that the pressure change across a node is negligible.

The K2 statement in Eq. (18) applies to all cycles in the graph, which would give us more equations than we need since the vectors of coefficients in Eq. (18) are linearly dependent. Instead, we need a minimal set of linearly independent K2 equations, corresponding to the cycles in a cycle basis of the graph, and we need only apply K2 to these cycles.

To construct a cycle basis of the graph G, we use standard results from graph theory [24]. We first construct a spanning tree T, defined as a connected subgraph which contains all the nodes of G and no cycle, as shown in Supplementary Fig. 11a on an example. Any tree T has N − 1 edges. Therefore there are E − (N − 1) edges in the graph G but not in the spanning tree T ; for each such edge e, we denote by C_e the unique cycle in the graph created by adding the edge e to the tree T (see Fig. 11b). The set of all such cycles C_e is then a cycle basis of G, and thus there are E − (N − 1) cycles in this set. There are therefore E − N + 1 independent cycles in the cycle basis.

FIG. 11:

To compute the spanning tree T we use a breadth-first search (BFS) algorithm [25]. We start from an arbitrary node and explore its neighbouring nodes. We add to T any previously unexplored node (and the corresponding edge) which does not result in the creation of a cycle in T . We then repeat this (in an arbitrary order) on the neighbours of the previous generation of nodes added to T, until no more nodes are left to explore. This is illustrated on an example in Supplementary Fig. 11c.

A similar algorithm is also used to compute the cycles C_e in the cycle basis. This time, however, denoting e = (i, j), the algorithm is started from i and set to terminate as soon as j is visited, yielding a path in the tree from i to j, which, together with the original edge (i, j), completes a cycle.

3. Pressure boundary conditions at exit nodes

At this point in the modelling, we have M + E independent variables (the flow rates in each edge and at the exit nodes), N equations from K1 (i.e. Eq. (17)), and E − N + 1 equations from K2 (i.e. Eq. (18)). The remaining M − 1 equations follow from the pressure conditions at the exit nodes. Specifically, the pressure difference between any pair of exit nodes i and j is zero

Here the sum Σ_i→j is defined to be over any path P_ij from i to j. Thanks to Kirchhoff’s second law, this quantity is path-independent. This indeed gives us equations, which we pick to be Σ _1→j Δp = 0 for j = 2, …, M .

We then solve the resulting linear system of M + E equations numerically. Note that we do not need to specify the value of the fluid viscosity μ in our algorithm because it cancels out in the K2 equation, Eq. (18), and in the pressure boundary condition equation, Eq. (19).

E. Simulating particle transport

With our solution for the flows in the active network at each instant of time, we now proceed to track the motion of Brownian particles inside the network using a discretisation of their stochastic equations of motion, as a model for the transport of proteins in the ER network.

We use the simplest approach where we superimpose Brownian motion onto advection by the flow inside each tubule. Let x(t) denote the position of a Brownian particle in a tubule and x_n the finite-difference approximation of x(nΔt), where Δt is a discrete time step. The displacement of the particle at each time step can be obtained approximately using an explicit first-order Euler scheme

where U is the instantaneous flow velocity, and the random noise term X(Δt) is drawn from a zero-mean Gaussian with variance (X(Δt)X(Δt)) = 2DΔtI [26], where D is the Brownian diffusivity of the particle. In our simulations, we take the diffusion constant to be the mean intranode diffusivity measured in Ref. [7], D ≈ 0.6 μm²s⁻¹. We include interactions between particles and walls by assuming that particles perfectly reflect off walls (i.e. elastic collisions). The particles are modelled as rigid spheres of diameter 5 nm, and the size of the particle matters only during elastic collisions with the walls of the tubules. As relevant in the limit of low volume fraction, we neglect hydrodynamic interactions between particles and perform ensemble averaging of the trajectories of many independent particles.

When a particle enters a node, we model its dynamics as follows. We consider a particle to have entered a node only if it has reached the end of a tubule, say of length l, at which instance we assign the particle to the node-point (i.e. the single point associated with the node in the graph description of the ER network). Although nodes contain a three-dimensional volume, their typical nodal length scale is of the order R ≪ l, and thus approximating them by point nodes is appropriate on the scale of the whole network. To decide towards which of the connected tubules the particle leaves the node, we estimate the values of the Péclet number Pe in each of the tubules. We define a local Péclet Pe_i = U_iR/D where U_i is the mean flow velocity through tubule i, with U_i > 0 for flow out of the node and U_i ≤ 0 otherwise. We then assume that the particle enters a neighbouring tubule i connected to the node with a probability proportional to max(Pe_i + 1, 0). This ensures that we have the expected behaviour in both limits of Pe: at high (positive) Péclet numbers, the probability is proportional to the flow speeds in each of the connected tubules, while at low Péclet the exit of the node is limited by diffusion and thus the exit is equally likely in each tubule.

F. Data processing: instantaneous speeds and average edge traversal speeds

During each simulation, we compute the edge traversal speeds as follows. A particle is defined to traverse an edge (i, j) if it travels from node (i.e. junction) i to node j, or from j to i, without visiting i or j in between. The corresponding edge traversal time is then the time between the arrival at the target node and the most recent departure from the node of origin. The edge traversal speed is naturally defined as the length of the tubule (i, j) divided by the edge traversal time.

The average edge traversal speed associated with an edge (i, j) is then defined as the mean over all edge traversal events across (i, j) of the edge traversal speeds.

In addition we also compute for each particle the ‘instantaneous’ speeds defined by V_n = |X(t_n+1) − X(t_n)|/Δt, where t_n = nΔt with Δt = 18 ms, which is the same temporal resolution as in the particle tracking carried out in Ref. [7].

G. Incorporating slip boundary conditions

The methodology we have detailed thus far assumes no-slip boundary conditions at the tubule walls for the fluid flow. However, the membrane-bound lipids themselves could also flow in response to the nanoscale luminal flows. This may be modelled by introducing a finite slip boundary condition on the tubule wall. The slip boundary conditions with a slip length λ ≥ 0 at the wall r = a(z, t) are given by

and n = 0, with , which simplifies to

We may then derive the long-wavelength solution for the flow field inside an axisymmetric deforming tubule as follows. Using an ansatz for the axial component u that is motivated by the uniform-radius Poiseuille flow with slip,

we may solve the incompressibility condition . Regularity at r = 0 constrains the integration constant to be zero, yielding

which automatically satisfies the boundary condition in Eq. (22). Note that the mass conservation equations are not affected by the introduction of a slip length.

Using the new solution for u, the modified Hagen-Poiseuille expression with slip may be derived as before to be

which does recover the no-slip result as λ → 0.

These results may then be used to simulate particle transport with slip boundary conditions. In Fig. 12a we plot the distributions of average edge traveral speeds obtained from simulations of a C1 network pinching with the original pinch parameters from Ref. [7] for four different values of the boundary slip length. In Fig. 12b we further display the profiles of the longitudinal flow (Eq. (23)) for different slip lengths with the volume flux fixed.

FIG. 12:

H. Estimate of forces required for pinches

In this section we derive an order-of-magnitude estimate for the forces required to pinch an ER tubule.

We first estimate the difference in the membrane’s elastic energy, ΔE = E_pinched − E_unpinched, between the pinched (E_pinched) and unpinched configurations (E_unpinched). In the absence of spontaneous curvature, the Helfrich free energy h per unit area of a membrane is given by

where k_c and are bending rigidities, H is the mean curvature, and K is the Gaussian curvature [27]. The mean and Gaussian curvatures may be expressed in terms of the principle curvatures κ₁, κ₂ as H = (κ₁ + κ₂)/2 and K = κ₁κ₂. We take κ₁ and κ₂ to be the principle curvatures in the directions normal and parallel respectively to the tubule’s longitudinal axis.

The dominant contribution to ΔE is from E_pinched, specifically from the region near the pinch site i.e. the centre of the pinch, where the tubule radius is smallest, taken to be at z = 0. We have κ₁ = 1/a(z, t) and κ₂ = 𝒪 (R/L²), yielding H = 1/2a + 𝒪 (R/L²) and K = 𝒪 (R/aL²). The dominant contribution to the membrane energy is therefore from H, in a neighbourhood of z = 0. Since real pinches do not have kinks, a linear term, for the purposes of estimating membrane energy, is unphysical, and we expect a(z, t₀) = b₀ + 𝒪 (z²) near z = 0, where b₀ is the pinch radius in the centre of the pinch in the maximally contracted state and t₀ denotes a time at which the tubule is in a pinched state. We may therefore write a(z) = b₀ + Rz²/L² to obtain an order-of-magnitude estimate of the membrane energy in the pinched state.

Evaluating the integral for E_pinched then yields the leading-order estimate

We may take k_c, (k_B is the Boltzmann constant times and T the room temperature) [28], with estimate values R = 30 nm, and L = 70 nm. The minimum pinch radius b₀ may be estimated to be 10 nm allowing for membrane thickness and incomplete squeeze, yielding

The hydrodynamic dissipation during a pinch may be calculated using the flow field in a pinch (see §IV C for intraluminal flow field), and it turns out to be much smaller than the elastic contribution to the work done to pinch a tubule.

An estimate for the force required to pinch the tubule may then finally be obtained as F ∼ ΔE/R ∼ 30 pN.

I. Modelling of alternative flow generation mechanisms

We have described in detail our modelling of a network driven by the pinching of tubules. In this study we explore two other flow generation mechanisms: the contraction of tubular junctions and the contraction of ER sheets. Our model for pinching tubules is readily generalised to account for these mechanisms, as we describe below.

1. Experimental estimates of junction volumes

From fluorescence microscopy images of ER networks, we measure the fluorescence intensity of junctions (i.e. the number of pixels in a junction multiplied by the mean intensity per pixel). We can then translate this intensity into an estimate of junction volume, assuming that intensity is directly proportional to volume. In order to calibrate the fluorescence intensity, we use our measurements for tubules. Specifically, we use the measured intensities of tubules and their known volumes, obtained by measurements of tubule lengths and the assumption that they are cylinders of radius 30 nm, in order to determine the proportionality constant to convert between pixel intensity and volume. Our new measurement of the intensities of the junctions then allows us to obtain estimates of their volumes; we obtain twelve values, ranging between 0.0020 μm³ and 0.0081 μm³, with a mean of 0.045 μm³ and a standard deviation of 0.0021 μm³.

2. Mathematical modelling of contractions of tubular junctions

To include the contribution of tubular junctions into the model (Fig. 1c), we assume that in addition to the pinch sites along tubules, each tubular junction pinches independently of other tubules and other tubular junctions. Given a junction, we assume that it expels the same volume ΔV of fluid during a contraction for all its pinches (and takes in the same volume when it relaxes). Each pinch is assumed to create a sinusoidal flow source, so a pinch lasting for a time 2T produces a flow rate S(t) = ΔV sin(πt/T)π/2T at a time t measured from the beginning of the pinch, where the numerical factors in S(t) are chosen such that the volume expelled during a contraction is indeed . We accommodate the flows from the junctions mathematically by modifying our K1 equations, Eq. (17), to allow the normal nodes to also carry non-zero sources (as opposed to just the pinch nodes, as was the case before); the other equations in the model remain unchanged.

3. Mathematical modelling of contractions of perinuclear sheets

In the tubules-pinching model, we included M exit nodes located towards the exterior of the network and through which flow could enter and exit the system in order to conserve mass. To account for the connection to a perinuclear sheet (Fig. 1f), we now assign randomly a number of these exit nodes, denoted M₂ < M, to be ‘sheet nodes’ i.e. nodes which are directly connected to a peripheral sheet, so that a number M₁ = M − M₂ > 0 of exit nodes remain. This is illustrated in Supplementary Fig. 13, where we show the C1 network from Fig. 8 with both sheet nodes (blue asterisks) and exit nodes (red squares).

FIG. 13:

A sheet contraction+relaxation lasting a time 2T produces a total source S_sheet(t) = V_sheetπ sin(πt/T)/2T at a time t from the beginning of the pinch, where again, the integral of S_sheet over a contraction gives a volume V_sheet of fluid expelled by the sheet.

Similarly to the mathematical model with tubular pinches only, our independent variables are the E tubule fluxes, the M₁ sources at the exit nodes, and the M₂ sources at the sheet nodes. As before, the K1 (Eq. 17) and K2 (Eq. 18) equations give us E + 1 equations. The requirement that the exit nodes are at the same (exit) pressure give us equations, namely the equality ΔP_ij = 0 for i, j exit nodes. We make the additional assumption that the sheet nodes are all at the same mechanical pressure (i.e. that they are connected to the same reservoir), which provides an additional M₂ − 1 equations. Requiring that the sources at the sheet nodes sum to the prescribed flow rate S_sheet yields one additional equation, so we again have a system of E + M independent linear equations.

4. Experimental estimates of volumes of peripheral sheets

To model the contraction of peripheral sheets (Fig. 1d), we need to estimate their contained volumes. Using the open-source image analysis software Fiji [29], we identify nine regions roughly occupied by peripheral sheets in a microscopy image of an ER network; this is illustrated in yellow in Supplementary Fig. 14. We then measure their areas (in μm²) and convert them to volumes by multiplication with the diameter of a tubule; we take the diameter to be 60 nm i.e. we assume that the effective thickness of a sheet is equal to the tubular diameter. From our nine data values, we finally obtain a mean sheet volume of 0.12 μm³ and a standard deviation of 0.04 μm³.

FIG. 14:

J. Derivation of theoretical bounds for active flows driven by pinching tubules

1. Advection due to a single pinch

In this section we calculate an upper bound for the axial distance Δz a particle can be advected by the flow produced by an individual pinch. Recall the formula for the volume ‘source’ due to a pinch,

Among all possible ways for a particle to be transported by a pinch, an upper bound on the transport distance Δz ≤ Δz_max can be reached if all of the following conditions are satisfied: (i) All of the source flows to one side of the pinch (i.e. there is no leakage on the other side); (ii) The particle travels outside the pinching region (axial flows within the pinching region produce smaller advective displacements than those outside, as may be verified numerically); (iii) The particle travels along the centreline of the tubule (i.e. at twice the cross-section averaged flow velocity, a standard Poiseuille result); (iv) The minimum pinch radius b₀ (recall from Fig. 9) is 0.

Under conditions (i) and (ii), the cross-section averaged speed corresponding to maximal transport generated from a single pinch can then be computed using Eq. (29) as

Using condition (iii), the position of the particle along the centreline z(t) satisfies then the ordinary differential equation

Integrating this equation from t = 0 (start of contraction, with b(0) = R) to T (end of contraction, with b(T) = b₀) then leads to

The value of Δz is maximized when b₀ = 0 (i.e. when condition (iv) holds), yielding the upper bound

2. Extension to nonlinear interactions between two pinches

An isolated pinch is only capable of reciprocal motions. The simplest system capable of producing non-reciprocal motions is illustrated in Fig. 5 and consists of two pinch sites arranged in series near the midpoint of a long horizontal tubule.

We may calculate an upper bound on the net particle displacement that can be achieved after both pinch sites pinch exactly once. We make the assumptions (ii)-(iv) as above, and in addition, allow pinches to spend extended amounts of time in their completely closed state; any deviation from these assumptions will result in net transport that is further reduced.

How can we maximise the positive (i.e. rightward in Fig. 5) particle displacements induced by the contractions, and minimise the magnitude of the negative (leftward) displacement induced by relaxations? As in the figure, let us denote the pinch on the left ‘pinch 1’ and that on the right ‘pinch 2’. Pinch 1 produces maximum positive displacement (magnitude 4L/3, i.e. half of the optimal from Eq. (33)) when pinch 2 is completely open i.e. when the hydrodynamic resistance to its right is minimal. Pinch 2 then produces maximum positive displacement (of magnitude 8L/3, i.e. the optimal value in Eq. (33)) when pinch 1 is completely closed and all pinch-induced flow continues to be directed rightwards. Similarly, pinch 1 produces minimal negative displacement when pinch 2 is completely closed (zero average displacement, since no flow is then allowed to escape to the right of pinch 2), while pinch 2 then produces a negative displacement of magnitude 4L/3 to the left when pinch 1 is completely open. These optima can be achieved by the non-reciprocal sequence of motions illustrated in Fig. 5: close pinch 1, close pinch 2, open pinch 1, open pinch 2. The coordination between two pinches can therefore be used to generate the net displacement of 8L/3 equal to the theoretical upper bound from Eq. (33).

Supplementary figures

Acknowledgements

EL and PHH are supported by the European Research Council under the European Union’s Horizon 2020 research and innovation program (Grant No. 682754 to EL) and by the Cambridge Trust. EA is supported by the UK Dementia Research Institute [award number UK DRI-2004] which receives its funding from UK DRI Ltd, funded by the UK Medical Research Council, Alzheimer’s Society and Alzheimer’s Research UK