Neural control of body-plan axis in regenerating planaria

Worms were fixed in Carnoy’s fixative as described in [78], were bleached in methanol, and were subsequently stained using the mouse-anti-synapsin primary antibody 3C11 at 1:50 dilution. 3C11 (anti SYNORF1) was deposited to the DSHB by [82]. A standard immunohistochemistry protocol was used as outlined in [78], using a goat-anti-mouse-Alexa488 secondary antibody (Invitrogen) at 1:250. Samples were mounted using VetraShield (Vector Laboratories) and imaged on a Nikon AZ100 Multizoom Macroscope or on a Leica SP8.

The dsRNA synthesis was performed as described in [83]. The template used was Dugesia japonica mRNA for β-catenin-B, partial cds (GenBank: AB499795). Each dsRNA injection was performed in the pre-pharyngeal area of 1~1.5 cm long worms, consisting of three pulses of 32.2 nL each, for three consecutive days. In order to achieve only a partial block of β-Cat, which was required to test the model prediction, the amputation was performed on the third day of dsRNA injection, at least 6 hours after the injection. Regeneration was at 13°C. As a control, VenusGFP dsRNA was injected, and displayed no phenotype different from uninjected worms. N = 20 replicates were used for experimental and control groups, respectively.

Dynein transport function was inhibited by pre-soaking intact animals in 3 μM Ciliobrevin D (Sigma, stock 25 mM in ethanol) for 3 days before cutting and placing trunk, pharynx or pre-tail fragments into fresh 3 μM Ciliobrevin D solution for 3 days before removing the inhibitor and replacing with Poland Spring water. Regeneration was at 20°C and outcomes were scored after 14 days. N = 47 replicates were performed as well as N = 50 controls with matching concentration of ethanol.

For all other pharmacological treatments, excised pre-tail fragments were treated immediately after cutting for the first three days of regeneration, after which the regenerating fragments were washed in Poland Spring Water and allowed to regenerate for 10 days at 20 °C before scoring regenerative outcomes. All pharmacological stock solutions were aliquoted and stored at -20 °C. Bromocriptine (Tocris 0427) 50 mM stock solution in DMSO was used at a final concentration of 2 μM, MDL12330A (Tocris 1436) 50 mM DMSO stock solution was diluted to 20 μM, U0126 (Tocris 1144) was used at a concentration of 18 μM, and 3-isobutyl-1-methylxanthine (IBMX) (Sigma I5879) was dissolved in DMSO and used at a concentration of 200 μM. A referenced summary of projected downstream targets for these pharmacological agents is given in Table 3. Bromocriptine treatment was also used on N = 52 worms cut into 5 equal fragments each.

To detect the pattern of cilia driven flow, an experimental procedure similar to [56] was used. Planaria were placed ventral side up in a small dish so that they attached to the surface of the water and a small amount of Carmine powder (Sigma) was sprinkled onto their ventral side. The movement of the powder particles was recorded on a Nikon AZ 100M microscope with an Andor DL-604M camera at 100 ms frame rate for 2 mins at a time. Slime and powder build-up was periodically removed. The functionality of the method was controlled by applying the same assay to worms treated with 3% ethanol to block cilia function [85] and worms cooled to 4°C, both of which showed a complete lack of cilia-driven powder movement. To analyze the flow pattern, recordings of 1.5 to 2 s were selected in which the worm showed no muscle-driven movement and analyzed using the PIVlab application [86, 87] in MATLAB (MathWorks). Particle image velocimetry (PIV) was calculated in an interrogation area of 64 x 48 x 32 pixels for subsequent passes for all frames, data was smoothed, and all frames were averaged to result in the shown flow pattern.

Ciliary beat was used as an indicator of cellular polarity, as the establishment of long-range order by physiological mechanisms is a key implication of our model. As a polarized microtubule array within individual cells is a required prerequisite for establishment of planar cell polarity [88–92], and cilia beat direction is in turn controlled by planar cell polarity [93], the beat direction of planaria cilia are read-outs of a net microtubule array alignment in individual cells, and ciliary beat direction is widely used as a readout of planar polarization in tissues [40, 94–99].

It is crucial to develop tissue simulations of proposed signaling networks to fully test quantitative models of regeneration, ensuring that they produce the correct outcomes for known functional data in the field. Additionally, an existing gap in understanding planarian regeneration is how to describe self-regenerating morphogen gradients independent of size-scale, while maintaining adherence to realistic time-scales and conditions. Using realistic values for biological parameters (see S1 Text) and known regulatory interactions (e.g. the canonical Wnt/β-Cat pathway), we developed a model to address these knowledge gaps, and tested its predictive abilities by simulating diverse experimental conditions and comparing the outcomes, either against published results, or new experimental data, in a quantitatively rigorous manner. The main predictions of the full model, and their state of validation by previously-reported or presently-reported experimental evidence, are summarized in Table 2.

The core of the model is a regulatory network including major molecular signals known to be involved in regulation of anterior-posterior axis regeneration in planaria, namely Wnt, β-Cat, ERK, Notum, and Hh (see Table 1). Predicted morphogen signaling gradients were found to be consistent with experimentally-observed planaria body-plan outcomes for homeostasis (Fig 2A) and normal wild-type regeneration (Fig 2A‘o’).

We first tested the model’s ability to correctly predict outcomes of experimental manipulations known to alter regenerative outcomes in planaria by simulating 10 different RNAi and pharmacological interventions as perturbations to key components of the regulatory network (Fig 2). Specifically, we modeled RNAi against Wnt, Hh, β-Cat, APC, and Notum, as well as chemical manipulations of ERK, cAMP, dopamine, and serotonin (Fig 2 and Tables 1 and 3). Of these 10 different interventions, cAMP inhibition by MDL12330A and cAMP enhancement by the phosphodiesterase inhibitor IBMX have not been previously reported. Overall, the model generated stable morphogen concentration patterns, which were predictive of body-plan outcomes observed in experiments (Fig 2) for interventions associated with both abnormal posterior regeneration (e.g. head on the posterior to regenerate 2H, Fig 2Ai to 2Aiii), abnormal anterior regeneration (e.g. loss of head or tail on the anterior, Fig 2Av to 2Ax), or loss of tail (Fig 2Aiv). These model predictions, where regenerated body-plan is interpreted in terms of stable morphogen gradients after cutting the model, matched corresponding experimental outcomes, which we compiled from published work from a number of different labs, as well as from new experiments (Table 1 and 3). Details of predicted morphogen gradients and morphological outcomes for all interventions considered in the model are summarized in S3 Fig.

The model also explains puzzling observations that have been reported experimentally with regards to Notum [32, 34]. In the model, Notum is transcribed at anterior-facing wounds (S4 Fig) and therefore appears as an anterior-polarized gradient under normal regeneration conditions, as is observed experimentally [28, 29, 32] (Fig 4A). Consistent with this anterior localization, experimental RNAi data shows that Notum is required for head induction [32], which our model reflects as a loss of anterior signals following simulated RNAi against Notum (Fig 2Aix). Paradoxically, the same experimental data [32] shows loss of Notum signal in 2H heteromorphoses, which is also recapitulated in our model (Fig 4B and 4C). Conversely, treatments leading to loss of anterior development (e.g. APC RNAi, S3L Fig) are predicted by the model to show very high Notum levels in spite of producing 0H regenerates—this is further consistent with the previously reported data [34].

The model provides insight into these otherwise paradoxical observations with respect to Notum, showing how they occur due to a confluence of four factors (Fig 1E): (1) expression of Notum is induced by β-Cat (indirectly via β-Cat’s induction of expression of a predicted Notum Regulating Factor (NRF), where NRF induces Notum expression); (2) NRF is subjected to vector transport by dynein (towards the anterior of fragments, and in the opposite direction to transport of Hh); (3) Notum inactivates Wnt at the anterior of regenerating fragments; and (4) β-Cat is ultimately what inhibits ERK signaling and therefore anterior regeneration. Taken together, while Notum is required to decrease levels of Wnt [32, 100, 101], and thereby induce anterior regeneration under normal circumstances, interventions that are downstream of Wnt signaling will show non-intuitive outcomes with respect to Notum. Manipulations leading to enhanced destruction of β-Cat (including decreases in Wnt and cAMP levels), decrease β-Cat, and therefore increase ERK and anterior regeneration signaling, while reducing levels of Notum. In contrast, as APC (and other elements of the β-Cat destruction complex) signal downstream of Wnt and Notum, RNAi to APC leads to significant increases in β-Cat levels and therefore to a strong inhibition of head formation with concurrent strong increases in Notum levels due to stimulated expression by increased β-Cat. Finally, the model correctly predicts that Notum works in conjunction with Patched (Ptc) to inhibit β-Cat signaling, such that inhibition of Notum or Ptc both lead to 0H or 2T heteromorphoses [39, 40] (Fig 2Aviii and 2Aix). Please see S3 Fig, which presents modeling results illustrating these trends.

Due to the spatialized nature of our model and the ability to flexibly simulate different physical perturbations, the model can predict regenerative outcomes for fragments of different sizes, shapes and locations within the original animal. This enabled us to extract predictions of differential regenerative outcomes for fragments along the anterior-posterior axis (Fig 4). In the standard case, the model predicts polarized gradients of ERK, β-Cat, and Notum consistent with normal regeneration of single-headed worms in each fragment (Fig 4A). For a partial RNAi block of β-Cat, such that β-Cat is no longer being synthesized but levels have not fully decayed in the system, the model predicts higher levels of morphogen gradient polarity disruption for more anterior fragments (Fig 4B and 4D), which translates to a prediction of higher numbers of 2H regenerates closer to the head. This prediction fundamentally depends on the premise that Hh is synthesized and transported in neural tissue. To test this prediction, we performed a novel investigation of the effects of partial RNAi against β-Cat in changing the incidence of regenerative outcomes of worms cut into 5 fragments along the head-tail axis immediately after 3 consecutive days of RNAi injections. Our experimental observations show the same trend as the prediction of the model, with 100% 2H organisms regenerating from head fragments, and progressively lower incidence of 2H heteromorphoses forming in fragments cut closer to the tail (N = 20), see Fig 4B and 4D.

In contrast to the observation of high 2H regeneration frequencies in more anterior fragments for a partial β-Cat RNAi, the model predicts the reverse pattern for perturbations inhibiting cAMP (or any perturbation enhancing degradation of β-Cat in the canonical Wnt/β-Cat signaling pathway) where a higher level of morphogen polarity disruption, and therefore more 2H incidence in a population of regenerates, is predicted in fragments cut closer to the tail (Fig 4C and 4E). We tested this prediction by chemically inhibiting cAMP signaling through the application of Bromocriptine [25, 42], and did indeed observe the highest number of 2H in the pre-tail fragment, with very low 2H numbers in the pre-pharynx and head fragments (Fig 4C and 4E), with the experimental results showing the same trend as model predictions.

Note that these contrasting trends are not limited to a small sub-set of parameters, but appear to be a natural component of the planaria model, which can also easily be seen in the simpler 1D simulations (see S2C and S2D Fig).

In summary, these results indicate that our proposed model is capable of recapturing and accurately reflecting past and novel experimental data, as well as offering explanations for previously puzzling observations. At the same time, the model can make new, non-obvious predictions for the outcomes of genetic and chemical manipulations, especially in combination with complex physical manipulations, some of which we validated here (see Table 2).

In addition to the signaling network within the individual cell, our model relies on vector transport of certain signaling factors, specifically Hh and NRF, which were identified as necessary in order to create a model capable of describing outcomes of RNAi experiments and the variety of new experiments reported herein (see Table 2). While the vector transport aspect of the model is compatible with several mechanisms of directional transport involving an array of potential cell types, including transcytosis [102], cytonemes [103], and electrophoretic transfer of small, charged morphogens via gap junctions [48, 104], experimental data presented here suggests a role for the nervous system anatomy of planaria in directing regeneration outcomes. We therefore proposed that the directional transport of morphogenic factors occurs with respect to aligned microtubules of the neuronal axons [105–107], where vesicles of Hh are hypothesized to be transported by the motor protein kinesin, while the here-proposed Notum Regulating Factor (NRF) is hypothesized to be moved by dynein. Neuronal transport or expression of these factors has previously been shown for Hh [30, 31, 39]. To test the importance of cellular motor protein activity for the patterning of the AP axis, we inhibited dynein-based transport through the application of low doses of Ciliobrevin D to regenerating fragments, which did not impact ciliary motion as assessed by gliding behavior. Blocking dynein-based transport led to inhibition of head formation during regeneration (N = 23/47 0H), consistent with the phenotype observed with Notum RNAi (see S3I, S3K and S5 Figs). Cell division was not disrupted by Ciliobrevin D treatment, as the tail blastema regenerated at a normal rate.

Although we cannot rule out the possibility that other morphogen-regulating factors are also being transported by dynein, the consistency between the dynein inhibition and Notum RNAi supports the hypothesis that NRF is being transported via dynein towards the anterior blastema where we predict it is responsible for inducing Notum expression and subsequent head formation.

As vector transport of morphogens is proposed to be determined by nerve axon polarity, and morphogen concentration is instructive for regenerated body-plan, we predicted that the average (net) polarity of nerves contained in a cut fragment would drive the patterning of the primary axis in head/tail regeneration. We sought to test this prediction by analyzing regeneration outcomes in a variety of different cutting scenarios involving the monopolar and bipolar nerve polarity maps of 1H and 2H heteromorphoses (Figs 5, 6 and 7).

The high nerve axon density and strong anterior-posterior alignment of axons in the VNC [57] are predicted to lead to VNC dominance of the transport direction in fragments containing the VNC (Fig 5). Therefore, when a 1H worm is cut perpendicular to its head-tail axis, transport along the VNC is predicted to lead to re-establishment of morphogen gradients with their original anterior-posterior polarity in each cut fragment (Fig 1G and 1F), a prediction that has been previously observed in experiments [29, 34]. To further test the involvement of VNC, L-shaped cuts were made into the margins of 1H worm trunk fragments, such that the VNC was diverted towards the induced third blastema via cutting, which induced a side-outgrowth containing either a forward-facing or a backward-facing VNC (Fig 5). Dependent on this directionality, the side outgrowths formed either a head when the forward-facing VNC was diverted (N = 46/46), or a tail when the backward-facing VNC was diverted outwards (N = 43/43), see Fig 5. When the L-shaped cuts were performed so that only tissue outside of the VNC was diverted towards a side outgrowth, no head regeneration was observed (N = 25/25), see S6 Fig. Taken together this data indicates that the VNC is important in generating head or tail regenerative outcomes, depending on cut orientation, and that other tissues, such as muscle, which are known to play an important role in patterning in planaria [18, 20], are not alone sufficient to determine new head/tail identity in tissues. These data further indicate that directional transport on the VNC has the ability to determine head/tail identity, as predicted by the steady-state morphogen gradients of ERK and β-Cat produced by the model for these same two cutting scenarios (Fig 5G, 5H, 5O and 5P).

To be able to functionally test the hypothesis that transport along the vector field formed by the net distribution of nerve axons distributed throughout the planarian tissue is the driver of subsequent patterning events, we sought a test case involving a cut fragment where the average nerve direction was not aligned with the prior head-tail axis and the gradient of morphogens that existed before the injury. Specifically we tested small fragments cut from the margins of the worm between the VNC and the body edge, where synapsin stains indicated a net nerve polarity oriented 90° to the anterior-posterior axis and thereby perpendicular to the original morphogen gradient (Fig 6B and 6C). The regeneration of rectangular fragments cut from this margin area (with cuts on all 4 sides so as to not bias blastema formation, and the long edge of the rectangle corresponding to the original anterior-posterior axis), was tracked over time, as well as observed in detail via a time course of synapsin staining. This revealed that the newly forming head in these fragments developed perpendicular to the original body axis and morphogen gradient and in alignment with the nerve fiber direction (Fig 6C–6F), exactly as predicted by the model (N = 94/94 regenerated 1H). These side pieces very likely contain nerve bodies in addition to the axons oriented 90° to the anterior-posterior axis, as seen by PC2 staining [63] and persistent synapsin staining over multiple days. These nerve bodies are hypothesized, by our model, to be the site of production of Hh and NRF, which are transported along the axons to the long edge of the fragments, perpendicular to the previous axis.

In contrast, control fragments of the same size and shape, but containing a piece of the VNC, regenerated following the orientation of the VNC and in accordance with the head-tail orientation of the original tissue (N = 103/103 regenerated 1H), see Fig 6G–6K. Further depictions and quantification of nerve directionality in regenerating fragments for the VNC-free and VNC-containing fragments are shown in S7 Fig, as well as the S1 Dataset. This confirms a highly non-intuitive prediction of our model, given that it is, to our knowledge, the first case where regenerative polarity of a fragment is fundamentally different from that of the original tissue without any applied treatment. This indicates that the morphogen gradient expressed in planaria tissue for various position control genes (PCGs) is not deterministic for regenerative outcomes but rather is set up by the underlying neuronal polarity which determines patterning.

The 2H worms feature a non-trivial nerve polarity map (Fig 7G and 7H), which is proposed to exhibit bipolarity associated with each head and the location of the midpoint of the animal (Fig 7, where white arrow shows midpoint). This bipolarity was directly observable using cilia flow assay (Fig 7A–7F and S5 and S6 Videos). Cutting 2H worms with respect to the nerve polarity map lead to experimental regeneration outcomes consistent with modeled morphogen gradients resulting from vector transport on the field derived from nerve polarity for a variety of cutting scenarios (Fig 7I–7N). Experiments showed 0% 2H (0/401) for short fragments containing the head (Fig 7I), 8% 2H (11/140) for longer fragments with head but posterior cut not crossing the symmetry line (Fig 7J), 94% 2H (85/90) for cuts containing the head but posterior cut passing the symmetry line (Fig 7K), 99% 2H (89/90) for amputations of both heads in fragments containing the symmetry line (Fig 7L), and 7% 2H (3/41) for amputation of both heads in fragments that did not contain the symmetry line (Fig 7M). Cutting 2H worms diagonally at angles greater than 45° from the symmetry line regenerated 88% 2H (85/97) in fragments that showed regeneration (Fig 7N). Experimental regeneration outcomes showed statistically significant corroboration with model-predicted outcomes (based on Markov-model body-plan inference from predicted morphogen gradients) for the same fragments (Chi2 test, p < 0.05). This modeling and experimental data in 2H worms further indicates that transport of morphogens occurs with respect to nervous system polarity, and demonstrates further correct predictions of a wide range of non-intuitive regeneration outcomes.

Models were run with planarian body shapes ranging from 0.3 cm to 2.4 cm in length, and were found to produce morphogen gradients with similar pattern form and of similar magnitude in models of very different size (S8 Fig). Moreover, cutting these models into 3 fragments showed that all models, except for the smallest, were able to fully re-polarize morphogen gradients (S8 Fig). Small fragments (below approximately 1 mm in length) begin failing to regenerate normal gradient polarities, with 0T and 2H heteromorphoses predicted as fragments become smaller (S8 Fig), an outcome which is experimentally observed for small or thin fragments [108–110].

Taken together, the results of the above-described analyses and new experiments show that our model correctly predicts the outcomes of a range of diverse amputations that probe the relationship between positional information in the intact worm and mechanistic transport processes that set the axial polarity of the resulting fragments. We conclude that regenerated anterior-posterior axial polarity depends on the average polarity of nerve axons contained in a regenerating fragment, and that these are dominant drivers in cases where their direction diverges from the pre-existing directionality of the prior AP axis.

Despite major advances in understanding cellular and molecular mechanisms required for regeneration and patterning, achieving a deep understanding of underlying large-scale pattern control during development and regeneration remains an outstanding challenge, crucial to achieving the goals of biomedical regeneration and the mitigation of birth defects. Here we utilize a tight integration between computational modeling and experiment to develop a context for understanding fundamental anterior-posterior body-plan control mechanisms in regenerating planaria, where instructive cues are ultimately described in terms of molecular positional information patterns. We present a comprehensive, spatialized model that correctly predicts and explains existing data, and also makes a number of surprising, novel predictions that were confirmed by experiments.

This work provides evidence that augmenting reaction-diffusion mechanisms with vector transport on hypothesized, though realistic, nerve polarity maps improves on limitations of reaction-diffusion models in a manner that is consistent with a variety of experimental observations for planaria. While reaction-diffusion mechanisms can elegantly generate various pre-patterns of molecular information to account for self-assembly of morphogens in a variety of systems [6, 17, 111], their patterns are often dependent on the size-scale of a tissue, and are therefore not sufficient for use in regenerating planaria [48, 49]. While scale-free models have been developed to address the limitations of traditional reaction-diffusion schemes [49], they have been found to have several limitations for planaria, including [48]: (1) having steady-state concentration magnitudes highly dependent on the size of the model, which presents difficulties for downstream biological signaling where consistent concentrations are required in fragments of very different sizes; (2) requiring many weeks to form a gradient for a model sized on the order of an average planaria length (~1 cm), whereas in reality, morphogen gradients reach a stable state at 72 hours or less [28, 30]; (3) being unable to describe gradients for both 1H and 2H worms with the same model [48, 49]; and (4) as gradients of the scale-free model form with respect to the longest length of a fragment, the model may produce non-realistic gradients in cut fragments [48, 49]. Therefore, after extensive assessment of a variety of reaction diffusion models, we concluded that reaction-diffusion in isolation appears unable to appropriately account for self-assembly of accurately instructive morphogen gradients in planaria fragments of a large variety of sizes spanning several orders of magnitude (e.g. 0.05 to 5.0 cm) and contexts (e.g.

2H worm regeneration).

In contrast, vector transport schemes described herein create self-assembling morphogen gradients that are virtually independent of scale, form with respect to the polarity of the transport field rather than with respect to the longest dimension of a fragment, and have relatively consistent concentration magnitudes independent of size [48, 50, 51, 62]. Previous reports have proposed vector transport of neurally-produced Hh on the VNC in an anterior to posterior direction as a key mechanism for maintaining axial polarity in regeneration [30, 39], however, this hypothesis has not been formally explored nor modeled. Our model expands on this original concept, predicting that in addition to Hh, a Notum-inducing substance such as NRF must be subjected to vector transport, and also, that the transport field includes all nerves, not only the VNC.

A foundational hypothesis of our model is that certain morphogens (NRF and Hh) are subjected to vector transport, and that this transport mechanism is crucial to account for the scale-free self-assembly of morphogen gradients directing the body-plan in regeneration (Fig 1). We herein show experiments, both of regenerative outcomes following various cuts, as well as chemical inhibition of dynein transport, that strongly support the concept that vector transport occurs with respect to the planaria nervous system. It is important to note that the vector transport aspect of the model is compatible with several mechanisms of directional transport involving an array of potential cell types, including transcytosis [102], cytonemes (which have been observed to transport both Wnt and Hh) [103, 112–114], and electrophoretic transfer of small, charged morphogens via gap junctions [48, 104]. However, given our experiments have provided multiple lines of evidence for a role of planarian nerve orientation in directing regeneration outcomes, we have highlighted axoplasmic transport, a well-known phenomenon mediated by motor proteins operating on the aligned microtubules of neuronal axons [115], as the most promising candidate for vector transport mechanism underlying planarian regenerative axial control. Importantly, evidence exists that morphogens such as FGF [116] and β-Cat [117] are subjected to axoplasmic transport in rat neurons, where morphogens are transported from source to target cells via synapses [118]. For consistency with published experimental observations involving Notum, it was necessary to introduce NRF as a β-Cat induced, neurally-transported factor, which in turn induces Notum expression.

While the molecular identity of NRF will be explored in future work, one possibility that corresponds with all existing experimental information to date is that NRF is Notum mRNA transported in the nerves, where axoplasmic transport of mRNA is a fairly well-established phenomenon [65, 119–121].