Synaptic polarity of the interneuron circuit controlling C. elegans locomotion

Introduction

Caenorhabditis elegans nematode worms possess a very small nervous system composed of only 302 neurons connected by about 5000 chemical synapses and 3000 gap junctions (White et al., 1986). Because of its smallness a precise map of neural connections was possible to construct (White et al., 1986; Chen et al., 2006). This places C. elegans in a unique position among all other animals (Varshney et al., 2011), for which we have at best only rudimentary connectivity data to test various concepts regarding neural wiring and function (Cherniak, 1994; Karbowski, 2001, 2003; Chklovskii, 2004; Chen et al., 2006; Kaiser and Hilgetag, 2006). However, despite this achievement we still have a very limited knowledge about the nature of most of the worm's synaptic connections, i.e., whether they are excitatory or inhibitory.

The simplicity of the C. elegans nervous system does not preclude these worms from executing various non-trivial behaviors such as locomotion, feeding, mating, chemotaxis, etc. (Hobert, 2003; de Bono and Maricq, 2005). To understand the neural basis of these behaviors requires some information not only about the pattern and strength of the connections but also about the type of their synapses. The same neural circuit can perform different functions depending on the signs of synaptic polarities it contains. Specifically, circuits in which excitatory synapses dominate can sometimes become epileptic. On the other hand, networks with only inhibitory connections could be silent, and therefore in many situations useless. Thus, it may seem that some sort of an intermediate regime is necessary for a proper functioning of the nervous system (van Vreeswijk and Sompolinsky, 1996). For example, it was proposed that mammalian cortical networks operate in a dynamic state in which excitation is effectively balanced by inhibition (Haider et al., 2006; Vogels et al., 2011), although anatomical number of excitatory connections dominates over inhibitory in the cerebral cortex (DeFelipe et al., 2002). For nematode worms, a similar issue has been addressed only sporadically. On a modeling level, in the context of a tap withdrawal circuit (Wicks et al., 1996), and experimentally, studying genes that influence the ratio of excitatory to inhibitory signaling (Jospin et al., 2009). We think that this topic deserves more attention both theoretical and experimental, if we are to understand the functioning of worm's circuits (Gray et al., 2005; Sengupta and Samuel, 2009; Ha et al., 2010).

Movement direction in C. elegans is governed by five distinct locomotory command interneurons (AVB, PVC, AVA, AVD, AVE), each in two copies (left and right). All of these 10 interneurons directly connect a downstream group of dorsal and ventral body wall excitatory motor neurons (Chalfie et al., 1985). The topology of connections between the command interneurons is well known (White et al., 1986; Chen et al., 2006), however, their synaptic polarities are not. Conventional thinking is that AVB and PVC control forward motion, while AVA, AVD, AVE control backward motion (Chalfie et al., 1985). This reasoning is based on the fact that the former interneurons connect mainly motor neurons of type B [experimentally shown to be critical for forward locomotion (Haspel et al., 2010)], whereas the latter interneurons connect mostly type A motor neurons [ruling the backward motion (Haspel et al., 2010)]. However, this simple locomotory picture, relying on a single neuron function doctrine may turn out to be too simplistic. Indeed, many laser ablation experiments show that removal of both AVB and PVC reduces forward motion, but does not abolish it completely (see below). Similarly, worms lacking the “backward” interneurons AVA, AVD, and AVE exhibit a comparable frequency of reversals as intact wild type (WT) worms (Piggott et al., 2011). Moreover, the major backward interneuron AVA makes also connections (both synaptic and gap junctions) with the forward B motor neurons (Chen et al., 2006). Thus, perhaps the decision to move in a particular direction is generated by a collective activity of all command interneurons, rather than by an activity of a particular interneuron or a particular connection.

Our aim is to investigate the problem of decision making for movement direction in C. elegans on the level of its interneuron network. The main question we pose is how two antagonistic behaviors, i.e., forward and backward motions, can be controlled by the same circuit of mutually coupled pre-motor interneurons. A strictly related to this question is the problem of synaptic polarities of these interneurons and the input pattern they receive. Specifically, by applying structural perturbations to the circuit we want to determine, using mathematical modeling, which combination of synaptic polarities (together with an input pattern) gives the best match to the experimentally observed locomotory output of C. elegans. This knowledge allows us to answer a question about a relative influence of inhibition and excitation in the command interneuron circuit. Moreover, this approach provides an insight about a degree of interneuron collectiveness in choosing the direction of motion.

Results

The Command Interneuron Circuit for C. elegans Locomotion

To simplify data analysis and mathematical modeling we grouped the left and right members of each locomotor circuit neuron as one model neuron. Thus, in our circuit controlling worm's motion there are five command interneurons, one distinguished upstream polymodal sensory interneuron called ASH, and a modulatory neuron DVA (Figure 1A). A recent study (Li et al., 2006) indicates that DVA plays a role of a sensory neuron in locomotion. We included this neuron explicitly in the circuit, since it has direct connections with body wall motor neurons (Chen et al., 2006), similar to five command interneurons. Because of this similarity, we want to investigate whether DVA can serve additionally as a command interneuron.

The neurons in the circuit are modeled as a single passive compartment with leak conductance. Connections between neurons are either by chemical synapses (of unknown polarity) or by electric synapses known as gap junctions. Chemical synapses transmit signals using graded transmission. The strength of the connection between two arbitrary neurons is proportional to the number of anatomical contacts between them determined from the empirical data (Chen et al., 2006). Additionally, each pre-motor interneuron receives a constant in time excitatory input from upstream (mostly head) neurons, which can be either weak or strong (this is variable in the model). Overall, our model captures long-term averages in neural activities that are associated with average locomotory output in C. elegans. All the assumptions made in the model and equations describing activities of all neurons are presented in Materials and Methods and in Supplementary Information. The main parameters used in the model are described in Table 1.

In an intact circuit for WT worms fluctuations in interneuron activities control forward and backward motion, and the distribution of these activities determines the relative proportion of forward and backward motion (roughly three to one). We wanted to perturb the system and investigate its corresponding output by performing laser ablations of selected neurons in the locomotory circuit. We reasoned that a gradual removal of neurons from this circuit (Figure 1B; see Materials and Methods) would not only affect its physical structure, but also would redistribute the remaining neurons activities, which in turn, should modify the worm's locomotory behavior. In particular, the pattern of interneuron activities in the circuit should change, altering the ratio of times spent in forward and backward motion. Thus, associating experimental average times of forward and backward locomotion for every ablation type with changes in the average activity levels of the circuit model for fixed sensory inputs can allow us to predict synaptic polarities of the interneurons.

Experimental Results

We performed single, double, triple, and quadruple ablations in the interneuron circuit. In total, we generated 17 types of ablation and recorded corresponding mean times the worms spent in forward (T_f) and backward (T_b) motion, as well as in stopped phase (T_s). These experimental results are presented in Table 2. From all the ablations executed, only removal of ASH and PVC neurons increase the time spent in forward locomotion in relation to WT. This is an indication that these particular interneurons have a definitely negative influence on the forward direction, and for that reason are likely to be inhibitory (see below). All other removals have a detrimental effect on the duration of forward motion, even those traditionally associated with backward motion (AVA, AVD). In particular, the ablation of AVA has the most dramatic effect on T_f, leading to its 10-fold reduction in comparison to WT, although the reversal frequency increases only mildly by a factor of 2. Moreover, the AVA ablated worms, including their multiple ablations, spent a lot of time not moving (stopped phase), far more than WT and worms with other types of ablation. For instance, for the combined ablation ASH + AVA + AVB, we find the largest stopped mean time of 1.16 s.

Worms with multiple ablations reverse roughly as frequently as worms with single ablations (Table 2). Generally, the frequency of direction reversals does not correlate well with the average time worms spent in forward motion (Table 2). For example, worms with killed ASH reverse approximately as often as worms with removed AVD, despite the fact that ASH worms spent three times more time in forward motion. Similarly, worms with almost equal T_f (about 0.9 s), i.e., AVB + PVC and AVA + AVB + PVC, differ in reversals by a factor of 2.5.

An interesting result is that ablating the modulatory neuron DVA causes a sharp decline in the forward motion timing in comparison to WT, and a two-fold increase in reversals (Table 2). Also, the combined ablations of DVA with PVC and AVB show a similar property. This clearly suggests that this neuron has a significant influence on the interneuron circuit output, which could be more than just its sensory modulation.

Theoretical Results

In our circuit model there are seven neurons (five pre-motor interneurons, DVA, and ASH), each of them can be either excitatory or inhibitory. Thus, there are 2⁷ = 128 possible copies of the circuit associated with synaptic polarities, i.e., the sign of ϵ_i (for i = 1, …, 7). Two examples of the polarity copies are: (1) ϵ_ASH = −1, ϵ_AVB = −1, ϵ_PVC = −1, ϵ_DVA = −1, ϵ_AVA = −1, ϵ_AVD = −1, ϵ_AVE = −1 and (2) ϵ_ASH = 1, ϵ_AVB = 1, ϵ_PVC = 1, ϵ_DVA = 1, ϵ_AVA = 1, ϵ_AVD = 1, ϵ_AVE = 1, which correspond to all inhibitory and all excitatory neurons, respectively. Additionally, each of the six pre-motor neurons (excluding ASH) receives an upstream sensory input coming mostly form the head, except for PVC for which it comes predominantly from the tail. The sensory input is represented by a binary variable z_i that can have two values: 0 for a weak input, and 1 for a strong input (see Equation 7 in Materials and Methods). Consequently, every polarity copy can be found in additional 2⁶ = 64 activity or input configurations. This implies that, in total, our circuit model has 2⁷ · 2⁶ = 8192 distinct polarity-input configurations.

For each possible configuration (i.e., synaptic polarity and an upstream input) of the circuit we performed 17 “computer ablations” analogous to the experimental ablations shown in Table 2, by setting ϵ_i = 0 if the neuron with an index i was removed. This procedure removes all the connections (synaptic and electric) coming out of this ablated neuron. For each ablation we computed the fraction of time corresponding to forward motion, i.e., T_f/(T_f + T_b), using Equation (8). Thus, we generated 18 theoretical fractional times associated with every circuit configuration (17 types of ablation plus WT) and computed their Euclidean distance (ED) to the experimental values in Table 2 (see Materials and Methods for a more detailed description). The configuration with the smallest ED value corresponds to the optimal solution that predicts synaptic polarities and the pattern of the upstream input. Our strategy was to vary the level of this input (σ, κ), and for each level to find optimal values of synaptic and gap junction conductances (q_s, q_e), together with the system noise amplitude (η), which give the best fit of the theoretical T_f/(T_f + T_b) to the experimental data.

Winning configurations

We find that the best match to the experimental data is achieved in the case when all seven neurons in the circuit are inhibitory (Tables 3–5). The winning synaptic polarity configuration is associated with the combination number 1, which gives the best (lowest) ED = 0.363, and the largest correlation with the data points, which is 0.743 (Table 3). The distribution of winning values of the upstream input z_i is non-homogeneous, and non-zero only for AVB and PVC neurons, implying that much sensory excitation comes to the neurons controlling directly forward motion (Table 3). In Figure 2 we display a comparison of theoretical and experimental values of T_f/(T_f + T_b) across different ablations for this winning combination of synaptic polarities. In general, there is a good fit of the theoretical points to the data (correlation about 0.74).

The second place among all synaptic polarities is taken by the combination # 17, for which all neurons, except AVB, are inhibitory. This combination has ED value very close to that obtained by the winning combination # 1 (Table 3). Moreover, these two combinations appear the most often among the winners, also for other choices of parameters describing the sensory input (Tables 3–6). The third place is taken by a combination with the number 11, in which only two neurons are excitatory: AVD and DVA. This combination also appears quite often among the leading polarities, but its ED is in some separation from the two winning synaptic combinations 1 and 17.

Optimal solution depends both on synaptic and input configurations

The degree of the match between theory and experiment, i.e., ED, depends not only on the pattern of synaptic polarities but also on the pattern of incoming sensory input (z_i) to the network (Figure 3). We noticed that we obtained better fits if we allowed the interneurons to receive a heterogeneous input from upstream neurons. A given synaptic polarity configuration can produce a slightly different locomotory output (slightly variable ED) depending on how many, and which, interneurons receive a strong input (Figure 3). However, despite this subtlety the overall emerging picture is such that the smallest ED values are associated with configurations dominated by inhibitory connections, while the largests ED correspond to networks dominated by excitatory synapses, regardless of the input pattern. Among the configurations with the lowest ED, the most optimal are those with a moderate number (typically 2 or 3) of interneurons receiving strong input. The highest value of ED is about 2.0 and it occurs for synaptic configurations with prevalent excitatory connections, with only a minor influence of the input pattern (Figure 3). These results suggest that the ED is much more affected by the pattern of synaptic polarities than by the input pattern.

Likelihood estimates of interneuron synaptic polarities

To quantify the likelihood of a given synaptic polarity among leading combinations we associate with each neuron a probability that it is inhibitory, for each input value σ and ASH activity level (see last columns in Tables 3–6). This probability is defined here as a fraction of times the ϵ = −1 appears in the row for each neuron. One can notice that some of the locomotory interneurons, such as AVE, AVA, and PVC, are inhibitory with probabilities that are close to 1, regardless of the input magnitude. The polymodal neuron ASH is also in this category. For the rest of the pre-motor interneurons these probabilities are not that high, but nevertheless are at least ≥ 0.5. We also computed an average probability that a given neuron is inhibitory across different input levels coming to the locomotory circuit (average values for all Tables 3–6). These average probabilities are: 1 for ASH, 0.875 for AVA, 0.5 for AVB, 0.5 for AVD, 0.938 for AVE, 0.656 for DVA, and 0.719 for PVC. Thus, from the whole group of seven neurons investigated here and implicated in locomotion control, only AVB and AVD have about equal chances to be excitatory and inhibitory.

One might wonder if these probability estimates hold if we include more winning combinations, not just eight as in Tables 3–6. Including 20 leading combinations for the winning input values σ = 8 mV and κ = 0.6 (corresponding to Table 3), gives qualitatively a similar picture. Again, the neurons AVB and AVD have the highest likelihood of being excitatory, although inhibitory polarities of these interneurons have the lowest ED values (Figure 4). Taken together, the above results strongly suggest that the majority of interneuron connections are inhibitory.

Dependence of ED on model free parameters

The ED between theoretical and experimental fractions of time spent in forward locomotion depends in a non-monotonic manner on model free parameters: q_s, q_e, σ, κ, and η. The first four parameters are neurophysiological in nature and determine the overall balance between currents flowing in the interneuron network. The last parameter, η characterizes a shape of the transfer function between neural activities and behavioral locomotory output, or in other words it characterizes a network noise level.

In Figure 5 we show a dependence of ED on synaptic and electric conductances (q_s and q_e) in the system. In general, ED has a minimum for a narrow range of these conductances, indicating that an optimal solution exists for this range. A qualitatively similar picture emerges also for the rest of the free parameters. Typically, the parameters controlling the strength of the upstream input, i.e., σ and κ, should be in some intermediate range to reach a minimal value of ED. For instance, the winning polarity combinations for σ = 4 and σ = 12 mV have larger ED than that for σ = 6 or 8 mV (Tables 3–6). The same is true for the value of κ, characterizing ASH activity (see Materials and Methods), with the optimal κ = 0.6 for each σ level.

In Figure 6 we show ED as a function of noise level η. Again, either too strong or too weak η increases ED, and there is an optimal value of this parameter for which ED has a minimum.

Discussion

The Main Findings

Using a combination of experimental (laser ablations) and computational (circuit model and optimization) tools we were able to decipher the likely synaptic signs of the interneurons composing the small network commanding C. elegans locomotion (Chalfie et al., 1985; White et al., 1986). It turns out that probably most of these neurons, i.e., synapses they sent, are inhibitory. In particular, the average probabilities that a particular interneuron is inhibitory are: 0.875 for AVA, 0.5 for AVB, 0.5 for AVD, 0.938 for AVE, 0.719 for PVC, and 0.656 for DVA. These numbers suggest that although some of the connections coming out of these neurons might be excitatory, the majority of the connections are clearly inhibitory.

Because of a suppressing nature, the command interneuron circuit must receive a sufficient amount of excitation from upstream (in large part sensory) neurons to be functional, i.e., to appropriately activate downstream motor neurons. Our computational results indicate that the behavioral data are best explained if the command circuit receives a mixed, heterogeneous input [denoted by X_i in our model; Equations (6) and (7)]. The best fit to the data is obtained if the largest excitation comes to forward interneurons AVB and PVC (Table 3; Figure 3). In this sense, the existence of sensory stimulation is an important factor for directional motion generation, which is in general agreement with the experimental findings (Zheng et al., 1999; Piggott et al., 2011).

Role of Interneuron Gap Junctions in Locomotion

How are downstream motor neurons activated given an inhibitory nature of synapses in the pre-motor interneuron circuit? This probably occurs due to strong gap junction couplings between two major interneurons, AVB and AVA, and the downstream motor neurons of type A and B. Our results suggest that the worm moves forward because AVB receives a stronger upstream input than the “backward” interneurons. AVB then excites downstream B motor neurons via strong gap junctions. The chemical synapses from AVB to B are weak and thus non-significant (Table 7). Therefore the sign of these synapses is irrelevant for forward motion, and AVB does not necessarily have to be excitatory.

The issue with backward motion is more subtle. Recent calcium imaging studies show that AVA is active during backward movement (Chronis et al., 2007; Faumont et al., 2011). How can one explain this? The likely answer lies again in the strong gap junction connections between AVA and A motor neurons, and between AVB and B motor neurons (Table 7). Specifically, during backward motion AVA, due to its large sensory input, synaptically inhibits other interneurons including AVB, but at the same time excite downstream A motor neurons via strong electric coupling. Thus AVB sends less excitation to B motor neurons via its strong gap junctions than does AVA to A neurons. Consequently, the activity of A prevails over the activity of B neurons (i.e., E_b > E_f in our model), and the worm moves backward, even with inhibitory synapses in the locomotor circuit.

Relative Importance of the Patterns of Synaptic Polarity and Sensory Input on the Results

Our results indicate that locomotory output of the interneuron circuit depends on its synaptic polarities as well as on sensory input pattern. This can be seen in Figure 3, where the ED characterizing the degree of the match between the theory and the data is displayed as a function of synaptic and input configurations. From these two main factors affecting ED, the synaptic polarity pattern seems to be a more important determinant. This is because the lowest values of ED are generally associated with circuit configurations in which synaptic inhibition dominates, whereas the highest ED's are related to mostly excitatory connections (Figure 3). On the other hand, the sensory input pattern does not exhibit such a simple general trend. More precisely, increasing the amount of strong input in the network does not lead to an explicit decrease (or increase) of ED, but rather to its fluctuations. Instead, ED assumes the smallest values for a moderate number of interneurons receiving strong input, which is usually 2 or 3 (Figure 3).

Properties of ASH and AVB Neurons Strongly Affect the Results

From all neurons present in our locomotory circuit, two neurons ASH and AVB are particularly important for the network performance. This is strongly related to the question of specificity in synaptic polarity and input pattern. The map in Figure 3 shows that the synaptic polarity of ASH is strongly correlated with the level of ED. In most cases, synaptic configurations in which ASH is inhibitory decrease ED, and increase it if ASH is excitatory. The exceptions are configurations with small inputs, where the opposite takes place. Thus, for the most part the synaptic polarity of ASH should be inhibitory in terms of the network optimality.

AVB neuron, which provides the main signal to the downstream motor neurons controlling forward locomotion, is also critical for the optimal solutions. The strength of the sensory input this neuron receives strongly correlates with ED (Figure 3). Generally, strong input in AVB decreases ED, whereas weak input in AVB increases ED, and this effect is essentially independent on synaptic polarity configurations (Figure 3). This suggest that on average AVB should receive a strong excitatory input coming from upstream neurons.

If we combine the action of both ASH and AVB we obtain an interesting picture. ED values are the smallest for cases in which ASH is inhibitory and AVB receives a strong input (Figure 3). In an opposite case, when ASH is excitatory and AVB gets a weak input, ED values are the largests.

Dependence of Network Performance on Model Free Parameters

Performance (ED) of the interneuron circuit depends also on several neurophysiological parameters and on the level of noise in the system. Generally, these parameters require fine-tuning to obtain the best fit between the theory and the data. Synaptic conductance q_s and gap junction conductance q_e determine the relative strength of synaptic and electric connections. Optimal values of these parameters that give the best model performance are always in a physiological range, and yield very similar values ~0.1 nS. However, ED is more sensitive to changes in synaptic conductance q_s than to changes in q_e (ED has a broad minimum as a function of q_e; Figure 5), which suggests that precise values of synaptic conductance are more important.

ED has always a minimum as a function of noise amplitude, regardless of the values of other parameters (Figure 6). This optimal value of noise is in the range 0.7–1.1 mV (Tables 3–6), which seems to agree with the magnitude for voltage fluctuations obtained from electrophysiology (Piggott et al., 2011). The optimality of the system performance for some finite level of noise resembles a phenomenon known as a stochastic resonance, which is a ubiquitous mechanism in many physical and biological systems (McDonnell and Abbott, 2009).

The circuit performance not only depends on the pattern of the sensory input (see above), but also on the strength of this input (Tables 3–6). There exists some optimal level of upstream excitation for which ED has a minimum. Either increase or decrease of this level has a negative influence on the system performance. This non-monotonic dependence can be explained in the following way. When the incoming sensory input is too strong its contribution to an interneuron voltage is much bigger than the contributions coming from synaptic and gap junction couplings. Thus manipulation of parameters associated with connectivity does not change the network output, and one cannot improve the system performance. On the other hand some minimal input is required to stimulate the network, so the signal flows to the motor neurons.

The level of excitation in ASH neuron, controlled by κ, can be also viewed as a measure of an additional input coming to the interneuron circuit. Therefore, a non-monotonic dependence of ED on κ should not be a surprise, which is what we observe in our computational analysis.

Asymmetry in Forward and Backward Interneuron Activities Determine the Likely Direction of Motion

Our results indicate that a decision to move in a particular direction can be made on a small circuit level composed of the six pre-motor interneurons (including DVA). Specifically, the output from these interneurons is fed to the two types of body wall motor neurons (B and A) controlling forward and backward motions, whose relative asymmetric activities (E_f and E_b in our model) determine the likely direction of worm's movement. This simple decision making mechanism can explain 74% of the correlations between the experimental data and computational results (see Table 3). Moreover, this behavioral picture is consistent with the findings in a recent experimental study (Kawano et al., 2011), in which it is shown that the imbalance between activities of A and B motor neurons is a likely scenario for the selection of worm's motion direction. However, in contrast to these authors the imbalance between A and B neurons in our model is caused not so much by a strong AVA-A gap junction coupling, but by the asymmetric upstream excitatory input to command interneurons (in the winning combinations, AVB and PVC neurons receive much stronger upstream excitation than the rest of the circuit).

Qualitative Interpretation of Ablation Data

The interesting experimental result concerning single-neuron ablations is that removal of certain interneurons causes an increase in forward motion timing, while removal of others leads to its dramatic decrease. Specifically, only killings of ASH and PVC neurons increase significantly the time T_f spent in forward locomotion. In relation to ASH, this suggests an important role of the sensory input. There are two surprises here. First, PVC was thought as promoting forward locomotion (Chalfie et al., 1985). Second, given that the sensory ASH neuron makes only weak or intermediate synaptic connections with the command interneurons (with all “backward” interneurons and AVB; see Table 7), it should not have such a strong influence on the motion characteristics. The solution of these puzzles is that PVC and ASH are probably inhibitory, i.e., the strongest synapses connecting these neurons with other postsynaptic targets should be inhibitory. Additionally, ASH should be highly depolarized in order to significantly downregulate the locomotory (mostly backward) interneurons via its synapses. Recent calcium imaging indicates that ASH to AVA synapse is likely excitatory (Chronis et al., 2007). This experimental result does not necessarily contradict our result regarding the negative ASH polarity. Indeed, our results concern only an average total polarity of a given neuron, and it is possible that some of its weak synapses can have a reverse polarity in relation to the strongest. In fact, the ASH to AVA synapse is relatively weak (third in strength out of four coming out of ASH; see Table 7), and it could be dominated by stronger inhibitory synapses to other neurons.

Removal of AVA interneuron causes a large reduction in T_f, despite the fact that this neuron belongs to the “backward” locomotory circuit, and one might naively expect that it effectively prohibits forward motion. Moreover, single and multiple ablations associated with AVA cause an increase in stopped time (Table 2). This suggests that removal of AVA decreases a difference between activities of forward and backward motor neurons, i.e., E_f − E_b may become comparable with a threshold for movement initiation Δ (see Supplementary Information). This in turn may suggest that when AVA is absent, backward motor neurons are more active. Taken together, these results imply that the overall synaptic effect of AVA is most likely inhibitory.

The ablation results for the DVA neuron indicate that it plays a more significant role in the locomotory circuit than just its passive modulation. From Table 2 it is evident that killing DVA has one of the biggest impacts on T_f. Based on this, we hypothesize that DVA might serve also as a command neuron in generation of forward locomotion, which is a novel function for this neuron.

The case with multiple ablations is more complicated. These type of ablations do not have an additive property, i.e., removal of more neurons does not necessarily lead to a progressive drop in the forward motion timings. For example, double AVB and PVC ablation has T_f = 0.91 s, but additional removal of AVD actually increases T_f to 1.33 s (Table 2). The latter may seem paradoxical, however, one has to remember that backward neurons do not act in isolation, but participate in the whole interneuron network activity, and thus indirectly also influence the forward motor neurons. Apart from that, interneurons interact among themselves both synaptically (non-linear in nature) and via gap junctions (bidirectional in nature). This may additionally mask a single interneuron contribution to the locomotory output of the circuit. As a result it is very difficult to predict in advance the effect of any particular ablation on worm's locomotory characteristics in the case of multiple ablations. For this, one needs to perform detailed computations on a network level, as was executed in this study.

Collective, Mutually Inhibitory Interactions Between Command Interneurons Underlie C. elegans Direction of Locomotion

These ablation results suggest that a picture in which a single neuron or a single connection control a specific behavior, advocated in several former studies (Chalfie et al., 1985; Gray et al., 2005), may be oversimplified. Instead, our findings support the idea that behavioral (locomotory) output depends to a large degree on a collective activity of neurons comprising the “functional circuit” (Zheng et al., 1999). That is, the same neuron can participate to some extent in opposite behaviors. Obviously, some neurons or connections in the functional circuits may be more dominant than others for a particular behavior, but the presence or absence of a particular neuron in the circuit is generally not critical for its operation. This partial redundancy in neural function is probably evolutionary driven to ensure a robust circuit performance.

Similarly, none of the interneuron ablations, either single or multiple, abolishes the worm's movement or body oscillations completely. This clearly indicates that none of the interneurons alone is a Central Pattern Generator, which again speaks in favor of collective rather than individual interneuron dynamics as a determinant of locomotion.

Our main result that the pre-motor interneuron circuit has mainly inhibitory synapses is qualitatively similar to two earlier (Wicks et al., 1996; Zheng et al., 1999) and one recent (Qi et al., 2012) study. Wicks et al. (1996) investigated computationally a tap withdrawal circuit in C. elegans and predicted that most interneurons composing it were inhibitory (Wicks et al., 1996). That study concluded that PVC and AVD interneurons were probably excitatory. Our results suggest that AVD is equally likely to be inhibitory as excitatory, whereas PVC with high probability should be inhibitory. The possible sources of the discrepancy can be that Wicks et al. (1996) used an older incomplete connectivity data for the pre-motor interneurons (White et al., 1986), did not include AVE neuron, and used a little different set of neurons in their analysis. In another study, Zheng et al. (1999) hypothesized that the locomotory interneuron circuit should act as an inhibitory switch in order to explain qualitatively data on motion direction transitions. A recent experimental work (Qi et al., 2012) also suggests that the pre-motor interneurons should use inhibition as a main synaptic signaling.

An interesting question is which neurotransmitters mediate inhibitory interactions between interneurons. The most likely neurotransmitter between interneurons is glutamate. In mammalian brains, it is known to be exclusively an excitatory signal, since AMPA and NMDA postsynaptic receptors conduct mostly Na⁺ and K⁺ with an effective reversal potential around 0 mV. However, in C. elegans the situation is more complicated because locomotor interneurons contain apart from these receptors, also GluCl postsynaptic receptors (Brockie and Maricq, 2006). These channels are gated by Cl⁻ (with large negative reversal potentials), and therefore mediate inhibition to postsynaptic cells (Brockie and Maricq, 2006). Specifically, the currents associated with GluCl receptor have been observed in the AVA interneuron (Mellem et al., 2002), and they may also exist in other interneurons.

Given these two types of postsynaptic receptors, it is possible that a single interneuron can have both excitatory and inhibitory synapses on distinct postsynaptic targets. In this case, the synaptic polarities associated with each neuron in our study have an average character. More precisely, the determined probabilities that a given neuron is inhibitory are the fractions of inhibitory connections that the neuron makes with other postsynaptic neurons. Thus, for example, the inhibitory probability 0.5 found for AVB indicates that this neuron sends out about 50% inhibitory and 50% excitatory synapses.

Our Computational Model and its Extension

The theoretical approach in this paper blends a traditional neural network modeling with a probabilistic method for relating network activity to behavioral data. In particular, we envision the nematode locomotion as a three state system, one state for forward, second for backward motion, and third state for no motion. In this system there are transitions between the states that are caused by intrinsic relative activities of A and B motor neurons (E_b and E_f), as well as by the system noise (η in Equation 8). Note that many ablations in Table 2 have the ratio T_f/T_b close to unity, which in terms of our model implies that (E_f − E_b)/η ≪ 1, i.e., for these ablations a stochastic influence of the environment is bigger than the relative activities of A and B motor neurons. This is interesting and shows that sensory noise gains in importance as we remove more neurons from our circuit. This may also suggest that some of the interneurons act as filters for the environmental noise.

Our combined approach allows us to achieve a concrete goal, which is the prediction of synaptic polarities for the well defined locomotory interneuron circuit and the determination of the likely sensory input pattern. This prediction does not depend on a precise form of the transfer function between neural activation and locomotory output (compare Equations 15 and 16 in the Supplementary Information). In this sense our results are robust.

The knowledge of the probable synaptic polarities of the command interneurons may have a positive impact on future modeling studies of C. elegans locomotion. We hope, that this will enable more realistic simulations of the neuronal dynamics that can extend the scope of testable predictions of the current locomotory models (Karbowski et al., 2006, 2008; Bryden and Cohen, 2008). Our method of determining synaptic polarities, which combines structural perturbations with the computational modeling, is sufficiently general that can be also applied to other small functional neural systems in which synaptic polarities are unknown. However, it is important to keep in mind that our model, as every model describing biology, is subject to several assumptions (see the list in the Materials and Methods), and clearly has some limitations. The model does not include several subtle neurophysiological features. For example, a possibility that an individual neuron might release multiple neurotransmitters of similar importance, or that neuromodulators might provide an extra synaptic input throughout the network. These features can be addressed in future more detailed investigations.

Materials and Methods

The ethics statement does not apply to this study.

Collection of Experimental Data

Strain maintenance

For our automated locomotion analysis, we cultured C. elegans at 20°C on NGM plates seeded with the OP50 using standard methods (Brenner, 1974).

Automated worm tracking and data extraction

Worms tested by automated tracking were continuously cultured on E. coli OP50. For assaying various parameters of locomotion, 10 cm non-seeded NGM plates were used. These NGM plates used for recordings were equilibrated to 20°C for 18–20 h. After ablations of individual neurons, the worms were placed on plates with E.coli as a food source to recover. Ablated worms and mock controls were tested within 72 h of the L4 molt. They were then transferred to assay plates containing no food. After 5 min of acclimatization on these plates, the worms were video taped for 5 min. Data presented in this paper represent the locomotory behavior of worms when they were exhibiting “area restricted search behavior”. Data extraction and processing was done using image processing and analysis software as previously described (Cronin et al., 2005). From each video recording of 5 min, we used the middle 4 min, and used the software to derive values for times in forward and backward locomotion, as well as reversal frequencies. In our software, we used a velocity threshold for motion detection. Specifically, if the nematode velocity was below 0.05 mm/s, we classified this as stopped time or “no motion”. Every change of velocity direction that was above this value was classified as a reversal. The average time spent in the stopped phase is T_s. The time T_f the worm spent in forward motion is defined as an average duration of time counted from a moment of moving forward to stopping. Similarly, the time T_b spent in backward motion is an average of times from the initiation of backward movement to stopping. Generally, because of many reversals, the sum of the times T_f, T_b, and T_s is much smaller than the recording time of 4 min. The numerical values of T_f, T_b, T_s provided in this paper are population averages. All incubations and recordings were done in a constant temperature room at 20°C.

Laser ablations

For all species tested, we used the L1 larva stage for our ablations.

Description of the Command Interneuron Circuit Model

List of the assumptions used to construct the model

We make the following major assumptions in the theoretical model:

(1) In the interneuron circuit left and right members of each interneuron are grouped as one interneuron.

(2) Synaptic and gap junction strengths between any two neurons are proportional to the anatomical number of synapses and gap junctions between them.

(3) Pre-motor neurons do not generate sodium-type action potentials but their activities are graded, as C. elegans genome lacks molecules coding for voltage-activated sodium channels (Bargmann, 1998). This assumption is also consistent with electrophysiological observations in C. elegans and related nematodes (Davis and Stretton, 1989; Goodman et al., 1998).

(4) Each neuron releases a single neurotransmitter, or equivalently, there exists a dominant neurotransmitter type for each neuron. Thus, with each neuron in the command circuit we can associate a single dominant synaptic polarity.

(5) Worm's movement direction is determined by a relative imbalance in the activities of excitatory motor neurons of type A and B, which is in agreement with recent experimental observations (Kawano et al., 2011).

(6) Behavioral output of the worm can be formally described in terms of a three-state model. The three states correspond to forward motion, backward motion, and stopped period. Each state has its probability of occurrence, which for forward and backward states is given by an exponential function of a difference between activities of type A and B motor neurons (see below).

Derivation of the model equations

Equations describing interneuron circuit responsible for forward–backward motion transitions are given below. This is a non-linear model based on synaptic connectivity data from www.wormatlas.org [updated version of White et al. (1986) wiring diagram from Chen et al. (2006)].

We start with a standard membrane equation describing the graded dynamics of neuron with an index i (Ermentrout and Terman, 2010):

where V_i is the voltage of neurons i, C is the membrane capacitance per surface area, S is the total surface area of neuron i, V_r is the resting voltage, g_L is the total membrane ionic conductance per surface area that is composed mainly of a constant leak current (typical K⁺ channel conductance is much smaller for voltages close to V_r), g_{e, ij} is the gap junction conductance between neurons i and j. The symbol g_{s, ij} denotes synaptic conductance coming from j presynaptic neuron with synaptic reversal potential V_{s, j}. The function H₀(V_j) is a non-linear sigmoidal function characterizing synaptic transmission, and is given by

where θ₀ is the voltage threshold for synaptic activation, and γ is a measure of steepness of the activation slope. Generally, the synaptic input strongly depends on γ. We set its value at 0.15 mV⁻¹ following an earlier analysis (Wicks et al., 1996) (Table 1). This particular value yields realistic synaptic currents that cause changes in voltage membrane by at most several mV, which agrees with known neurophysiology in other animals (Koch, 1998).

Our goal is to write Equation (1) in a more convenient form for the investigation of synaptic polarities. We assume that the resting potential V_r (when no synaptic or gap junction input is present) for C. elegans interneurons is −40 mV, which agrees with earlier suggestions (Wicks et al., 1996), and it is close to a recent measurement (≈ −50 mV) in AIB neuron (Piggott et al., 2011). We want to re-define the voltage in Equation (1) as a deviation from its resting value, i.e., we introduce ΔV_i ≡ V_i − V_r. Let's denote by V_ex the reversal potential for excitatory, and by V_in the reversal potential for inhibitory synapses. The value of V_ex is around 0 mV (the current in excitatory synapses is mediated by Na⁺, K⁺, and partly by Ca⁺⁺). The value of V_in was reported between −70 mV (Purves et al., 2008) and −90 mV (Koch, 1998) (the current in inhibitory synapses is mediated by Cl⁻). As an average for V_in we take −80 mV. Consequently, we obtain for excitatory synapses V_i − V_ex = ΔV_i + V_r − V_ex = ΔV_i − 40, and for inhibitory V_i − V_in = ΔV_i + V_r − V_in = ΔV_i + 40. The resulting average numerical factors in both expressions have identical absolute values. Thus we can use an approximation: V_i − V_s,j ≈ ΔV_i − ϵ_jA, where A = 40 mV, and ϵ_j characterizes the synaptic polarity of the presynaptic neuron j. The value of ϵ_j is either 1 for excitatory synapses or −1 for inhibitory. Taking the above into account and dividing both sides of Equation (1) by g_LS, we can rewrite this equation as

We can determine the strengths of synaptic and gap junction connections between any i and j interneurons by their anatomical numbers N_{s, ij}, N_{e, ij}, and maximal elementary conductances q_s, q_e, i.e., g_{s, ij} = N_{s, ij}q_s, and g_{e, ij} = N_{e, ij}q_e. The data for N_{s, ij} and N_{e, ij} are available from the data set in Chen et al. (2006) (see Table 7). A typical range of conductances for chemical and electrical synapses is known from neurophysiology of other animals (Koch, 1998). The leak conductance g_L is taken as g⁻¹_L = 150 kΩ· cm² (Wicks et al., 1996), which comes form the neurophysiological measurements in a related larger nematode Ascaris (Davis and Stretton, 1989). The surface area S of all interneurons is very similar and around 15 · 10⁻⁶ cm² (White et al., 1986; Wicks et al., 1996), so we obtain g_LS = 0.1 nS, of which the inverse (i.e., 10¹⁰ Ω) is comparable to an experimental measurement of the total input resistance ~ 0.5 · 10¹⁰ Ω (Goodman et al., 1998). Consequently, we can estimate the ratios present in Equation (3) as: Ag_{s, ij}/(g_LS) = 400 N_{s, ij}q_s, and g_{e, ij}/(g_LS) = 10 N_{e, ij}q_e, where q_s, q_e are expressed in nS. We checked that the term ∑_j g_{s, ij}H(V_j)/(g_LS), which is associated with ΔV_i is generally much smaller than 1, since H(V_j) ≪ 1 for voltages not far away from V_r. Consequently, this term is neglected, which simplifies significantly the resulting equations for interneurons (see below and the Supplementary Information). Thus, we can write Equation (3) in an approximate form as:

where τ = C/g_L is the membrane time constant, w_ij is the synaptic coupling w_ij = 400q_s N_{s, ij}, and g_ij is the gap junction coupling g_ij = 10q_e N_{e, ij}. The function H(ΔV_j) in Equation (4) differs from the original function H₀(V_j) only by a substitution θ₀ → θ, where θ = θ₀ − V_r, i.e.,

It is important to keep in mind that synaptic polarities in Equation (4) are determined simply by the signs of ϵ_j coefficients.

Activity equations for interneurons

Equations describing activities of the interneurons in the locomotory circuit are presented in the Supplementary Information. They are similar in form to Equation (4) with an additional inclusion of the heterogeneous excitatory sensory input X_i, i.e.,

where the subscripts i, j are labels for our circuit neurons. The parameter ϵ_i denotes synaptic polarity of the neuron i, and it assumes value 1 (if the neuron is excitatory), value −1 (if inhibitory), or 0 (if the neuron is absent because of the ablation). Ablations in the circuit remove also electric connections if a neuron on either side of this coupling is killed. To include this effect we rescale the gap junction couplings by ϵ²_iϵ²_j factors. The square in the epsilon assures that we do not get unphysical negative values for gap junction conductance. We assume that the excitatory sensory input X_i coming from the upstream neurons to the interneuron i is constant in time and represented by

where x₀ = 2.0 mV, σ is a variable parameter characterizing the strength of a strong input that is constant in time, and z_i is the binary variable either 0 (weak input) or 1 (strong input). Thus neurons can receive only two types of the input: either weak (x₀) or strong (x₀ + σ). The “input” parameter z_i, similar to ϵ_i, is unknown. We want to find their optimal values for each interneuron.

The above pre-motor interneurons make synaptic and gap junction connections with downstream excitatory motor neurons. Two separate groups of these motor neurons generating forward and backward motion, called B and A, respectively, directly connect locomotory muscles. The activities of excitatory motor neurons (E_f and E_b for forward and backward motion, respectively) have a similar form to that of the interneurons and are given in the Supplementary Information. Generally, activity equations for interneurons and motor neurons are of similar kind to those used before in Karbowski et al. (2008) for analyzing forward locomotion.

We solve Equations (6, 7) using a second order Runge–Kutta method. We are interested only in the steady-state activities of the interneurons and motor neurons. One can think about these activities as temporal averages over sufficiently long periods of time driven by a constant bias input. This simplifying step significantly enhances the feasibility of the analysis. The steady-state values of motor neuron activities, E_f and E_b, are inserted in Equation (8); see below. All the results presented were obtained by choosing initial conditions corresponding to resting voltages, i.e., ΔV_i = 0 for each neuron. We also tried random initial conditions in which neurons start with voltages uniformly distributed in the range −10 to 10 mV. In both cases the steady-state values are the same, which means that steady-state activities are independent of initial conditions, at least in that range.

Theoretical ablations

Ablations or removals of neurons in the model are performed by setting ϵ_neuron = 0. For example, if we remove neuron AVB, then we put ϵ_AVB = 0 in all equations for neural activities.

Values of the connectivity matrix

The strength w_ij of synaptic input to neuron i coming from neuron j is given according to Equation (4) by the expression w_ij = 400q_s N_{s, ij}, where N_{s, ij} is the number of synaptic contacts of neuron i with presynaptic neuron j. The number N_{s, ij} is determined as an arithmetic mean for the right- and left-hand side interneurons. As an example, the right AVB neuron receives an input from both right and left PVC neurons, of which we take an arithmetic mean. Similarly, the left AVB neuron receives an input from both right and left PVC, of which we again take an arithmetic mean. Next, we take an arithmetic mean of these two arithmetic means, and obtain a single value representing average number of synaptic contacts between presynaptic PVC and postsynaptic AVB (N_s,AVB,PVC). The strength of gap junction g_ij between neurons i and j is given by an analogous formula g_ij = 10q_e N_{e, ij}, where N_{e, ij} is the number of gap junction contacts between i and j (arithmetic mean of right and left interneurons). Empirical data for N_{s, ij} and N_{e, ij} were taken from the database in Chen et al. (2006) and are presented in Table 7. Parameters q_s and q_e were taken in the range: q_s = 0.03 − 0.6 nS and q_e = 0.03 − 0.5 nS (Koch, 1998).

ASH neuron

From all upstream neurons we selected explicitly ASH neuron because of its polymodal sensory role. Specifically, it has been implicated in avoidance responses (Kaplan and Horvitz, 1993), which are associated with the modulation of locomotion. We do not write an explicit equation for ASH dynamics, because it receives a massive input form many other head neurons, of which we have no knowledge. Instead, we make computations for four different levels of ASH activity that we set by hand. We choose ASH = κθ, where κ = 0.2, 0.4, 0.6, or 0.75. The value of the normalized threshold θ is set to 45 mV.

Relationship to the behavioral data

From the experimental part we have average times the worms spent in forward and backward locomotion, which we denote as T_f and T_b, respectively. These average times should be somehow related to the average activities of the two type of motor neurons, E_f and E_b. The precise relationship between them is unknown due to the lack of direct physiological data. However, one can expect that domination of E_f over E_b should favor forward motion and its duration, and vice versa.

We can make some progress by using an analogy with statistical physics (Gardiner, 2004), and treating worm's locomotory behavior as a three state system influenced by both deterministic and stochastic factors. These three states correspond to forward movement, backward movement, and stopped time (no motion). There could be transitions between the states driven by sensory input from the environment (either deterministic or stochastic). However, we do not model such transitions. We are interested only in long-term “average” or steady states activities of the system. Our model is motivated in large part by experimental results of Kawano et al. (2011). In that study it was shown that C. elegans motion direction is determined by a relative activity of A and B motor neurons. In particular, it was suggested (Kawano et al., 2011) that when forward motor neurons are much more active than backward (i.e., E_f ≫ E_b), then there should be a high probability of finding the worm in the forward motion. Conversely, if the activity of backward motor neurons dominates over the activities of their forward counterparts (i.e., E_b ≫ E_f), then there is a high chance that the worm moves backward. Thus, it appears that the sign of E_f − E_b plays a key role in the choice of worm's motion direction. Moreover, one can expect that when the activities of both types of motor neurons are comparable (equal or almost equal), then C. elegans likely does not move.

We relate the behavioral observables, i.e., the fraction of time spent in forward motion, with the motor neuron activities using a transfer function known as a sigmoidal logistic function. Specifically, we propose

where η is the noise in the system. This parameter also determines the shape of the transfer function. This type of function is a standard tool often used in many computational studies in neuroscience, and this is the primary transfer function used in this study. Note that for cases in which activities of A motor neurons dominate, i.e., (E_b − E_f)/η ≫ 1, the time spent in forward motion is relatively short, i.e., T_f/T_b ≪ 1. The derivation of Equation (8) is presented in the Supplementary Information.

The goal

We want to determine which combination of neuron polarities: ϵ_ASH, ϵ_AVB, ϵ_PVC, ϵ_DVA, ϵ_AVA, ϵ_AVD, ϵ_AVE, together with their corresponding upstream inputs z_i, yields the best fit to the experimental values of T_f/(T_f + T_b).

Comparison of the theory with the data

In our model there are 8192 distinct combinations of synaptic polarities ϵ_i and the upstream inputs z_i, i.e., different configurations in which the circuit can be found. For each circuit configuration, we perform 17 interneuron ablations in our computer model, and compute theoretical values of T_f/(T_f + T_b) for each ablation. Next, we compute an ED of these values to the experimental values given in Table 2. The ED serves as a system performance (the lower ED the better), and it is computed according to the expression:

where R = T_f/(T_f + T_b) ≡ (T_f/T_b)[1 + T_f/T_b]⁻¹, and the subscripts th and ex refer to theoretical and experimental values of R. The subscript a refers to the ablation number, in the same order as in Table 2. In particular, a = 1 corresponds to the mock ablation, i.e., WT. All the results presented in Tables 3–6 were generated by using Equation (8) for the parameter R_th.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank Cori Bargmann for comments on the manuscript. The work was supported by the EU grant POIG.02.01.00-14-122/09 “Physics as the basis for new technologies—development of modern research infrastructure at the Faculty of Physics of the University of Warsaw” (Franciszek Rakowski), by the Marie Curie Actions EU grant FP7-PEOPLE-2007-IRG-210538 (Jan Karbowski), and by the grant from the Polish Ministry of Science and Education NN 518 409238 (Jan Karbowski).

Supplementary Material

The Supplementary Material for this article can be found online at: http://www.frontiersin.org/Computational_Neuroscience/10.3389/fncom.2013.00128/abstract

References