1. Computational and Systems Biology
  2. Plant Biology
Download icon

Rapid transporter regulation prevents substrate flow traffic jams in boron transport

  1. Naoyuki Sotta
  2. Susan Duncan
  3. Mayuki Tanaka
  4. Takafumi Sato
  5. Athanasius FM Marée  Is a corresponding author
  6. Toru Fujiwara  Is a corresponding author
  7. Verônica A Grieneisen  Is a corresponding author
  1. The University of Tokyo, Japan
  2. John Innes Centre, United Kingdom
Research Article
  • Cited 2
  • Views 1,397
  • Annotations
Cite as: eLife 2017;6:e27038 doi: 10.7554/eLife.27038


Nutrient uptake by roots often involves substrate-dependent regulated nutrient transporters. For robust uptake, the system requires a regulatory circuit within cells and a collective, coordinated behaviour across the tissue. A paradigm for such systems is boron uptake, known for its directional transport and homeostasis, as boron is essential for plant growth but toxic at high concentrations. In Arabidopsis thaliana, boron uptake occurs via diffusion facilitators (NIPs) and exporters (BORs), each presenting distinct polarity. Intriguingly, although boron soil concentrations are homogenous and stable, both transporters manifest strikingly swift boron-dependent regulation. Through mathematical modelling, we demonstrate that slower regulation of these transporters leads to physiologically detrimental oscillatory behaviour. Cells become periodically exposed to potentially cytotoxic boron levels, and nutrient throughput to the xylem becomes hampered. We conclude that, while maintaining homeostasis, swift transporter regulation within a polarised tissue context is critical to prevent intrinsic traffic-jam like behaviour of nutrient flow.


eLife digest

Every multicellular organism, including all plants and animals, faces the challenge of taking up the nutrients it needs and distributing them throughout its body. Plants absorb many nutrients including nitrogen and boron from the soil into their roots, often using tightly controlled processes that require energy to work. Plant roots contain several distinct layers of cells and the nutrients need to cross these layers to reach a channel at the centre of the root known as the xylem, which transports the nutrients to other parts of the plant.

Plants need boron to grow. However, high levels of this nutrient are toxic so plants have evolved to change the rate at which they absorb boron to optimize growth in different environments. When there is little boron in the soil, certain transporter proteins move to the surface of root cells to bring boron into the root more effectively. On the other hand, when plants grow in soils with high boron, their root cells have fewer of these transporters on their surfaces to prevent too much boron entering the plant.

This regulation of boron uptake appears logical, except for one detail: at any given location, the amount of boron in the soil is relatively stable and changes only very slowly. Why do plants invest energy in responding rapidly to the supply of a nutrient that changes so slowly in nature?

Sotta et al. used mathematics and experimental approaches to study boron uptake in a plant known as Arabidopsis. The work reveals that the plants ability to rapidly alter how efficiently boron moves into root cells actually serves to avoid internal “traffic jams” in boron transport. If the numbers of transporter proteins on the surface of root cells changed more slowly, individual cells would occasionally experience high levels of boron that would interfere with the movement of boron further into the root, causing a jam. Furthermore, these ‘peaks’ of boron could damage the individual cells they affect.

The findings of Sotta et al. reveal that, by being able to rapidly change the numbers of certain transporter proteins on the surface of root cells, plants can ensure they receive a steady supply of boron. This work suggests that to develop artificial systems that can adapt to changing surroundings, researchers will need to engineer solutions like those found in plants in order to avoid similar traffic jams in the systems. Along with considering how plants interact with their environment, studying how they avoid internal traffic jams in nutrient uptake may help efforts to alter plants, including crops, so that they grow better in harsh environments.



Robust growth and functioning of all living organisms, including plants, requires a well-balanced uptake of nutrients from the environment. However, nutrients are not always readily available in the environment at the organisms’ optimal levels. This issue is exacerbated for sessile plants which take up essential minerals, required for their growth and functioning, from the soil in which the nutrient concentrations are mainly determined by elemental contents of igneous rocks. Most soils are either deficient or in excess of essential elements (Roy et al., 2006). Plants have therefore evolved to take up these essential minerals from the soil in a regulated manner over a wide range of concentrations through their roots, conveying the nutrients to the rest of the plant body via the vascular network. For a well-balanced uptake of nutrients from soils by roots, nutrient transport processes need to be regulated in a nutrient-dependent manner, which constitutes one of the essential processes for plant growth and crop production.

For roots to transport soil-derived nutrients into its central vascular tissues, which will ultimately transport it to other regions of the plant, a flow of nutrients has to be established which crosses several cell files, from the cells in contact with the medium (typically epidermal or lateral root cap cells) to the xylem. Boron is one of such essential nutrients for plants, critical for cell wall composition (O'Neill et al., 2004). Borate is cross-linking a pectic polysaccharide, Rhamnogalacturonan II, which makes it indispensable for tissue growth (Kobayashi et al., 1996). Only available in the soil, boron is transported to the xylem through an intricate system of polarly localised membrane transporters that bring boron from the apoplast into the cytosol, and then back again into the apoplast (Miwa and Fujiwara, 2010). The boron transport mainly takes place in the form of boric acid and borate. At neutral pH, boron is mostly present in the form of boric acid (Woods, 1996).

Two families of transporters, BORs and NIPs, are involved in boron transport through the root. BOR proteins export boron from within the cell to the cell wall, while NIP proteins enhance bidirectional permeability. Both BORs and NIPs present polar localisation, often complementary. In the root, BOR1 (At2g47160) and BOR2 (At3g62270) locate on the inwards-facing membranes of the cell (Takano et al., 2010; Miwa et al., 2013), while NIP5;1 (At4g10380) locates on the outwards-facing membranes (Takano et al., 2010). Being polarised in such a manner, their combined action and localisation allows for a highly efficient boron uptake and transport to the xylem, even when boron is only available at very low concentrations in the medium. This transport system ensures that within the root tissue much higher boron concentrations can be reached than what is available in the medium, allowing for significant xylem loading which then provides sufficient boron for the growing shoot. Indeed, both nip5;1 and bor1, bor2 mutants present severe growth defects under low boron concentrations (Noguchi et al., 1997; Takano et al., 2006-06; Miwa et al., 2013). Such a coordinated and directed polar transport is reminiscent of the intricate polar auxin transport networks within the plant, although their axes of polarity are often perpendicular to one another (Robert and Friml, 2009).

Likewise, the boron transport system faces constraints: While boron is required for cell wall strength and stability, too high intracellular boron concentrations elicit DNA damage, growth retardation and eventually cell death (Sakamoto et al., 2011-09). Unsurprisingly, therefore, is has been found that the protein levels of the boron transporters are regulated by boron itself, with both NIP5;1 and BOR1, 2 protein levels dropping at higher boron concentrations (Takano et al., 2005; Takano et al., 2006-06; Takano et al., 2010; Miwa et al., 2013). Such boron-dependent regulation is considered essential to allow for boron homoeostasis, preventing boron toxicity when boron availability is high, but allowing for efficient uptake when availability is low. Protein down-regulation takes place through two distinct mechanisms. In the case of NIP5;1 boron reduces protein levels via mRNA degradation (Tanaka et al., 2011; Tanaka et al., 2016), while in the case of BOR1, 2 the mechanism involves increased protein degradation (Takano et al., 2005; Miwa et al., 2013). Surprisingly, however, in both cases the down-regulation of boron transporters by boron occurs on a short time scale: when a plant is transferred from low to high levels of boron, swift downregulation of BORs (via protein degradation), and NIPs (via transcript degradation) is observed (Takano et al., 2005; Tanaka et al., 2011; Tanaka et al., 2016). The BOR1 degradation occurs through endocytosis, apparent 30 min after the transfer from low to high boron media. After two hours BOR1 has already mostly disappeared, suggesting that the half life of BOR1 is well below one hour (Takano et al., 2005). The half life of NIP5;1 mRNA after the transfer from low to high boron media is 10–15 min (Tanaka et al., 2011). Such rapid time scales seem at odds with the expected natural variations of boron a plant would experience, as there is no evidence supporting considerable fluctuations in soil boron concentrations, neither spatially nor temporally. This is due to boron being available to plants as boric acid, which is highly water-soluble (solubility: 0.92 mol/L at 25°C). Consequently, boron is very mobile in the soil (Nable et al., 1997), rendering patchy heterogeneous boron levels throughout the soil neither stable nor probable. The main process that presumably would allow a plant growing at a fixed location to experience rising boron levels is through drought, a phenomenon which fails to account for the necessity of transporter down-regulation occurring on the order of minutes. Also watering of plants is not expected to quickly change the boron levels. Keren and Bingham, 1985, analysing the effects of solution-to-soil ratio on the boron concentration in soil water, propose that soil adsorption plays a role in buffering the fluctuations in boron concentration in soil water as a consequence of fluctuations in the water-soil ratio. Mathematical modelling was used to predict the boron concentrations in the soil solution from the water-to-soil ratio (Keren, 1981), revealing that boron concentration in the soil is robust against fluctuations in the water-to-soil ratio if soil adsorption is considered. From these in-depth studies, the picture consolidates that rapid fluctuations in boron concentration are indeed unlikely, or rare events even if they could occur under very specific and unlikely conditions. Given that there is no apparent necessity for swift regulation, it is thus surprising to find cost-ineffective down-regulation mechanisms through degradation of mRNA and protein underpinning this system, instead of more cost-effective down-regulation processes via transcriptional repression.

In short, down-regulation of boron transporters can be readily understood as a natural adaptive mechanism for plants to optimise growth and function at different geographical locations with varying natural boron concentrations. Nevertheless, it remains intriguing as to why plants have evolved such a swiftness in the regulation of these transporters. Puzzled by this behaviour, we questioned if other dynamical constraints are operating on polarised root tissue that could explain the need for these rapid timescales. To address this, we explored possible implications of dynamical transporter regulation in a parsimonious model that captures the nutrient flux across a small cross-section of a generic root tissue, with an explicit focus on the swiftness of the temporal regulation.

A simple model for boron uptake to probe implications of boron transporter regulation swiftness

In Arabidopsis thaliana, boron is shuttled from the soil across the epidermis, cortex and endodermis into the stele, to finally be taken up into the xylem, mediated by the action of secondary active boron exporters (BORs) and the enhancement of permeability through NIPs. Previously, we analysed boron patterning by using a two-dimensional cross-sectional model of the entire Arabidopsis root meristem, considering all spatial nuances of transporter localisation, polarity and intensity, whilst neglecting transporter dynamics and regulation (Shimotohno et al., 2015-04). Instead, transporter levels and activity were static, fixed to the observed steady state transporter expression under constant boron medium condition of 0.3 μM. While neglecting regulation of transporter expression or activity, a characteristic boron pattern emerged on the level of the root, presenting highest concentrations at the stem cell niche. This boron profile gradually decreased longitudinally up to the start of the differentiated tissues. The computationally simulated characteristic boron distribution, which we confirmed experimentally through LA-ICP-MS, strongly indicates that the mature root tissue is involved mainly in transporting boron from the soil into the xylem, while the distal meristematic tissues have a differential boron transport function, namely, to locally provide higher boron levels for fuelling local growth (i.e., cell wall material) at those regions. In that study, however, we did not address the dynamics of the transporter regulation nor how nutrient homeostasis is achieved. Here, we investigate the dynamics of the transporters and their mutual feedback with the generated boron distribution. We do so by spatially focussing on the differentiated tissue region involved in xylem-uptake.

At the elongation and differentiation zone of the root, transporters possess a striking polar expression; BORs are located at the inside-facing membranes of the cells, while NIPs are concentrated to the outside facets (see Figure 1A,B). Although protein levels can present strong variation between cell files, the transcription of the transporters takes place throughout the entire root tissue, even for NIP5;1 (Figure 1—figure supplement 1). We therefore use the simplifying assumption that all cell files are intrinsically the same in respect to the potential of expressing the transporters, and variations arise solely as a consequence of the nutrient distributions.

Figure 1 with 2 supplements see all
A polarised tissue model for nutrient boron transport across roots.

(A) Confocal microscope images of BOR1-GFP and GFP-NIP5;1 localisation in A. thaliana roots, revealing polar localisation of both transporters, facing inwards and outwards, respectively. Cell file identities are shown in the top-left with numbers representing, in increasing order, cell files from the outermost to the innermost cells, such as ep, epidermis; co, cortex; en, endodermis; va, vasculature. (B) Detailed view along a transversal section of the tissue at the proximal meristem. lrc, lateral root cap. (C) Schematic diagram of the NIP-BOR boron transport model, consisting of a transversal root cross-section composed of n cells between the medium and the xylem. For a simplified A. thaliana model, we consider n=4. Depicted is a generic root model consisting of n cell files, not showing all cells and cell walls in between W2 and Wn2 for illustration purposes. The model includes intracellular and apoplastic compartments, as well as membrane-based properties such as transporter activity and background permeability rates. lw: cell wall width; lc: cell width; hc: cell height. For a full description of the parameters, see Table 1.


Our analysis focuses on the transversal nutrient flow through the root that results from boron entering and leaving the different cell files transversally. Effectively, in our simple model, we only consider a single row of cells and, for simplicity, only four cell files over which boron is transported from the soil into the xylem (Figure 1C). To capture the dynamics of boron transport and transporter activities in root tissue, the model’s variables are the boron concentration in cells (C) and cell walls (W), and transporter activities of NIPs (N) and BORs (B) for each individual cell (n). The mutual dependency between these variables is described using ordinary differential equations (ODEs). For each cell file of the model, the transporters are produced at an equal rate, have the same transport potential, and are regulated in the same way. The modelled cells loosely map to the outermost epidermis (cell 1), the cortex (cell 2), the endodermis (cell 3), and finally to the pericycle/vasculature (cell 4). The last cell of the row is endowed with an upward convective flow, which captures the shootwards convective flow of the xylem. For reasons of simplicity, we consider these four cells to have equal dimensions, and make the same assumption for all intermediate cell walls (Figure 1C).

Given that transporter regulation can be observed in response to varying the boron conditions in the medium (Takano et al., 2005, Takano et al., 2006-06Takano et al., 2006-06, 2010; Miwa et al., 2013), (Figure 1—figure supplement 2), boron sensing underpinning those responses could either involve measurements of the intracellular concentration or involve measurements of the intercellular concentration (i.e., the levels within the apoplast). In the case of NIP5;1, we recently revealed that NIP5;1 mRNA degradation is triggered by boron-dependent ribosome stalling during the translation process (Tanaka et al., 2016). The boron-dependent degradation was also observed in an in vitro system without cell wall fraction, which suggests that the sensing mechanism depends on the cytoplasmic boron concentration. For the BORs the precise subcellular location of the boron sensing has not yet been determined, but it is reasonable to assume that the regulation is also responding to cytosolic boron levels, given that the purpose of the nutrient homeostasis is to keep the nutrient levels within the cell within bounds, and to this end directly sensing the cytoplasmic boron concentration is beneficial. In short, given the current knowledge, cytosolic boron sensing is the most parsimonious assumption.

Hence, our model assumes that all the boron transporters of a cell are regulated by the cytosolic boron concentration within that respective cell. Note that we loosely use the term ‘boron’ in our model to refer to both boric acid as well as borate. They can be found in a dynamic chemical equilibrium, B(OH)3+H2OB(OH)4+H+, with a pK value of 9.24 (Woods, 1996). Given that this process does not involve any enzyme kinetics, this distinction is not explicitly considered in our model. Moreover, we solely focus on the import of soluble boron (both boric acid and borate), not considering any boron that is bound to the cell wall.

The temporal dynamics of the soluble boron concentration within the cytosol (Ci) and in the cell walls (Wi, being the cell wall adjacent to the inward facing facet of cell i) is determined by their inflow and outflow from and to neighbouring compartments:

(1) Ci˙={(p+Ni)(Wi-1-Ci)-BiCi+p(Wi-Ci)}1lc,
(2) Wi˙={p(CiWi)+(p+Ni+1)(Ci+1Wi)+BiCi}1lw.

Here, Ni and Bi represent boron permeabilities due to, respectively, the bi-directional transport by NIPs and the directional transport by BORs within cell i. Transporter-independent boron permeability through the plasma-membrane is also taken into account, described by parameter p. To capture the tissue context, boundary conditions are set at both extremities of the transversal cell series. The outer cell wall of the first cell (cell 1) is in contact with medium. Given the rapid diffusion rate of boric acid in water, we assume that the boron concentration in the outermost cell wall (W0) is the same as in the medium, constant over time. The xylem transport that occurs from roots to upper tissue is represented by a removal term, with rate a, attributed to the innermost cell (i=n) of the cell row:

(3) Cn˙=(p+Nn)(Wn-1-Cn)1lc-aCn1hc.

Assuming that transport activities are proportional to protein concentration, it follows that the dynamical regulation of NIP and BOR permeability (Ni and Bi) are in direct proportion to the dynamics of their proteins (we will later show that our main insights do not depend on the assumption that the BOR transporter does not saturate under the range of concentrations here treated). It was therefore not necessary to introduce explicit variables for protein concentrations. Instead we directly use permeability, representing the protein concentrations and their regulation. The protein concentration, and thus the effective permeability, is determined by production and degradation rates.

As mentioned, accumulation of NIP5;1 mRNA is regulated through boron-dependent mRNA degradation. This suggest that the production of NIP5;1 protein is boron dependent. On the other hand, degradation of NIP5;1 protein is not found to be boron dependent (Takano et al., 2010; Tanaka et al., 2016). Thus, the levels of Ni vary due to a production term that is dependent on cytosolic boron concentration and a constant degradation term:

(4) Ni˙=αNkNnNkNnN+CinN-ξNNi,

where αN is the production rate, ξN is degradation rate, and nN is the Hill coefficient capturing the boron-dependent inhibitory regulation.

Accordingly, given that BOR protein levels are regulated through boron-dependent degradation (Takano et al., 2005), but there is no evidence of boron dependent production, Bi varies over time due to constant production and boron-dependent degradation terms:

(5) Bi˙=αB-ξB(1+dBCinBkBnB+CinB)Bi,

where αB is the production rate, ξB the basal degradation rate, dB the maximum boron-independent degradation rate, and nB the Hill coefficient capturing the boron-dependent inhibitory regulation.

We wish to use this model to understand the significance, if any, of the swift transporter regulation. To assess the effects of regulation swiftness on transport, we therefore introduce a time-scaling factor ϵ, which multiplies both the entire Ni˙ and Bi˙ terms (Equations 6 and 7). To align ourselves with the current experimental evidence that indicates swift regulation of NIP5;1 on the mRNA level and of BORs on the protein level, we modulate the production rate and degradation rate for NIP5;1 and BOR1, respectively. However, changing production rates of NIPs through a factor ϵ without changing the degradation rate would result in changes of the overall protein concentration which might influence the system’s behaviour, not due to the swiftness of the regulation but due to the changes in equilibrium protein concentration. This also holds for the BOR regulation. Thus, to avoid such possible undesired effects, the time-scaling factor multiplies the entire right-hand side of the dynamics:

(6) Ni˙=ϵN(αNkNnNkNnN+CinN-ξNNi),
(7) Bi˙=ϵB{αB-ξB(1+dBCinBkBnB+CinB)Bi}.

We define ϵ=1 as representing the swiftness of regulation in the wild type scenario, while reducing ϵ allows us to simulate behaviours that would occur with less rapid regulatory dynamics, and increasing ϵ allows us to consider even higher swiftness of regulation than experimentally observed. Note that the regulation of the transporter activity by cytosolic boron concentrations occurs through a sigmoidal relationship, in which kB and kN are the boron concentrations at the half maximum of the Hill function, defining the sensitivity of the response to the boron concentration. We set kB30kN, consistent with our observations that the sensitivity of BOR expression is much less than that of NIP5:1 (Figure 1–Figure supplement 2), but will also show that our results do not require large differences between these two parameters.

We have assessed what is the biologically acceptable range for each parameter, based on the available literature (see Material and methods for details). Default parameter values were chosen to lie within these valid ranges, and are given in Table 1.

From steady state flows to oscillatory dynamics

When regulation swiftness lies within the experimentally observed regime (ϵB=ϵN=ϵ=1), a constant flow of boron that is transported across the root ensues (Figure 2). Boron concentrations within the cells (Cn) remain constant over time. This steady state behaviour, however, is disrupted when the transporter regulation swiftness is decreased ten-fold (ϵ=0.1). Oscillations in the boron concentration arise in all cells, but most strikingly in the endodermis (cell 3) and cortex (cell 2) (Figure 2A), as well as in the cell walls (Figure 2B). Decreasing the transporters’ regulation swiftness even further (ϵ=0.01, 0.001) consistently enhances the amplitude of such oscillations, and enlarges their periods (Figure 2A,B). We conclude that even when environmental boron levels (C0) and xylem activity (a) remain constant, simply deviating from the rapid experimentally observed transporter’s regulation swiftness is sufficient to push the boron transport system into an oscillatory dynamics.

Figure 2 with 1 supplement see all
Oscillatory behaviour arises due to decreased transporter regulation swiftness.

Time development in the four-cell NIP-BOR model of the boron concentration in the cells (A) and cell walls (B), and of the transporter activity of BOR (C) and NIP (D). The simulations are started using the default parameter setting (Table 1). At the time points 48, 72 and 96 hr, the parameters determining the transporter regulation dynamics (ϵN and ϵB) were each time reduced to one tenth of their previous value. C1 to C4 are the boron concentrations in the outermost cell 1 up to the innnermost cell 4. Wi represents boron concentration in cell wall fraction between cell i and cell i+1 and Ni are transporter activity of BOR and NIP in cell i, respectively.


The reason the cells undergo such dramatic concentration variations at lower regulation swiftness can be intuitively linked to the accompanying changes in transporter levels which also ensue (Figure 2C,D). They are initially triggered by minute fluctuations in the concentrations, but then cause increasingly larger changes in cytosolic and apoplastic concentrations. It is not immediately clear, however, how it depends on either the NIP or BOR transporter regulation, why regulation swiftness exacerbates oscillation amplitudes and how their interdependencies either produce stable flows or unstable oscillations that propagate throughout the tissue. What we can already conclude though, is that within a polarised tissue rapid transport dynamics are necessary to ensure flux homeostasis and avoid instabilities in concentration levels. But to fully understand the process we will first observe how in the model the steady flow regime behaves in terms of transporter expression levels when we assume that the regulatory dynamics are as swift as experimentally observed.

Parsimonious model, under wild type dynamical settings, generates stable transporter expression

Given that our question is focused on the possible advantages of swift transporter regulation within a polarised tissue context, our model purposely ignored details of tissue-specific regulation, using the simplifying assumption that all cells are endowed with the same potential of expressing transporters. This assumption may initially seem at odds with the actual distinct GFP levels of the GFP-tagged proteins (Figure 1B), which qualitatively indicates stronger BOR1 expression in the epidermis and endodermal files and pronounced NIP5;1 expression in the outermost cell file (lateral root cap or epidermis) only. It was therefore interesting to observe that our parsimonious model could already account for general features of the the observed transporters’ expression patterns (Figure 2C,D) solely as a consequence of the relative positioning of the cells within the cell row: in the model, NIP only becomes highly expressed in the outermost cell, the epidermis (N1, representing NIP5;1 in cell 1, is by far the highest variable in the steady flow case). These relative levels are in qualitative agreement with the experimental findings regarding NIP5;1. Moreover, in contrast to this strong expression in the epidermis, in the model NIPs are present at much lower amounts in all the other inner cell files (cell 2, cell 3 and cell 4), again qualitatively similar to the experimental observations. These model results can be understood as follows. Although all in silico cells share the same intrinsic properties, as a consequence of NIPs acting as permeability facilitators rather than directional transporters, it is impossible for the cell file in contact with the medium (cell 1, the outermost cell file, that is, epidermal cell) to reach higher concentrations than the level in the medium. Only the next cell inwards, cell 2, is capable of reaching higher levels, due to the directed polar action of the BORs. The inevitable low boron levels in the epidermis give rise to the observed lack of NIP downregulation within the first cell file, and hence to its preferential expression.

The expression patterns of the active exporters, BORs, are quite different from NIP5;1, but again there are similarities between the model results (B values) and confocal images, with BORs distinctly expressed in the epidermis (cell 1) at the highest levels, in the cortex (cell 2) at high levels, whilst in the endodermis it is relatively weaker (cell 3) (Figure 1B). The model’s qualitative correspondence between the transporter expression patterns and those experimentally observed suggests that concentration-dependent regulation is at least sufficient for such expression patterns to arise. Corroborating this, GUS expression under the NIP5;1 promoter is indeed present in the inner root tissues, albeit the NIP5;1 protein levels in those cells are very low, supporting the notion that these files are likely responsive and capable of expressing NIP5;1, but inhibited to do so by the boron levels that arise due to the cells’ positioning within the tissue (Figure 1—figure supplement 1). However, these results, while suggesting sufficiency of our parsimonious assumptions in relation to observed transporter levels, do not offer evidence that tissue-specific regulation is not occurring. At most, the qualitative matches justify, at first instance, that the model is kept simple, without the need of introducing ad hoc tissue-specific dependencies, when exploring the effects transporter regulation dynamics. Furthermore, these results indicate that under our experimentally observed parameter values, the system presents steady nutrient throughput and constant cellular concentrations and transporter levels over time.

The spatial nature of the unstable flows: traffic-jam-like behaviour

We next investigated more carefully the spatial-temporal characteristics of the oscillatory behaviour that arises when transporter swiftness is decreased (Figure 3A). Kymographs make it visually clear that oscillations in the individual cells are in fact spatially coupled, manifested as a boron wave that propagates backwards over the tissue, that is, the pulse moves in the direction contrary to the nutrient flow itself. The initiation of the instability, as seen by the first peak in concentration, occurs close to the xylem (i.e., the last cell, here, cell 4), followed by a strong rise in endodermal concentrations (in cell 3) around 90 min later. Due to boron’s inhibitory influence on the transporters, this intracellular rise leads to a decrease in both NIPs and BORs within that cell. This transporter downregulation causes a drop in boron throughput across that cell, leading to an accumulation of boron in the adjacent, outward facing, cell wall. Background permeability rates along the plasma membrane allow this apoplastic rise in concentrations to trigger an increase boron concentrations in the next outermost neighbouring cell, in this case, the cortex cell (cell 2). Cortex concentration levels thus rise, again triggering a shutdown of transporters, causing the same process to ensue in a spatially coupled manner in the outermost cell, the epidermis. In the epidermis, however, boron levels can never rise above the soil concentration levels, as discussed previously, such that the backward travelling wave loses its amplitude and terminates. The overall process is observed to be cyclical, with levels rising again close to the xylem, until a next wave is triggered in the innermost cells (Video 2). Moreover, we found that such traffic-jam-like behaviour is robust when smaller differences between kB and kN are considered (Figure 2—figure supplement 1A–C), as long as kB is sufficiently large as not to too strongly interrupt boron throughput through the tissue altogether (Figure 2—figure supplement 1D). Obviously, if kB is too low, and traffic-jam-like behaviour are prevented, the plant would also not be able to take up boron, thus, this would be a parameter setting that is biologically irrelevant. Therefore, although our interpretation of the experimental results (Figure 1—figure supplement 2) suggest larger differences between kB and kN, given that these are indirect estimates that rely on underlying assumptions, this robustness analysis shows that, within a biologically reasonable window, our results hold even when kB were to be smaller.

Figure 3 with 1 supplement see all
Traffic jam-like behaviour occurs for various tissue sizes.

(A–G) Kymographs of boron concentration in NIP-BOR models consisting of four (A) up to ten (G) cells. Boron concentrations in cells and cell walls are represented through a heat map colouring with development over time displayed vertically. Boron concentrations exceeding 3000 μM are displayed using the same red colour. Regions of high concentration can be observed to displace to the left, contrary to the direction of boron flow due to the transporter activity. Default parameters were used (Table 1), with exception of the transporter regulation swiftness, ϵN and ϵB, which were both set to 0.01. Above each kymograph a schematic diagram of the tissue structure is shown, corresponding to that of Figure 1C.


We next wondered if this traffic-jam-like effect dissipates over larger tissues. Our four-cell model is a simplified representation of A. thaliana, which has less cell files between the epidermis and the xylem than most other roots of experimental plant models. However, if one distinguishes the pericycle from the vasculature, a five-cell model would be more appropriate. We show that the variables of the five-cell model present similar dynamics as our default model when BOR transporter regulation is lowered to 400 μM (Figure 4—figure supplement 2A–D). Lowering kB slightly stabilises the dynamics (as shown in Figure 2—figure supplement 1), which counteracts increasing destabilisation of the constant flow regime when the number of cell files in the tissue increases. Conversely, given that this renders smaller tissues (such as our default four-cell model), more robust against traffic-jam-like behaviour, our default setting should be regarded as a worse-case scenario for such dynamics to occur, rather than being a special case. To further appreciate tissue-size dependencies, and explore the generality of these effects for other plants which might have highly deviating cell file numbers, we gradually increased the in silico tissue from 4 to 10 cell files. Figure 3B–G shows that the phenomenon is conserved irrespective of the number of cells between the medium and the convective xylem flux. Larger transversal tissue segments are able to generate higher wave amplitudes, but in all cases the peak in boron propagates contrary to the nutrient flow with approximately the same speed. The same phenomenon was still observed even when the model contained 100 cell files (Figure 3—figure supplement 1), suggesting that large tissues do not prevent these effects from occurring, but rather increase its likelihood.

In short, when transporter regulation swiftness is sufficiently slow, due to the downregulation of transporters under the control of increased boron concentrations, a traffic jam-like behaviour emerges: a high boron concentration peak appears, correlating with locally lower transport rates, triggering a low-transport high-boron wave that back-propagates from cell to cell in the direction contrary to the transport polarisation direction within the tissue. This occurs independently of the size of the plant tissue under consideration.

We next studied the importance of the relative dynamics of BOR and NIP regulation to trigger these phenomena, that is, which of the processes have to be sufficiently slow.

Traffic jam-like behaviour depends on swift regulation of both transporters

Under what conditions and parameter values for transporter regulation rates does the traffic jam phenomenon manifest itself? In our previous simulations, we simultaneously decreased the swiftness of both the NIP and the BOR regulation. Oscillations in boron concentration then arose due to a change in stability of both the flow and steady state of boron concentrations. We next probed the system’s stability while varying the regulation swiftness of the NIPs and BORs independently. This was done by analysing the equilibrium of the simplified, A. thaliana inspired, four-cell model. We linearised the ODEs around the equilibrium using a Taylor expansion, and then evaluated the stability of the equilibrium in terms of the largest eigenvalue, for 40,000 different parameter combinations. The phase-portrait that emerges (Figure 4A) shows a large region in which the system becomes unstable and oscillations arise (indicated in red), as well as a region in which sufficient regulation swiftness ensures that the system is stable and oscillations do not develop (indicated in blue). For the stable (blue) parameter settings, any perturbation dampens out (Figure 4B), whereas for parameter settings within the red region, oscillations dominate. The frequency of the oscillations is determined by a relative lack of swift regulation for both transporters, as shown by the colour map in Figure 4. This phase diagram is qualitatively unaffected for the five-cell model at kB=400 μM (Figure 4—figure supplement 2F). Furthermore, we also checked if saturation in the activity of the BOR transporter might obstruct or stabilise the system against the appearance of the oscillations. (The NIPs, being a channel, should be less prone to undergo concentration-dependent saturation effects in permeability.) Contrary to this, we found that such a saturation in permeability of the BOR transporter does not change the likelihood or transporter-response dependency to present traffic-jam-like behaviour, but rather, when they arise due to slow dynamics, the oscillations present much higher amplitudes in concentration variations, with peak boron levels becoming more than twice as high due to the transport saturation (Figure 4—figure supplement 3).

Figure 4 with 3 supplements see all
Stable- and unstable-flow regimes and their parameter dependencies.

(A) Diagram showing combinations of transporter regulation swiftness (ϵNϵB plane) for the four-cell NIP-BOR model for which stable flows (blue) or oscillations, that is, ‘traffic jams’, (red) ensue. (B) Oscillatory periods of the boron concentration variation. The white curved line represents the boundary between the stable and the unstable region, equivalent to the boundary between the blue and red area in (A). Four insets plot C2 [μM] over time [hour] for the indicated parameter values (white dots), illustrating the behaviour of the model at those points within the ϵNϵB plane.


From this combined analysis it becomes apparent that the parameter values we derived to capture wild type dynamics of transporter regulation (i.e., point 1,1) fall into the stable regime. Hence, the wild type system is robust against perturbations and does not present traffic jam-like behaviour. However, independent small reductions from these values in either NIP or BOR regulation rate can cause the tissue to experience oscillations. Not only does the analysis reveal that both transporters are required to present a rapid response, and that this is independent of the BOR-dependent fluxes becoming saturated, but their combined swiftness is synergistic: if one transporter is extremely fast, the other can be bit slower than what is otherwise needed. However, for any parameter setting, a certain speed in responsiveness needs to present to prevent oscillations.

Physiological implications of flow instabilities for the root system

Traffic jams are commonly experienced by vehicle drivers in urban areas. In terms of human transport, traffic jams are an undesirable outcome as they reduce the throughput through the highway and increase the time to reach the destination. Traffic jams displace backwards in space (contrary to the flow) while giving rise to increased car densities at the specific location of the traffic jam (Kesting and Treiber, 2013). All those properties can also observed in our plant tissue model, namely diminished throughput of boron from the soil to the xylem and backwards propagating pulses of locally increased boron concentrations within cells. Moreover, again analogous to real-life traffic jams, traffic-jam-like behaviour in plant nutrient uptake has undesirable implications in the biological context as well. High intracellular boron is detrimental to plant health as it causes DNA damage and ultimately causes cell death (Sakamoto et al., 2011-09). Plants are therefore likely to have evolved mechanisms to avoid high intracellular boron levels to limit boron-induced damage. Lower throughput across the root implies reduced xylem loading and hence a reduction in the boron absorption efficiency, critical for plant growth at low boron conditions.

To better quantify the magnitude of these effects, we performed a simulation screen to evaluate the expected physiological impact of the transporter swiftness in the form of increased boron levels within the cells due to the instabilities, as well reduction in throughput (see Materials and methods for details). We observe increasing maximum levels in boron concentration for slower regulation swiftness (Figure 5A), as well as higher fold-changes in boron concentration (Figure 5B). This suggests that sufficiently rapid transport regulation is important to reduce the risk of DNA damage when either considering absolute concentrations or the increase over the basal equilibrium boron concentrations. Moreover, under the conditions in which traffic jams and large-amplitude variations occur, we measured a considerable reduction in total throughput through the tissue (Figure 5C). In our default model we assume the BOR transporter-driven flux to be linear with the cellular boron concentration. If we instead consider Michaelis-Menten saturated BOR transport (using a biologically reasonable Michaelis constant — the boron concentration at which the flux becomes half as large due to saturation, as well as the concentration at which the flux reaches its half-maximum value — of cB=1000 μM) the system’s capacity of uptake and throughput only marginally decreases. In stark contrast, the detrimental effects of the instabilities that arise due to slowing down the regulatory dynamics are amplified when BOR saturation is considered (see details of saturation implementation in the caption of Figure 4—figure supplement 3). The impact of the traffic jams becomes much more severe, with peak levels during the concentration fluctuations more than twice as high as the default scenario (compare Figure 2 to Figure 4—figure supplement 3). This is a direct consequence of the cell not being able to effectively efflux boron when intracellular levels spike. Taken all into consideration, we infer that selective pressures operating on the root’s capacity to absorb boron optimally and robustly could result in the system evolving to a regime of rapid regulation of transporters.

Physiological detrimental impact of the oscillations for the root.

In the BOR-NIP model with 4 cells, physiological output for the simulation period between 24 and 48 hr was displayed on the ϵNϵB plane. (A) Maximum boron concentration in any cell. (B) Maximum amplitude variation in boron concentration for any cell (highest boron concentration/lowest boron concentration). (C) Average throughput, measured as the output from the last cell, which represents the xylem, between 24 and 48 hr of simulation.


Traffic jam-like behaviour is independent of bottlenecks in the tissue context

Our model thus far consisted of a segment of polarised tissue, of a variable number of cell files, linking the soil to the xylem via a polar flow of nutrients. Consequently, both soil and xylem, although presenting stable characteristics over time, do constitute abrupt boundary conditions: the first sets a constant medium concentration, and the second presents a constant convective flow away from the transversal section. We thus questioned if these boundary conditions are responsible for triggering the traffic jam-like regime that arises under lower rates of transport gene regulation. During the last decades a similar discussion has been taking place regarding what triggers the analogous traffic jams in discrete macroscopic systems such as road networks (Kesting and Treiber, 2013). Although traffic jams are readily triggered by bottlenecks (such as road obstructions, roundabouts, etc.), theoretical models of traffic jam behaviour predicted that the collective systems parameters alone could be sufficient to render the free flow state unstable. It depends on the car density exceeding a critical value as well as how speed relates to the local car density. To prove this theoretical insight, Sugiyama et al., 2008 developed an experiment in which cars were confined to move on a homogeneous circular road. Above a critical car density they observed a transition from a free flow to a traffic jam state, due to the collective effect of the vehicles. They concluded that, although no bottleneck was present, the intrinsic parameters of the system rendered it unstable, such that the smallest of naturally occurring fluctuations was enhanced to disrupt the free flow. This experiment illustrates the non-intuitive notion of traffic jams being generated spontaneously under certain critical parameter regimes.

Inspired by this strategy, and for us to be able to rule out that the origin of the nutrient-flow instabilities within the plant tissue is provided by the soil boron concentration or the xylem flow, we constructed an in silico polarised tissue composed of 6 cells, wrapped up into a ring, that is, the first cell is connected to the last. The ring structure avoids boundary conditions (i.e., bottlenecks) and offers a scenario in which uninterrupted flows (either stable or unstable) can be studied. We allow the dynamics of the polarised tissue to evolve from an initial situation in which concentrations are homogeneous over the whole ring-tissue. When transporter regulation is as swift as in wild type (Figure 6A), a steady flow with no oscillation in the concentrations is observed. When transporter regulation swiftness is decreased to sufficiently low values as to cause unstable flows (as derived from the analysis shown in Figure 4), oscillations arise. After an initial period in which the flows seem constant, random numerical fluctuations in the simulation bring forth a rise of boron concentration in a random cell, which, due to the mechanism of inhibition by boron of the boron transporters, results in local interruption of flows and the back propagation of the boron concentration pulse (Figure 6C,D, Video 1). The frequency of these oscillations is decreased, and the wave broadened, as the regulation speed is decreased (Figure 6C,D). These simulations demonstrate that the unstable flow regime — which we now can conclusively eliminate as being triggered by boundary conditions — continuously self-propagates.

Traffic jam-like behaviour in the absence of bottlenecks: the tissue positioned in a ring.

Kymographs of the boron concentration in the cells and cell walls of the six-cell NIP-BOR ring model, represented by heat map colours and time development shown progressing vertically. The default parameters were used (Table 1), with ϵN and ϵB varied as specified for each kymograph.

Video 1
Simulation of six identical polarised cells endowed with BORs and NIPs orientated as to generate clockwise boron flows.

Time is indicated above the simulation, and dynamics show an initially constant flow, with low boron concentrations in all cells. Due to the intrinsic instability of the system, minute numerical noise generates boron peaks in the cells, a pattern which propagates counterclockwise (i.e., against the flow) through the tissue. Cytosolic boron concentrations, as well as those in the cell walls between cells, are depicted through the colour bar shown in the movie. Transporter levels (of both BORs and NIPs) are characterised through the colour changes along the radial membranes, again depicted through the colour bars shown in the movie.


Mechanisms underlying rapid regulation

Our model analysis highlights the requirement for rapid regulation of both BOR1 and NIP5;1 transporters to stabilise boron flux through root cells and minimise risks associated with transient high levels of boron in cells, both considered to be important constraints for the plant. Current experimental evidence shows rapid downregulation of the BOR1 protein (Takano et al., 2005; Takano et al., 2010), and NIP5;1 transcript (Tanaka et al., 2011; Tanaka et al., 2016), although the fine details regarding the molecular mechanisms remain to be elucidated. The modelling however predicts permanently fast regulation, not only when boron levels drop, but also when they remain constant. This is however difficult to prove with techniques that we have used previously.

To experimentally explore the regulatory swiftness of the bidirectional transporter NIP5;1 at the transcript level, we used single molecule RNA FISH (smFISH). In this method ~40 complementary fluorescently labeled oligonucleotide probes are used to visualise RNA at the cellular level (Raj and Tyagi, 2010; Duncan et al., 2017). Initially, we combined probe sets to detect mRNA for both NIP5;1 and the housekeeping gene PP2A (specifically A3 scaffolding subunit of PP2A, At1g13320) and compared abundance and cellular distribution (Figure 7A) (Czechowski et al., 2005; Duncan et al., 2016). Consistent with published smFISH experiments, PP2A mRNA appeared evenly diffused throughout all cells (Figure 7B) (Duncan et al., 2016). In contrast, we observed many nuclear accumulations of NIP5;1 RNA that were accompanied by few cytoplasmic mRNA copies (125 nuclear accumulations in 149 cells) (Figure 7B). Large nuclear smFISH RNA signals have been reported for highly transcribed genes. They are typically accompanied by many copies of cytoplasmic mRNA and considered to represent bursts of transcription, where multiple Pol II associate with a locus (Battich et al., 2015; Duncan et al., 2016).

Figure 7 with 1 supplement see all
Nascent NIP5;1 RNA accumulates at sites of transcription.

(A) Schematic of probes used to detect NIP5;1 intron 1 (red), NIP5;1 exonic (green) and PP2A exonic (red) RNA. (B) Representative images of cells showing NIP5;1 (green) and PP2A (red) exonic RNA distributions. (C) Representative images showing co localisation of exons (green) and intron 1 (red) regions from NIP5;1 nascent transcripts. DNA labelled with DAPI (blue). Scale bar: 10 μm in (B) and 1 μm in (C).


Arabidopsis genes that contain introns must undergo pre-mRNA splicing before mature transcripts are translated (Morello and Breviario, 2008-06). Co-transcriptional splicing is common for plant genes and this ensures that introns are removed from the pre-mRNA and rapidly degraded close to the site of transcription (Reddy et al., 2013-10). This system allows intronic smFISH RNA labelling to be used as an effective method to identify loci actively engaged in transcription (Levesque and Raj, 2013; Duncan et al., 2016; Rosa et al., 2016). We used this approach to further investigate NIP5;1 transcription and found 72% of cells with at least one NIP5;1 intron signal. As expected, this was lower than the 84% previously reported for the more highly expressed gene PP2A (Duncan et al., 2016). When we combined NIP5;1 exonic and intronic RNA smFISH probe sets we observed 92.5% of nuclear accumulations co-localised with intron signals (n=111) (Figure 7C). This is consistent with the nuclear accumulations representing sites of ongoing RNA production and degradation, rather than separate nuclear storage compartments.

The suggested rapid degradation of NIP5;1 mRNA after transcription was further supported by an independent approach, in which we performed qRT-PCR using probes specific to pre- or mature mRNA to measure mature mRNA/pre-mRNA ratio in root cells. The average mature mRNA/pre-mRNA ratio of NIP5;1 was less than 25% of that of PP2A (Figure 7—figure supplement 1). In accordance with our interpretation of the smFISH results, this again indicates a high degradation rate of NIP5;1 mRNA.

In light of our model predictions, these observations combined point to a highly dynamic sensing system where RNA is turned over at or near the site of transcription, to limit NIP5;1 levels when boron is above a threshold level. We recently demonstrated that NIP5;1 transcript degradation occurs through ribosome stalling, triggering degradation in the cytoplasm (Tanaka et al., 2016). We also identified an upstream sequence that promotes mRNA degradation, but does not trigger ribosome stalling. The data presented in Figure 7 provides evidence of mRNA degradation occurring close to the site of transcription. Together these results suggest that mRNA degradation could occur through two coordinated B-dependent processes; one where cytoplasmic degradation is triggered by ribosome stalling and another localised in the nucleus where RNA is turned over at the site of transcription. Relief of this repression could then ensure rapid boron uptake when required. Our model suggests that such a dynamic mRNA production/decay process could support the necessary swift boron-mediated transporter regulation critical for the evolution of controlled, stable boric acid flows through polarised tissue.

These results also support the notion of a constant turnover (both producing and breakdown the RNA) occurring for regulating NIP5;1 transcript levels. We propose that this ‘wasteful’ process can now be understood in light of constraints of stable flows through polarised tissue. Interestingly, the regulation of the transporters is of a different nature, one predominantly on the level of protein degradation (BOR1), the other occurring on the level of transcript accumulation (NIP5;1).


We conclude that the experimentally observed swiftness of boron transporters’ boron-mediated downregulation can be understood as a necessary mechanism to maintain optimal constant xylem loading rates over time, while also avoiding DNA damage due to oscillations within cells. One could argue that without a substrate-mediated downregulation of the transporters there would not be the need to avoid traffic jam-like behaviour through swift regulatory dynamics in the first place. However, were the system not to present the regulation, the plant would be unable to shield itself from toxic soil boron levels (Figure 8A) and be much less efficient in growing at low soil boron levels (Figure 8B). Thus, some kind of inhibitory regulation necessarily has to be present, to ensure plant plasticity under different environmental conditions. As a consequence, the speed of the transporter downregulation needs to be sufficiently high. We showed that this is an intrinsic requirement, simply because boron transporters are polarly located: The requirement itself does not depend on external boron fluctuations in the medium, number of cell files composing the root tissue, possible saturable activity of the transporter itself nor the magnitude, location or strength of the xylem flow.

Maximum boron concentration and throughput with and without transporter regulation.

(A) Maximum boron concentration; and (B) throughput, for the four-cell NIP-BOR model, plotted for varying boron concentrations in the medium. The default parameters are used (Table 1). When transporter regulation is not considered, NIP and BOR activity were fixed to their unregulated equilibrium value in the wild type model, namely αNξN and αBξB.


Our study has been based on the well-characterised system of directed boron transport, owing to the richness in quantitative experimental measurements regarding BOR and NIP regulation. However, the implications of our derived insights apply to any system in which the dynamics of polar transport of a given substance relies on an inhibitory feedback between the concentration of that substance and its local transport rate. For example, the phosphate and iron transport systems present both polar transporters as well as substrate-dependent regulation, serving as candidates for such phenomena and dynamical constraint avoidance. In the case of the phosphate transport system, transporters localised to the outer side of the epidermis are responsible for uptake (LePT1 and LePT2, Liu et al., 1998), and they have recently been found to carry a domain possibly involved in phosphate sensing which leads to rapid degradation of phosphate transporter in a phosphate-dependent manner (Gu et al., 2016). Similar regulation is also reported for the iron (Fe) transport system, where polar accumulation of IRT1, a major Fe transporter for Fe uptake into roots, is regulated in an Fe-dependent manner (Barberon et al., 2014). Even other transport systems which do not shuttle nutrients might display this sort of instability, such as polar auxin transport, for which stable flows might also require fast response dynamics regarding the activity of auxin importers and efflux facilitators (Robert and Friml, 2009). Indeed, for any polarised tissue in which adaptive regulation of transport levels is necessary we can expect to find intrinsic temporal constraints operating to ensure that steady state flows are maintained and large fluctuations avoided.

In molecular terms, we discovered that although both transporters are predicted to work synergistically and are both required to be effectively similar in regard to rapid regulation, their regulations are biologically encoded in unique ways. We here experimentally focused on the regulatory swiftness of the bidirectional transporter NIP5;1, and found strong evidence that it constitutes a system in which production is maintained high but degradation is rapidly controlled transcriptionally, allowing for the necessarily rapid boron-dependent response. It will be interesting to speculate what underlying reason, if any, has caused the directed exporter, BOR1, which is similarly known and required to be swiftly regulated, to be controlled on a biologically distinct level.

For this work, the analogy with motorway traffic jams was helpful to build an understanding of principles of polarised transport dynamics, albeit both systems obviously deviate in important aspects. Boron concentrations are continuous and can reach arbitrarily high levels, as opposed to (incompressible) cars. Yet, the discrete nature of cells and their regulated polar transporters, combined with the boron-dependent inhibition of the transporters themselves, resembles cars slowing down in response to busy regions along the motorway. Moreover, in both systems highly sensitive responses can prevent jamming (Sugiyama et al., 2008). We therefore propose that the notion of traffic jams can be universally instructive and helpful for biologists considering physiology and regulated growth through polarised transport in plants, as well as for future efforts to enhance plant growth under diverse environmental conditions.

materials and methods

Observation of protein localisation of boron transporters

For plant culture, MGRL medium was used (Fujiwara et al., 1992) supplemented with 1% sucrose and solidified with 1.5% gellan gum. For Figure 1, A. thaliana transgenic lines carrying BOR1-GFP (Kasai et al., 2011) and GFP-NIP5;1 #8 (Tanaka et al., 2011) were germinated and pre-incubated for 5 days on normal MGRL medium plate, and then grown for 2 days on MGRL plates with 0.3 μM boric acid. For Figure 1—figure supplement 2, plants were grown for seven days on MGRL medium with the indicated boric acid concentrations. Images were captured with a confocal microscope (FV1000, Olympus, Japan). GFP fluorescence was observed with 473 nm excitation and 510 nm emission. Cell wall was stained with 10 μg/mL propidium iodide aqueous solution for 3 min, and observed with 559 nm excitation and 619 nm emission.

NIP5;1 promoter activity was observed using a transgenic A. thaliana strain carrying NIP5;1 promoter-GUS (strain −2,180ΔUTR312 #8 in Tanaka et al., 2011). Seedlings were germinated and cultured on MGRL medium plate. Four-day-old seedlings were stained with GUS staining solution: 100 mM Na2HPO4 (pH 7.0), 0.1% Triton X-100, 2 mM K4Fe(CN)6, 2 mM K3Fe(CN)6, 0.5 mg/mL X-Gluc (5-bromo-4-chloro-3-indolyl β-D-glucuronide cyclohexylammonium Salt, Glycosynth, UK) in a decompressed desiccator for 20 min at room temperature. After rinsed with phosphate buffer (pH 7.0), stained seedlings were mounted in 45°C molten 5% agar (gelling temperature 30–31°C, Nacalai, Japan) aqueous solution and solidified at 4°C for 3 hr. The samples were trimmed into blocks and sliced with a vibrating microtome with thickness of 100 μm. The sliced sections were observed with microscope.

Model parameter choices

Although we have constructed our model as parsimonious as possible to qualitatively explore what kind of behaviours could emerge from a polarised tissue with transporter regulation, our model requires a few important parameters to be set. We have sought to explore the possible dynamics of our system while staying well within the expected limits of what is experimentally considered plausible in terms of the boron transport system. Below we give a brief description of the data our parameter choices were based on, and, when possible, how reasonable ranges could be established.

The plasmamembrane boron permeability rate (p) is set to a maximum rate of 8·10-2 s-1 and a minimum rate of 4.4·10-3 s-1, based on the estimation of boric acid membrane permeability by Raven, 1980 and on more recent experimental measurements performed on internodal cells of the charophyte alga Chara corallina (Stangoulis et al., 2001). Note that the plasmamembrane’s permeability rate for boric acid is relatively high compared to other nutrient elements, due to boric acid existing in a non-charged form under physiological pH.

The degradation rates ξB, ξN represent the effective decay rate of the actual transporter activity, not just of the transporter protein. This includes processes such as its removal from the membrane, its inactivation due to usage, etc. As it amalgamates a range of inactivation/degradation processes, we set it higher than typical protein degradation rates, but lower than typical membrane residence times, using a half life of 9 s. The valid ranges for the production rate of the transporter activities, αB and αN, were derived from the degradation rates ξB, ξN, combined with qualitative conditions for the transporter equilibrium levels. It is as yet not possible to obtain these values directly through experiments, as several intermediate steps are involved. Instead we impose that the transporter activity steady state level does not exceed the plasmamembrane’s permeability by a thousand fold:

(8) p<B*,N*<1000p.

Combining these conditions we could then extract the valid ranges for αB and αN. As explained in the main text, we take kB30kN, based on our observations of the sensitivity of BOR and NIP5;1 expression under different boron concentrations (Figure 1—figure supplement 2), further corroborated by earlier data on NIP5;1, as displayed in Figure 1B in Tanaka et al., 2011. We also ran simulations to show robustness of the observed behaviours for smaller variations between kB and kN by assuming lower values of kB (Figure 2—figure supplement 1, Figure 4—figure supplement 2). The parameter dB determines the maximum fold increase in BOR degradation at very high intracellular boron concentrations. It was set to a 50-fold increase, based on rough assessments of experiments in which plants were transferred to high boron conditions, with the caveat that assessment of the temporal change in functional, membrane-located BOR as a function of the intracellular boron concentration is very challenging. Typical cell sizes (height, hc; width, lc) and cell wall width (lw) were estimated from confocal microscope images of A. thaliana roots. The xylem loading rate a was calculated as follows: Hosy et al., 2003 reported that in mature Arabidopsis the rate of transpiration is 0.1 g water/g FW tissue/hour at 70% humidity under daytime conditions. Given an Arabidopsis FW at 18 days of around 0.2 g/plant, transpiration is 20 mg/hour/plant, or roughly 5·106 μm3/s. Considering a stem xylem cross-sectional surface area being roughly 107 μm2, the convective velocity a should therefore be around 5·106 μm3/s / 107 μm2 = 0.5 μm/s.

Numerical analysis

Numerical and analytical calculations were computed using Wolfram Mathematica 10.2. To calculate time development of the model, ODEs were solved numerically using the NDSolve[] function, setting the initial conditions for each variable to zero. The parameters used in the models are shown in Table 1. The solutions of the ODEs were visualised as a line plot or kymograph, using the Plot[] or ListDensityPlot[] functions.

Table 1
Model Parameters.
ParameterUnitDescriptionAcceptable rangeDefault
BORϵB-Time constant for transporter regulation--1
αBμm s-2Production rate of transporter activity4.9×1023.7×10-92.0×10-1
ξBs-1Basal degradation rate1.6×10-67.6×10-2
kBμMBoron concentration for half-maximum in Hill’s function10001600
dB-Amplitude of increased degradation rate by boron100050
nB-Hill’s coefficient--2
cBμMBoron concentration at which the flux reaches its half-maximum value5001000
NIPϵN-Time constant for transporter regulation--1
αNμm s-2Production rate of transporter activity4.9×1023.7×10-92.0×10-1
ξNs-1Basal degradation rate1.6×10-67.6×10-2
kNμMBoron concentration for half-maximum in Hill’s function1000120
nN-Hill’s coefficient--2
Cell sizelcμmCell width20510
lwμmCell wall width20.20.5
hcμmCell height150520
Otherpμm s-1Membrane background permeability of boron8×10-22.3×10-33×10-2
aμm s-1Xylem loading rate (in the last cell)50000.5
c0μmBoron concentration in medium50000300

To evaluate the physiological impact (maximum concentration and amplitude of cytosolic boron, and throughput towards shoots) of transporter regulation swiftness, numerical simulations were repeated for combinations of ϵB and ϵN, which were independently varied between 0.01 and 2, at intervals of 0.01. For each simulation, minimum and maximum boron concentration for each cell and average throughput through the cell row between 24 to 48 h simulation time were recorded using the ‘EventLocator’ method of the NDSolve[] function.

For the ring model, the boundary conditions in the original model were removed by connecting the last cell in the tissue-segment to the first one. This was done by using a tissue composed of six cells and then introducing a cell wall between cell 6 and cell 1. The boron concentration in the novel cell wall is referred to as W6-1. This generates a six-cell ring model. Because within a ring-model the total amount of boron is conserved, introducing an additional ODE for this cell wall component effectively means an overdetermination of variables, causing numerical issues. Hence the value of the new variable was instead defined as follows:

(9) W6-1=(Btotal-lchci=16Ci-lwhci=15Wi)1lwhc,

where Btotal is the total amount of boron in the whole system, which is constant over time as there is no input into nor outflow from the system. The initial boron concentration of all cells and cell walls in the simulation were set to 300 μM, and Btotal was determined accordingly.

The local stability analysis was performed as follows: First, all possible equilibria for the set of ODEs were found using the NSolve[] function. Next, analytical descriptions of the Jacobian matrix of the ODEs around the equilibria were derived using the D[] function. Finally, the stability of the system for varying values of ϵB and ϵN was determined according to the signs of the eigenvalues of the Jacobian matrix. (Please note that changing ϵB and/or ϵN does not change the number of equilibria or their value, only their stability can change.) The system is stable when the real parts of all the eigenvalues are negative, and unstable in any other case. The imaginary part of the largest eigenvalues was then used to estimate the oscillatory period in the unstable region, based on the following conversion:

(10) Oscillatoryperiod=2πλim.

To compare the predicted oscillatory period close to the unstable equilibrium with the actual oscillatory period of the system, the oscillatory period was also calculated numerically (Figure 4—figure supplement 1). For certain parameters individual variables can present rather complicated patterns, for example temporarily maximum values, followed by a slight decrease and then further rise before decreasing again. We therefore used the following algorithm to determine the oscillatory period: The ODEs were numerically solved using the NDSolve[] function. The ‘EventLocator’ method was then used to collect all extrema for three variables, C2, C3 and N2, for the period between 24 and 72 hr. Using the first set of extrema as a reference point, the second and potentially later sets of extrema were then compared against the first set. The first case for which all differences were within ±0.1% for all three variables was then labeled as being identical to the reference point, and the corresponding time interval reported as the oscillatory period.

smFISH probe design

The probes were designed using the online program Stellaris™ Probe Designer version 2.0 and ordered from LGC Biosearch Technologies (Petaluma, CA). For probe sequences see Table 2.

Table 2
Probes used in the smFISH experiments.
NIP5;1 exon probesNIP5;1 intron 1 probesPP2A mRNA

Observation of NIP5;1 RNA

Col-0 seeds were stratified for 2 d at 4°C and then germinated and grown under 16/8 hr light conditions at 20°C on MGRL plates (30μM boric acid). smFISH was carried out on seedlings as described by Duncan et al., 2017. Briefly, 5d old seedlings were fixed for 30 min in 4% paraformaldehyde. Root squash samples were prepared on slides. Tissue permeabilisation was achieved by immersing the samples in 70% ethanol for a minimum of one hour. Two washes were carried out with wash buffer (10% formamide and 2x SSC). 100 mL of hybridisation solution (0.3 M Sodium chloride and 30 mM tri-Sodium citrate dihydrate, referred to as ’2x SSC’, 10% dextran sulfate and 10% formamide) containing probes (final concentration 25 nM) was then left to hybridise at 37°C in the dark overnight. Each sample was washed twice with wash buffer following hybridisation buffer removal with the second wash left to incubate for 30 min at 37°C. Each slide was then incubated with 100mL 4’ 6-diamidino-2-phenylindole dihydrochloride (DAPI, Sigma, St. Luis, MO) (100 ng/mL) at 37°C for 30 min. The DAPI solution was removed before 100 μL 2x SSC was added and then removed. Samples were then incubated for 2 min at room temperature with 100 μL GLOX buffer (0.4% glucose in 10 mM Tris, 2x SSC) before being mounted in 100 μLof GLOX buffer containing 1 μL of glucose oxidase (#G0543, Sigma) and 1 μL catalase (Sigma) under 22 mm x 22 mm No. 1 coverslips (Slaughter, Uppminster, UK) and sealed by nail varnish. A Zeiss Elyra PS1 inverted microscope was used for imaging. A x100 oil-immersion objective (1.46 NA) and cooled EM-CCD Andor iXon 897 camera (512×512 QE>90%) were used to obtain all images. (Quasar®570 probes: excitation by 561 nm laser, signals detected between 570–640 nm. Quasar®670 probes: excitation by 642 nm laser, signals detected between 655–710 nm. DAPI: excitation by 405 nm laser, signals detected between 420–480 nm.) For all experiments a series of optical sections were set up with z-steps of 0.2 μm. All smFISH images were de convolved with Huygens Professional version 16.05 (Scientific Volume Imaging, the Netherlands, http://svi.nl) using the QMLE algorithm with Good signal to noise settings and 50 iterations. Maximum Z projections in Figure 7 were performed using Fiji (an implementation of ImageJ, a public domain program by W. Rasband available from http://rsb.info.nih.gov/ij/).

Determination of mature mRNA/pre-mRNA abundance ratio

Col-0 seedlings were germinated and grown on MGRL plates containing 30 μM boric acid for eight days. Total RNA was extracted from whole roots of about 25 seedlings with four replicates, using NucleoSpin RNA Plant (MACHEREY-NAGEL). To reduce contamination of genome DNA, 500 ng of prepared RNA was subjected to DNase treatment with RQ1 RNase-Free DNase (Promega), following the manufacture’s protocol. Reverse transcription was performed with Super Script III reverse transcriptase (Thermo Fisher Scientific) with random primers. To compare the abundance of cDNAs with different sequences, absolute quantification by realtime PCR was conducted, using PCR-synthesised DNA fragments with known concentration as standard templates. To detect mature or pre-mRNA specifically, exon-exon or intron-exon probes were designed (Figure 7—figure supplement 1 and Table 3). Realtime PCR was performed with Thermal Cycler Dice® Real Time System and SYBR Premix Ex Taq II (Tli RNaseH Plus) (Takara), following the supplied standard shuttle 2-step PCR protocol. As standard DNA templates, DNA fragments corresponding to cDNA of mature- or pre-mRNA of NIP5;1 and PP2A were amplified by PCR with Prime Star GXL (Takara), using Col-0 cDNA or genome DNA as a template. The primers used are shown in Table 3. The PCR products were subjected to electrophoresis and the target fragment was purified from the gel slices with QIAquick Gel Extraction Kit (Qiagen). The concentration of the purified DNA was estimated by NanoDrop 1000 (Thermo Fisher Scientific). Tenth dilution series of 1 pg/μL~1 ag/μL template DNA were applied for standard curve. The copy number of each standard DNA solution was normalized by realtime PCR Ct values of the exon-exon proves which amplify both mature- and pre-mRNA sequences (Figure 7—figure supplement 1 and Table 3). As a negative control, in parallel, after the total RNA extraction, the same samples were processed with the same procedure except reverse transcription, confirming that realtime PCR signals from contaminated genome DNA does not affect the final results.

Table 3
Primers used in qRT-PCR experiments.
NIP5;1_qRT-PCR-1caccgattttccctctcctgatProbes for realtime PCR
NIP5;1_mRNA_FatttaggtgacactatagtaagctcaaagactaaccaaacFor amplification of the standard DNA fragments
Video 2
Dynamics of transporters and fluxes.

Simulation output revealing the details of boron flows over the membranes, due to transporter dynamics and background membrane permeability rates.



  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
  7. 7
  8. 8
  9. 9
  10. 10
    Boron in water, soils, and plants
    1. R Keren
    2. FT Bingham
    Advances in Soil Science 1:229–276.
  11. 11
  12. 12
    Traffic flow dynamics: data, models and simulation
    1. A Kesting
    2. M Treiber
    Springer, Berlin.
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
    FAO Fertilizer and Plant Nutrition Bulletin
    1. RN Roy
    2. A Finck
    3. GJ Blair
    4. HLS Tandon
    Plant Nutrition for Food Security. A Guide for Integrated Nutrient Management, FAO Fertilizer and Plant Nutrition Bulletin, 16, Rome, Food and Agriculture Organization of the United Nations.
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37

Decision letter

  1. Maria J Harrison
    Reviewing Editor; Boyce Thompson Institute for Plant Research, United States

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your article "Rapid transporter regulation prevents substrate flow traffic jams in boron transport" for consideration by eLife. Your article has been reviewed by three peer reviewers, one of whom is a member of our Board of Reviewing Editors, and the evaluation has been overseen by Naama Barkai as the Senior Editor. The following individuals involved in review of your submission have agreed to reveal their identity: Christophe Maurel (Reviewer #2).

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.


Previously, co-authors of this manuscript described the coordinated, polar locations of the NIP5 and BOR1 transporters and the surprising finding of their rapid boron-dependent regulation. Here, the authors have used a mathematical modelling approach to explore the significance of this rapid regulation and its necessity. Their model predicts oscillatory patterns of boron accumulation within a cell file and also a spectacular traffic jam behavior when the swiftness of transporter regulations falls below a certain threshold. Thus the authors argue that rapid regulation of transporters is needed in order to ensure unobstructed trafficking across cells and to avoid potentially toxic spikes of intracellular boron. Their findings provide an intriguing new way of thinking about nutrient transport in plants and the concepts they develop are potentially of significance for other polarized transport systems. The reviewers agree on these very positive aspects and also agree that the following points need to be addressed prior to publication.

Essential revisions:

1) The authors hypothesize complex relationships between B concentration and NIP or BOR protein abundance and thus activity (Equation 5). However, they may oversimplify the relationship between B concentration and transport itself. Proportionality of transport to B concentration gradients (Equations 1 and 2) is consistent for channels such as NIPs. Is it always applicable to BOR1, assuming that this transporter may saturate at high Ci? The authors should know the likely Km of BOR1 (In any case, the bor1-1 mutant lacks a saturating B uptake component): Under steady state conditions (without oscillation), estimated Ci seems to go up to 250 μM and in oscillating conditions, it goes up to 1.5 mM. In the latter conditions, BOR1 must clearly saturate? In our opinion, the authors should explore this important point.

2) We wondered why the authors hypothesize such a huge difference (30-fold) between kB and kN. Is it consistent with actual variations of intracellular B concentration? Data from Figure 1—figure supplement 2 refer to external B concentration and intracellular B may not vary so much. Could the authors present a sensitivity analysis of their model to kB/kNratio?

3) While the authors have chosen a rather perfect example with the boron transport system, we think the concept they develop might be more widely applicable to a range of other nutrients. The authors should discuss other elements for which such a rapid reporter regulation might be necessary (those where coupled-transcellular transport might occur and where toxicity is an issue)

4) Figures 1 and 2. The labeling is exceedingly poor. There is nothing in Figure 2 to allow the reader to know what C1-C4 etc. are (eventually, we work this out). Please add this to the figure legend. It seems that some of the cells have an identity (C3 is endodermis?) but it is not possible to determine this from the figures. What is the gap between W2 and Wn-2? Please add suitable explanation.

5) Although we think that the reduction of radial transport into a cell file is acceptable in the context of this work, we do not understand why the authors chose a four cell model in which the endodermis is the last cell before the xylem (if we understood it correctly). This is clearly not the case, as the endodermis never touches xylem vessels directly, but is always at least separated from it by on layer of pericycle, which might be quite relevant when considering xylem loading. Please provide some explanation for this.

6) Connected to the above, we think the authors stretch it a bit far when they claim that the observed tissue-specific transporter accumulations can be largely explained by feedback regulation of otherwise constitutively expressed transporters. This might work for NIP5;1, but we question whether it can account for the complex accumulation pattern of BOR1 ? We don't think the statement in the text « in the cortex (cell 2) at high levels, whilst in the endodermis it is weaker (cell 3) » reflects what can be seen on the pictures in Figure 1 where BOR1 looks stronger in the endodermis than in the cortex. Please consider providing additional support for Figure 1 (other images?) and discuss these aspects.

7) The experiment in Figure 7 does not evaluate the outcome of the model. Instead, it further evaluates the mechanisms of NIP5 RNA regulation. This is also an interesting question but previous analyses suggested that NIP5 mRNA degradation involved ribosome-stalling controlled by uORFs in the 5' UTR. So how are the current data in Figure 7 reconciled with previously published data in Plant Cell (the Plant Cell data are also generated by some of the authors of this manuscript)?

8) This model has been reported previously (Shimotohno et al., 2015) although it tested different aspects of the system. We are surprised by this omission and request that the authors insert some discussion of the earlier work and its relationship to the current data.

9) The mathematical aspects are very well done with the exception of one small typo in Equation 1B where Wi-1 should be Wiprobably. Please check this.

10) The authors did a robustness analysis to show that the results are not specific to the parameter set they used, but this should be mentioned more explicitly in the text.


Author response

Essential revisions:

1) The authors hypothesize complex relationships between B concentration and NIP or BOR protein abundance and thus activity (Equation 5). However, they may oversimplify the relationship between B concentration and transport itself. Proportionality of transport to B concentration gradients (Equations 1 and 2) is consistent for channels such as NIPs. Is it always applicable to BOR1, assuming that this transporter may saturate at high Ci?

We here assume a linear relationship between the transmembrane flux due to the activity of the BOR1 transporter and the intracellular boron concentration. We agree with the referee that at a certain concentration this process should saturate and this assumption will not hold anymore. We have experimental reasons, however, to believe that this saturation is not occurring within the concentration range considered here. (The experimental justification is described below, answering the next reviewers’ comment). Nevertheless, we found it both important and interesting to investigate what effects such a saturation in the BOR1-driven transport would have on the system, in particular to its tendency to present traffic-jam like behaviour. We therefore altered the BOR-mediated flux (J) as a function of BOR protein levels at the membrane (B) and the intracellular boron concentration (C), from the standard linear regime (J = BC) into a Michaelis-Menten saturating form (J=BC1+C/cB, in which cB is the

Michaelis constant, the boron concentration at which the flux becomes half as large due to saturation, as well as the concentration at which the flux reaches its half-maximum value. We found that as long as cB is not unrealistically low, the traffic-jam-like behaviour occurs at very similar slowed-down transporter response regimes, but the impact of the traffic jams becomes much more severe, with the concentration fluctuations greatly increasing in amplitude (for cB = 1000µm, as now shown in the paper, peak values more than double in respect to the non-saturable peak values). We present the increased detrimental effects on the traffic jam when saturation is considered in a supplementary figure. We have however decided not to use such a term for all simulations presented throughout the paper. Firstly because the linear assumption over the relevant regime is not unrealistic (see below); secondly, to avoid introducing an ad-hoc parameter, its value likely to be too low; and thirdly, to prevent the false impression that the traffic-jam-like behaviour is due to or dependent on such a saturation.

The authors should know the likely Km of BOR1 (In any case, the bor1-1 mutant lacks a saturating B uptake component): Under steady state conditions (without oscillation), estimated Ci seems to go up to 250 μM and in oscillating conditions, it goes up to 1.5 mM. In the latter conditions, BOR1 must clearly saturate? In our opinion, the authors should explore this important point.

We agree that this is a very interesting point (see above). We believe, however, that a linear relationship is not unreasonable, given what we currently know regarding the saturation in BOR1 transport. Our view is based upon experimental data where expression of Arabidopsis BOR1-GFP fusion protein was determined for BY-2 tobacco cultured cells for a range of boric acid concentrations (Yamauchi et al., F1000Research 2013, 2:185). Although the following data was not included in that paper, we measured as part that research the concentration of boron in cells expressing BOR1-GFP fusion. When BY-2 cells expressing BOR1 were exposed to 3, 5 and 10 mM boric acid, intracellular boron concentrations were reduced by ˜10% compared to wild type BY-2 cells, irrespective of the boron concentration in the media. This suggests that the boron transport capacity of BOR1-GFP fusion in BY-2 cell is not saturated in that specific mM concentration range, and thus, we would not expect it to be saturated under the much lower concentrations that our in silico tissues are exposed to.

2) We wondered why the authors hypothesize such a huge difference (30-fold) between kB and kN. Is it consistent with actual variations of intracellular B concentration? Data from Figure 1—figure supplement 2 refer to external B concentration and intracellular B may not vary so much. Could the authors present a sensitivity analysis of their model to kB/kN ratio?

Firstly, we do believe these values are reasonable given the data. The referees note that “Data from Figure 1—figure supplement 2 refer to external B concentration and intracellular B may not vary so much”. It is reported that the boron concentration in root cell sap is almost linear to that in medium, at least within the range of medium boron concentration from 0 to 100µM (Takano et al., 2010).

For the medium boron concentration higher than 100µM, there are no reported direct measurements of cell sap boron. We therefore instead assume that for the higher external B concentrations presented in Figure 1—figure supplement 2, the concentration in the medium should not to be too far off from the concentrations within and sensed by the outermost cells, in its turn regulating their expression.

Nevertheless, whilst by default using the values for kB and kN that are in accordance with Figure 1—figure supplement 2 and the considerations above, we address the concern that the intracellular concentrations in the outermost cells might actually be lower than the concentrations in the medium, giving rise to an overestimation of kB, by strongly reducing the value of kB and analysing its impact on the dynamics. The traffic-jamlike behaviour is marginally affected when kB is reduced from 600 to 400µM, the value we also use for the five-cell-file scenario (see below). When kB is reduced further the amplitude of the oscillations get smaller, while for kB = 100µm no oscillations occur, driven by a dramatic (5-fold) reduction in the boron throughput through the tissue. Clearly, traffic-jam-like behaviour can be prevented by a dramatic reduction in default throughput, but this does not seem to be a very adaptive solution. (These results also presented in a figure supplement).

3) While the authors have chosen a rather perfect example with the boron transport system, we think the concept they develop might be more widely applicable to a range of other nutrients. The authors should discuss other elements for which such a rapid reporter regulation might be necessary (those where coupled-transcellular transport might occur and where toxicity is an issue)

We thank the reviewers of pointing this out. Now that we understand the intrinsic instabilities that the boron system is avoiding, we recognise a posteriori that similar temporal constraints are likely to be operating in a wider set of polar substrate-dependent transport systems, such as in phosphate and iron uptake, and even in auxin transport. In the case of phosphate transport, phosphate transporters LePT1 and LePT2 localise to the outer side of the epidermis to facilitate uptake (Liu et al., 1998), and it was recently suggested that phosphate transporters may carry a phosphate sensing domain that can trigger rapid degradation in a phosphate-dependent manner (Gu et al., 2016). Similar regulation has also been reported for the Fe transport system, where polar accumulation of IRT1 (a major Fe transporter for Fe uptake into roots) is regulated in an Fe-dependent manner (Barberon et al., 2014). We now more explicitly comment about this in the manuscript, as well as regarding the possibility that if under certain conditions auxin were to be exerting inhibitory regulation in the transport efficiency of AUX/LAX, PCP and/or PIN family members, it could also, within a polarised tissue, lead to intrinsic transport, and hence auxin, oscillations. Given that oscillations in auxin have received a significant amount of attention during the recent years, this might provide a potential alternative mechanism underlying these patterns.

4) Figures 1 and 2. The labeling is exceedingly poor. There is nothing in Figure 2 to allow the reader to know what C1-C4 etc. are (eventually, we work this out). Please add this to the figure legend. It seems that some of the cells have an identity (C3 is endodermis?) but it is not possible to determine this from the figures. What is the gap between W2 and Wn-2? Please add suitable explanation.

We apologise for this oversight. Our revised manuscript includes clearer figures and figure legends.

Figure 1: We have now added cell identity labels as well as a simple schematic that explains the mapping of each biological cell to those included in our model.

Figure 2: We have now added an explanation for each variable in the figure legend.

5) Although we think that the reduction of radial transport into a cell file is acceptable in the context of this work, we do not understand why the authors chose a four cell model in which the endodermis is the last cell before the xylem (if we understood it correctly). This is clearly not the case, as the endodermis never touches xylem vessels directly, but is always at least separated from it by on layer of pericycle, which might be quite relevant when considering xylem loading. Please provide some explanation for this.

Indeed, as suggested above, our four-cell model is meant to serve as a concise representation of a generic cross-section, and allows for a quicker visualisation of the dynamics of traffic jam behaviour in the concentration-time plots. Even this concise version of the model already contains 14 coupled ODEs, and, as the previous comment pointed out, this can already be overwhelming and requires careful attention to keep the figures readable and understandable. We are therefore reluctant to choose for the main story a more complex version.

Nevertheless, in response to this comment, we have now repeated the main simulations of the paper in a 5-cell context model and added our results as a supplementary figure to further support the robustness of the mechanism. As can be seen in this supplementary figure, 5-cell tissues manifest traffic-jam-like behaviour when transport regulation is lowered in a very comparable way. In fact, as we also now point out, a higher number of cell files actually enlarges the regime for which instabilities are found. Thus, showing the effect in the 4-cell model is not only parsimonious, but also a ‘worst-case’ scenario for these behaviours to manifest themselves.

Moreover, we now emphasise that while the 4-cell model is a “toy” model, and the 5-cell model a “schematic Arabidopsis” model, higher number of cell files, as displayed in Figure 3 and the supplementary fig, can represent other plant species.

6) Connected to the above, we think the authors stretch it a bit far when they claim that the observed tissue-specific transporter accumulations can be largely explained by feedback regulation of otherwise constitutively expressed transporters. This might work for NIP5;1, but we question whether it can account for the complex accumulation pattern of BOR1 ? We don't think the statement in the text « in the cortex (cell 2) at high levels, whilst in the endodermis it is weaker (cell 3) » reflects what can be seen on the pictures in Figure 1 where BOR1 looks stronger in the endodermis than in the cortex. Please consider providing additional support for Figure 1 (other images?) and discuss these aspects.

We agree, in contrast to the NIP5;1 patterning, which is a clear and robust result, the detailed quantification and interpretation of the BOR1 accumulation patterning is not a main focus of this work, but an interesting side observation that we made from the model. We have now altered the text so that this observation will not be over-interpreted and obscure the main findings of the model. These changes have been made to the relevant section, now entitled: “Parsimonious model, under wild type dynamical settings, generates stable transporter expression”.

7) The experiment in Figure 7 does not evaluate the outcome of the model. Instead, it further evaluates the mechanisms of NIP5 RNA regulation. This is also an interesting question but previous analyses suggested that NIP5 mRNA degradation involved ribosome-stalling controlled by uORFs in the 5' UTR. So how are the current data in Figure 7 reconciled with previously published data in Plant Cell (the Plant Cell data are also generated by some of the authors of this manuscript)?

We appreciate you raising this point. We recently demonstrated that NIP5;1 mRNA degradation occurs through ribosome stalling and that degradation occurs in the cytoplasm (Tanaka et al., 2016). However, we also report that an additional sequence exists upstream of the AUG-stops that also promotes mRNA degradation, but does not trigger ribosome stalling. The current data shown in Figure 7 (in our submitted manuscript), provides evidence of mRNA degradation occurring close to the site of transcription. Considering all evidence together, we propose that mRNA degradation could occur through two coordinated B-dependent processes; one localised in the nucleus, where RNA is turned over at the site of transcription, and in addition cytoplasmic degradation that is triggered by ribosome stalling. An alternative hypothesis could be that NIP5;1 mRNA is transported into the cytoplasm though the nuclear pores close to the site of transcription and the transported mRNA is rapidly degraded by the mechanism described in Tanaka et al., 2016. We have not further tested either of the hypotheses as we feel that it lies beyond the scope of this paper. But we now make this suggestion in the revised manuscript and explain how it can reconcile our current data with previously published work.

8) This model has been reported previously (Shimotohno et al., 2015) although it tested different aspects of the system. We are surprised by this omission and request that the authors insert some discussion of the earlier work and its relationship to the current data.

We now mention and explain the main insights of that previous model and thank the reviewers for this suggestion, as those findings complement this work greatly. We have added this explanation to the section that introduces the model and modelling choices.

9) The mathematical aspects are very well done with the exception of one small typo in Equation 1B where Wi-1 should be Wi probably. Please check this.

Thank you for spotting this typo!

10) The authors did a robustness analysis to show that the results are not specific to the parameter set they used, but this should be mentioned more explicitly in the text.

We now do so, and also have further increased our understanding of the robustness of this phenomenon through the new supplementary simulations that address the points raised above (concerning transport saturation, kB/kNrelationship, and choice of root file number). We now mention explicitly in the text that, taken this all together, there is a structural need of traffic jam avoidance in plants.


Article and author information

Author details

  1. Naoyuki Sotta

    Graduate School of Agricultural and Life Sciences, The University of Tokyo, Tokyo, Japan
    Data curation, Software, Formal analysis, Funding acquisition, Validation, Investigation, Visualization, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon 0000-0001-5558-5155
  2. Susan Duncan

    Department of Computational and Systems Biology, John Innes Centre, Norwich, United Kingdom
    Present address
    Earlham Institute, Norwich, United Kingdom
    Data curation, Funding acquisition, Validation, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon 0000-0001-9581-1145
  3. Mayuki Tanaka

    Graduate School of Agricultural and Life Sciences, The University of Tokyo, Tokyo, Japan
    Validation, Investigation
    Competing interests
    No competing interests declared
  4. Takafumi Sato

    Graduate School of Agricultural and Life Sciences, The University of Tokyo, Tokyo, Japan
    Formal analysis, Investigation
    Competing interests
    No competing interests declared
  5. Athanasius FM Marée

    Department of Computational and Systems Biology, John Innes Centre, Norwich, United Kingdom
    Conceptualization, Software, Formal analysis, Supervision, Investigation, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    Competing interests
    No competing interests declared
    ORCID icon 0000-0003-2689-2484
  6. Toru Fujiwara

    Graduate School of Agricultural and Life Sciences, The University of Tokyo, Tokyo, Japan
    Conceptualization, Resources, Formal analysis, Supervision, Funding acquisition, Investigation, Methodology, Project administration, Writing—review and editing
    For correspondence
    Competing interests
    No competing interests declared
    ORCID icon 0000-0002-5363-6040
  7. Verônica A Grieneisen

    Department of Computational and Systems Biology, John Innes Centre, Norwich, United Kingdom
    Conceptualization, Resources, Formal analysis, Supervision, Funding acquisition, Investigation, Methodology, Writing-original draft, Project administration, Writing-review and editing
    For correspondence
    Competing interests
    No competing interests declared
    ORCID icon 0000-0001-6780-8301


Biotechnology and Biological Sciences Research Council (BB/J004553/1)

  • Athanasius FM Marée
  • Verônica A Grieneisen

Japan Society for the Promotion of Science (25221202)

  • Toru Fujiwara

Engineering and Physical Sciences Research Council (BB/ L014130/1)

  • Susan Duncan

Japan Society for the Promotion of Science (15J11021)

  • Naoyuki Sotta

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.


This work has been supported by the UK Biological and Biotechnology Research Council (BBSRC) via grant BB/J004553/1 to the John Innes Centre and by the Japan Society for the Promotion of Science via Grants-in-Aid for Scientific Research (No. 15J11021 to NS and No. 25221202 to TF). Experiments were partly financed by the OpenPlant BBSRC/EPSCR Synthetic Biology Research Centre (Grant BB/L014130/1).

Reviewing Editor

  1. Maria J Harrison, Reviewing Editor, Boyce Thompson Institute for Plant Research, United States

Publication history

  1. Received: March 22, 2017
  2. Accepted: August 13, 2017
  3. Accepted Manuscript published: September 5, 2017 (version 1)
  4. Version of Record published: September 29, 2017 (version 2)


© 2017, Sotta et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.


  • 1,397
    Page views
  • 310
  • 2

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Computational and Systems Biology
    2. Immunology and Inflammation
    Richard R Stein et al.
    Research Article
    1. Computational and Systems Biology
    2. Structural Biology and Molecular Biophysics
    Wen Ma et al.
    Research Article