Front. Phys. Frontiers in Physics Front. Phys. 2296-424X Frontiers Media S.A. 731941 10.3389/fphy.2021.731941 Physics Original Research Resilience of Nematode Connectomes Based on Network Dimension-reduced Method Dongli et al. Resilience of Nematode Connectomes Dongli Duan 1 * Xixi Wu 1 Shubin Si 2 School of Information and Control Engineering, Xi’an University of Architecture and Technology, Xi’an, China School of Mechanical Engineering, Northwestern Polytechnical University, Xi’an, China

Edited by: Haroldo V. Ribeiro, State University of Maringá, Brazil

Reviewed by: Bradly J Alicea, Orthogonal Research and Education Laboratory, United States

Matjaž Perc, University of Maribor, Slovenia

 *Correspondence: Duan Dongli, mineduan@163.com

This article was submitted to Social Physics, a section of the journal Frontiers in Physics

 18 10 2021 2021 9 731941 28 06 2021 20 09 2021 Copyright © 2021 Dongli, Xixi and Shubin. 2021 Dongli, Xixi and Shubin

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

The whole map of nematode connectomes provides important structural data for exploring the behavioral mechanism of nematodes, but to further reveal the functional importance and resilience pattern of nematode neurons, it is necessary to effectively couple the regulatory relationship between neurons and their topology. Here, with a typical signal excitation function we propose a model to capture the interacting relationship between the neurons, because a differential equation depicts the activity of a neuron, n neurons mean we need high-D differential equations to capture the neural network. With mean-field theory, we decouple this N-dimension question into a one-dimension problem mathematically. In our framework, we emphatically analyze the characteristics, similarities and differences of the structure and dynamical behaviors of the neuronal system for Caenorhabditis elegans and Pristionchus pacificus. The comparing results of simulating method and theoretical approach show that the most important homologous neurons between C.elegans and P.pacificus are I2 and NSM, which may lead to their different behavior characteristics of predation and prey. At the same time, we expect that the x eff index can be used to reveal the importance of neurons for the functional evolution and degeneration of neural networks from a dynamic perspective. In the hermaphroditic and male C.elegans, we test the control level of the intermediate neuron groups over the output neuron groups and the single neuron. These results suggest that our theoretical approach can be used to reveal the effects of bio-connectivity groups, potentially enabling us to explore the interaction relationship of neural networks in humans and animals.

 resilience mean-field theory network dynamics nematodes connectomes neural networks 1 Introduction

Due to special physiological structure and easy modeling, nematode has become the primary model to reveal the neurons structure and functional mechanism for humans and animals [1, 2]. After four decades of exploration, nematodes are the first organisms, for which the wiring diagram of their entire nervous system has been mapped at the cellular level [3, 4]. Even Witvliet et al reconstructed full brain of eight isogenic Caenorhabditis elegans individuals across postnatal stages to investigate how it changes with age [5]. This contributes to a more detailed and accurate understanding and exploration of human and animal behavior at the neuronal level, as well as an analysis of the functional importance of individual neurons at the systemic level. Therefore, we coupled the structural information of the neuronal map of the nematode and the functional relationship between the neurons, predicted the functional involvement of specific neurons or synapses in defined behavioral responses [6, 7].

Resilience is system’s ability to retain its basic functionality while errors, failures and environmental changes occur, which universally presents in most dynamical systems [8]. Thanks to the rapid development of network theory, the description of the nonlinear dynamics that governs the interactions between the neurons have been extracted. According to external disturbance for neuron systems, deleting interneurons or synapses may disconnect the networks, or make a great impact on the individual behavior [9]. The methods to explore the nematode’s resilience mainly can be divided into two categories. 1) Biological experiments: such as the gene labeling and the neuron ablation. They can track some specific synapses or their function in detail. The labeling approach can capture the directionality of synaptic connections, and its quantitative analysis of synapse patterns display excellent concordance with electronic micrograph reconstructions [10, 11]. However, this approach takes a lot of energy and is not universal. 2) Theoretical experiments: using different mathematical methods and engineering mechanisms, which can abstractly study the function of outstanding nodes in nematodes [12–16]. For example, according to the theory of symmetry group, they found that the symmetry of the neural network has directly biological significance, and its correctness can be strictly proved by using the mathematical form of symmetry group. This form makes it possible to understand the importance of the structure-function relationship [17]. In addition, Yan also applied network control principles to the connectomes, which reveals both neurons with known importance and neurons which was previously unknown [18]. Within a network control principle, they add input signal to control the output of neurons by a linear framework. The results remain reliable with a small amount of disturbance to the reference connectors, but large disturbances are likely to cause distortion.

Here we construct a framework to explore the resilience of nematode connectomes: we use the signal excitation model to describe and solve quantitatively structure-function relationship of the nematode, which combines the biological behavior with dynamic behavior. The application of signal excitation mode is key to revealing how different behaviors, which contributes to understand the mapping from network structure to function [19]. Since the network contains a large number of neurons, each activity of synapse node or neuron node can be seen as a solution of the high-D differential equations. However, calculating the stable state of the neural system could be difficult or almost impossible, especially when the system is large-scale. With a network dimension-reduction method, we can derive an effective one-dimensional dynamic model which captures the system function of the neural networks. Finally, we use a weighted average activity x eff to measure the neural performance of the whole nematodes. Meanwhile, we introduce a weighted average connectivity β eff to express the structural strength of the nematode neuron system.

With the application of our framework into the neuron networks of P.pacificus and C.elegans, the differences between the predator and prey behaviors of these two nematodes can be clearly quantified at the neural network level, and the key neurons leading to the differences in nematode behaviors can be easily identified as well. Similarly, for hermaphroditic and male C.elegans, we apply our approach to discuss the core control of the intermediate neuron groups to output neuron groups. Core control means that the interneurons have great influence over the output neurons or all the neurons after neuronal perturbation.

Actually, our approach is a general mathematical framework, which couples the complex structure of the system with the nonlinear activation function, and enriches the research methods for exploring the functional behaviors of nematode. We can reveal the underlying principle of the neural network in this way. In the future, our theoretical methods can be applied to more specific linkage groups, to do targeted quantitative research [20].

 2 Neural Model and Network Dimension-reduced Method 2.1 Neural Dynamics Model

There are nonlinear regulation relationships between each synapse in neural networks. Taking hyperbolic tangent excitation function as an example in this paper, we give a neural network decoupling method with a general network structure. For the stimulation relationships between the trillions of neuron synapses in human brain or animals, it would enable the individuals show a much more complex nonlinear phenomenon of affine transformation. The activation function can be used to describe the coupling mechanism among synapses as d x i d t = I − x i R + J 2 ∑ j = 1 N A i j ( 1 + tanh ( n ( x j − α ) ) ) , Where, x i describes the activity of the neuron i, I is the basal activity of the neurons, R is the inverse of the death rate, α is the firing threshold, and J is the maximal interaction strength between pair of neurons. The adjacency matrix A depicts the topology of the neuron networks. The coefficient n governs the steepness of the sigmoid function, analogously to the Hill coefficient n in gene regulatory networks [18]. Actually, the Eq. 1 is made up of two parts: I − x i R describes the growth rate of the neuron itself, J 2 ∑ j = 1 N A i j ( 1 + t a n h n ( x j − α ) ) describes the positive excitation effect of its neighbor neurons on neuron i.

 2.2 Network Dimension-Reduced Method

When the system of Eq. 1 tends to a relatively stable state, which means that the activities of the N individuals in the network are constant values, so the N-dimensional nonlinear rate equations all equal to 0. However, calculating the stable state of the system requires to solve the N-dimensional nonlinear rate equations numerically or analytically, which could be very difficult, or almost impossible especially when the system is large-scale. We present a new dimensionality reduction to explore the interactive relationship between neuronal connectomes for nematodes with mean-field approximation. This approximation is exact in the limit where the node activities are uniformly distributed.

2.2.1 Theoretical Framework for Network Resilience

Once we get the value of x i for the neuron i from the solutions of Eq. 1, we can quantify the performance and the connectivity strength of the neural networks with a weighted average metrics on the network topology following the way in [21]. Defining an operator for the degree sequence of A as Γ y = 1 ⊤ A y 1 ⊤ A 1 = ⟨ y n n ⟩ = s j o u t y j s j o u t .

Here, y i is a variable of neuron i, and y is the weighted average of all the neurons. s j o u t is the outdegree of neuron j. In contrast, s i i n represents the indegree of neuron i. This operator is used to calculate the weighted average y over all nearest neighbor nodes of network A.

In order to simplify the N-dimension of Eq. 1 and approximate it to 1D formalism, we use a mean field hypothesis. Here we define I − x i R as F(x i ), define ( 1 + t a n h n ( x j − α ) ) as G(x i , x j ). Using Γ(G(x i , x)) ≈ G(x i , Γ(x)) approximate Eq. 1, we get d x i d t = I − x i R + s i i n J 2 ( 1 + tanh ( n ( Γ ( x ) − α ) ) ) .

For the states of all nodes in the system, Eq. 3 can be rewritten by the Hadamard product as d x d t = I − x R + s i n ◦ J 2 ( 1 + tanh ( n ( Γ ( x ) − α ) ) ) . Where we used the Hadamard product ◦, namely a ◦ b = ( a 1 b 1 , … , a N b N ) T . Eq. 4 is an equation composed of N dimensions. In order to reduce the dimensionality of this set of equations, we apply Γ ( x ) = 1 T A x 1 T A 1 to both sides of Eq. 4, namely d Γ ( x ) d t = Γ I − x R + s i n ◦ J 2 ( 1 + tanh ( n ( Γ ( x ) − α ) ) ) .

Using Γ(F(x)) ≈ F(Γ(x)) and Γ(G(x, Γ(x))) ≈ G(Γ(x), Γ(x)), we get the equation d Γ ( x ) d t = I − Γ ( x ) R + Γ ( s i n ) J 2 ( 1 + tanh ( n ( Γ ( x ) − α ) ) ) .

In the process of dimension reduction, we introduce two variables to describe the state of the system, where x eff measures the neural performance of the whole nematodes, β eff expresses the structural strength of the nematode neuron system, x e ff = Γ ( x ) = 1 T A x 1 T A 1 , β e ff = Γ ( s i n ) = 1 T A s i n 1 T A 1 = 〈 s o u t s i n 〉 〈 s 〉 .

Lastly, we have successfully decouple this N-dimension question into a one-dimension problem mathematically. Bring Eqs 7, 8 into Eq. 6, one gets d x eff d t = I − x eff R + J 2 β eff ( 1 + tanh ( n ( x eff − α ) ) ) .

Subsequently, when the system comes to a stable state, it is obvious that d x eff d t = 0 . Hence, we can further get that β eff ( x eff ) = 2 J x eff R − I 1 1 + tanh ( n ( x eff − α ) ) .

Equation 10 describes the theoretical relationship between the topology of the neuron networks and its dynamical performance. Actually, by Eq. 10 we have decoupled the complex system of Eq. 1 into a 1D problem.

 2.2.2 Simulation Method

From Eq. 10, we can obtain the theoretical curves of x eff for β eff. To verify the correctness of our dimension-reduced method for the system resilience of x eff, we use the ODE45 function to solve the neural dynamics of Eq. 1, which can ensure we get the numerical results of x i . Then, with the steady state vector x and Eq. 2, which can separate the system from the topology parameter β eff and the behavior parameter x eff. The steps to simulate the neural dynamics are following

1. use ODE45 function to solve Eq. 1. We can obtain the stable vector x. In our study we set I = n = 2, R = α = J = 1. The setting of this parameter we refer to [16]. Assuming the network has N neurons, the solution of the steady state vector is x = [x 1, x 2, …, x N ]. It is easy to know that in the steady state vector, the value of column i represents the steady state activity x i of neuron i.

2. compute the x eff and β eff from Eq. 9. For the vector x obtained in step 1 and the network indegree sequence s in , we can get the system topology parameter β eff and the behavior parameter x eff with the operator Γ y .

3. compare the simulation results with the theoretical solutions. From the theoretical curve of Eq. 10, we test whether the simulation values of β eff and x eff fall on the theoretical curves. In our experiments, the simulation points agree well with the theoretical curve, which shows the feasibility of our approach.

 2.3 Datasets for Neural Networks

In this study, we used two data sets: pharyngeal nervous systems of C.elegans and P.pacifica [22], all neuronal connections of hermaphrodite and male C.elegans [3]. The groundbreaking work of Cook and Bumbarger et al. provides a wealth data of neural networks, which can help explore the wiring problem from the perspective of structure. We want to couple the information of structure and dynamical function, so as to probe into the complex wiring problem of nematode from the network function and system resilience.

The connection matrix for pharyngeal nervous systems of C.elegans and P.pacifica was derived from the data set used by Bumbarger in his earlier work [22]. The matrix contains neuron connections and weights. The data and its experiment are used for homology comparison of two nematode in our study.

In fact, because of the complex wiring of the nematodes, their neurons have a precise functional grouping: input neuron, intermediate neuron, output neuron. So, for the dataset of hermaphrodite and male C.elegans, instead of focusing on the function of a single node in the nematode network, we focus on the influence of the intermediate neurons on the output neurons.

Connectome adjacency matrices of hermaphrodite and male C.elegans are gap junction, which were provided by Steven J. Cook [3]. Notably, they presented the complete wiring diagram of an animal nervous system for the first time in 2019, including adult C.elegans of both sexes. It is an important milestone in the field of neuron science. The data corresponding to one connectome has the following information: the node, the weight, neuron names in each group. e.g., InterNeurons contains all the neurons that belong to the “InterNeurons” group, neuron names in each subgroup, e.g., InterNeurons_1 contains all the neurons that belong to the “InterNeurons_1” subgroup of the “InterNeurons” group.

 3 Results 3.1 Resilience of <italic>C.elegans</italic> and <italic>P.pacificus</italic>

P.pacificus (Figure 1A) and C.elegans (Figure 1B) have highly similar synaptic structures, but quite different functions in the pharynx [23], we mainly explore the different connective structure and function in pharyngeal neurons.

Comparison of function and structure of synaptic connectivity in C.elegans and P.pacificus. (A) Neuron network for P.pacificus. (B) Neuron network for C.elegans. (C) Activity of 14 homologous synapses between C.elegans and P.pacificus. (D) System performance of C.elegans and P.pacificus while the neurons are disturbed. The red circles represent the system resilience of C.elegans while the homologous synapses are disturbed. The blue circles represent the system resilience of P.pacificus. The darker the color, the greater the effect of the neuronal perturbation.

We measured the functional differences of the connective structure for C.elegans and P.pacificus with x i . As shown in Figure 1C, by exploring the impact of the 14 pharyngeal homologous neurons, we found that the values of M1 for P.pacificus and C.elegans show the greatest difference, which means that M1 is the most functionally different synapse between the two nematodes. M1 is a motor neuron, and their structure centralities of M1 in C.elegans and P.pacificus are quite similar [22]. Our framework can help identify the potential critical neurons, which may be inconsistent with the traditional structural comparison.

Then, with the framework we explored the differences of the predation behavior for the two nematodes through the system resilience index x eff. As shown in Figure 1D, first of all, the activities of C.elegans and P.pacificus were in good agreement with our predicted theoretical curves. Secondly, we were awared that the x eff of the P.pacificus simulation point is much higher than C.elegans. That is to say, the total activity of the P.pacificus is much higher than C.elegans. Nematode feeding strategies differ greatly [24, 25]. P.pacificus have a necromenic association with scarab beetles [22]. Their dauer larvae rest on the insect and resume development after the beetle’s death to feed on microbes on the decaying carcass. P.pacificus can be easily cultured on bacteria [26], however, which is also predatory on other nematodes. In fact, the difference of predation strategy leads to the difference of neuron function. For example, differential but coordinated regulation of pm1 and pm3 in P.pacificus represents motor output that is specific to predatory feeding and does not exist in C.elegans [27, 28]. Actually, the key characteristic of P.pacificus is that it takes predation behavior actively for a long time [29, 30], while the C.elegans, which is preyed on, just eats the food around it. Thus, from the neural level, we can conclude that the system resilience index x eff could be used to reveal the essential differences of the predation behavior for the two nematodes that C.elegans are prey, and P.pacificus is more likely to be predators [31].

Thirdly, we tested 14 pharyngeal neurons, which showed that NSM and I2 are the neuron with the greatest reduction compared with the original simulation point in Figure 1D. It meant that NSM and I2 are the two most important homologous neurons. I2 is an interneuron, which play the role of coordinating corpus and tooth contractions during predation. I2 are more highly connected in P.pacificus and may function as network hubs. Although there are no more detailed studies showing the important function of NSM neurons, the high impact of which on the system resilience in our results is sufficient for further attention of biological experiments. Compared with Bumbarger’s earlier biological experiments [22], the difference is that they concluded that I1 and I2 are the most important candidate neurons from the view of the system structure, but we found NSM and I2 are the two most important homologous neurons while considering the synaptic structural and functional characteristics simultaneously.

In addition, a phenomenon occurs in our experiment from Figure 1D, after neuronal perturbation (node deleting) of a few homozygous neurons, the x eff of P.pacificus decreased, but x eff of C.elegans increased. M4 is one of these special neurons. With its necromenic beetle association, P.pacificus has an intermediate position between the microbivorous C.elegans and true parasites [32, 33]. In other words, C.elegans are prey and P.pacificus is more likely to be predators. As one of the important feeding neurons in the pharynx, M4 is closely related to predation strategy. So we made a bold guess that M4 may degenerate in P.pacificus, but evolve in C.elegans. More details in the Discussion.

 3.2 Quantitative Analysis of Control Model in Hermaphrodite and Male <italic>C.elegans</italic>

From the point of dynamics, the three-layer classification of neurons provides a good basis for the discussion of control. Here we see sensory neurons as input neurons, muscle neurons as output neurons, and other categories as intermediate neurons. We explore the control of the intermediate neurons over the output neurons.

At the beginning, we want to find a class of intermediate neurons that have core control over input neurons. Then, we can further get which of the intermediate neurons has the highest level of control. The control level is described by x i . The interneurons ultimately serve all the neurons in the system, so we also explore the control of the intermediate neurons over all the neurons. In our experiment, aiming at a group of output neurons or all neurons, if a group of interneurons are removed, we should firstly calculate the x i , the one changes greatest, indicating that the interneurons have core control level over the output neurons or all the neurons.

Originally, We found the intermediate neurons group with the greatest control over the output neuron group. For hermaphrodite C.elegans, there are two types of output neurons: SexSpecificCells-Muscle and BodyWallMuscles. As shown in Figure 2A, after perturbation of the OtherEndOrgan intermediate neurons, the x i of the SexSpecificCells-Muscle have dropped mostly, this is to say, the OtherEndOrgan intermediate neurons have the maximum control level over the SexSpecificCells-muscle output neurons. By contrast, after perturbing of the intermediate neurons, the x i of the BodyWallMuscles changes barely, which means that the intermediate neurons have no control over the BodyWallMuscles output neurons. Male C.elegans also has two output neurons: SexSpecificCells-BodyWallMuscles and SexSpecificCells-Muscle. As shown in Figure 2C, when perturbing the Sex-OtherOrgan intermediate neurons, the x i of the SexSpecificCells-BodyWallMuscles and the SexSpecificCells-Muscle both decease sharply, so the sex-OtherOrgan intermediate neurons have the maximum control over the SexSpecificCells-BodyWallMuscles output neurons, and also have the maximum control over the SexSpecificCells-Muscle output neurons. The final result is shown in a three-layer model, as shown in Figures 2B–D.

Control analysis of neuron group in hermaphrodite and male C.elegans. (A) Line chart of the change of the average activity by neuronal perturbation of the each interneurons group in hermaphroditic C.elegans. (pink as the input neuron, blue as the intermediate neuron, and yellow as the output neuron, same as below). (B) Model of the interneurons group with core control in hermaphroditic C.elegans. The groups of interneurons with maximum control were marked with yellow color, same as below. (C) Line chart of the change of the average activity by neuronal perturbation of the each interneurons group in male C.elegans. (D) Model of the interneurons group with core control in male C.elegans.

As shown in Figure 2A, in hermaphrodite C.elegans, no matter removing which group of interneurons, the line at the bottom of the Figure 2A is always stable, which means that the x i of the BodyWallMuscles is steady. When removing the OtherEndOrgan intermediate neurons, what stands out in our results is the x i of whole neurons (green line) decreases mostly. However, in Figure 2C, for the sexOtherorgan intermediate neurons in male C.elegans, all line corresponding to the sexOtherorgan drop sharply, which both have great effect on the neural network and the output neurons. As a result, we conclude that although these intermediate neurons have no direct control over the output neurons, they still play an important role in the overall functional mechanism of C.elegans. Hence, even though there is no decrease in the activity of the output neurons after the ablation experiment, the average activity of all neurons decreases significantly. In other words, the control mechanism over the output neurons and the whole neurons are not congruent.

Ultimately, we try to identify the intermediate neurons group with the greatest control level. We measure the neuron disturbance by degree in the intermediate group with core control, and calculate the effect of neuron control by the change of the mean value of x i . In Figure 3, the neurons in Y-coordinate are arranged in degrees. As can be seen from the Figure 3A, in the case of male C.elegans’ sexOtherorgan intermediate group of neurons, the value decreases mostly after the removal of CEPshVR neurons. So CEPshVR has the maximum control over the SBmuscle output neurons group. As shown in Figure 3C, the value in male C.elegans’ sexOtherorgan intermediate neurons decreases sharply with the disturbance of Int neurons. Consequently, the Int neurons have the maximum control over the Smuscle output neurons group. Similarly, as shown in Figure 4A, in the case of hermaphrodite C.elegans’ OtherEndOrgan intermediate neurons, the HSNR neurons have the maximum control over the sexspecificcells-muscle output neurons group. Briefly, for most of the interneuron groups, the higher the degree, the bigger circle in Figures 3, 4, the greater the controlling impact. To illustrate more clearly, we mark in Figures 3B,D, 4B, The gray area in the middle displays the specific neurons of the interneuron group, where the neurons of the core control level are marked in yellow.

Control analysis of single neuron in male C.elegans. (A) Activity line chart of SBmuscle, after neuronal perturbation of the Sexotherorgan intermediate neuronal group in the male C.elegans by degree. Degree refers to the degree of the node in this nematode network, we calculated degree and sorted them. The order is shown in figure: The size of the green circle represents the size of degree, same as below. (B) Control model of SBmuscle by the Interneuron group of Sexotherorgan. The neurons with maximum control are shown in yellow color, same as below. (C) Activity line chart of Smuscle, after neuronal perturbation of the Sexotherorgan intermediate neuronal group in the male C.elegans by degree. (D) Control model of Smuscle by the Interneuron group of Sexotherorgan. (SensoryNeurons:Sen; SexSpecificCells-SensoryNeuron:Sex-Sen; SexSpecificCells -InterNeurons:Sex-Inter; SexSpecificCells -OtherEndorgans: Sex- OtherEndOrgans. Same as the others).

Model exploration of single neuron hermaphroditism in C.elegans.(A) Activity line chart of SexSpecficalCells-muscle, after neuronal perturbation of the SexSpecficalCells-motor intermediate neuronal group in the hermaphroditism C.elegans by degree. (B) Control model of SexSpecficalCells-muscle by the Interneuron group of SexSpecficalCells-motor.

 4 Discussion

In our framework, the coupling method of network structure and dynamical interactions enables us to establish a dynamical model of neuron networks to explore the system performance. To quantitatively solve the high-dimension rate equations, we provide a way in which multi-dimensional features between neurons can be reduced into one-dimension equation. Our approach can help the realization of the mapping from structure to function, and the quantitative measurement from function to control for neuron systems.

Excepting this, we noticed that x eff of the P.pacificus simulation point is much higher than C.elegans. With its necromenic beetle association, P.pacificus has an intermediate position between the microbivorous C.elegans and true parasites [34, 35]. The omnivorous feeder P.pacificus should have the most complex metabolic pathways for nutrition and protection against defense and prey, comparing with the microbivorous C.elegans [36, 37]. So we infer that the essence of this phenomenon is the distinction of predation behavior caused the difference of synaptic connection, which changes the structure of nematode network.

Furthermore, we hope that our research could be used to illustrate evolution and degeneration of neuron: M4 most likely degenerates in C.elegans, evolve in P.pacificus. Combining the differences in the two nematode predation strategies, and the functions of M4, our evidence is as follows: M4 are known to be one of the major excitatory neurons in C.elegans, which is required for posterior isthmus peristalsis [38–40], it meant that M4 plays a crucial role in the act of eating or swallowing in C.elegans. However, the ablation experiments made by Edelman and Garry have shown that MC are the only neurons required for rapid eating in P.pacificus. But the ablation of other cholinergic pharyngeal neurons, such as M2S and M4, only slightly reduced feeding rates [41, 42]. In addition, Trojanowski’s previous results demonstrate that this robust and evolutionarily adaptable network is highly degenerate at both the neural and genetic levels, the same behavior can be stimulated by multiple neurons and different types of receptors [41]. Avery and Horvitz even found out, in the C.elegans gene ced-3 Mutant, probably MSpaaaaap (the sister of M4), can sometimes take over M4’s function. Using a laser microbeam to kill cells, they found that one of the extra cells in ted-3 worms, tentatively identified as MSpaaaaap, could become a functioning M4 neuron about half of the time, although it rarely or never fully replaced M4 [38]. Given the above, we can not be absolutely sure that M4 most likely degenerates in C.elegans, evolve in P.pacificus, but we hope our theoretical framework can help explore the function of neurons in the evolution and degeneration process.

To illustrate the power of our framework, we analyzed the control model of hermaphroditic and male C.elegans. Based on the data of Cook et al [3], our experiment abstracted the neuron connection group from the real network to the theory level by a nonlinear dynamical model, and obtained a measurable control index of the intermediate neurons to the output neurons. The control between different functional groups need to take into account both structural and functional characteristics. It means that our framework may help design more targeted and accurate biological ablation experiments. To understand vividly, we mapped the control model of the three-layer neurons, which is accurate to the individual neurons [43, 44].

In sum, our theoretical research perfects the exploration of the functional structure of nematodes, and breaks the limitation of analyzing the function and control ability from the structure. Comparing with the linear control model to predict motor neuron controllability, Yan et al. use network control principles, we emphasize the growth rate of the node itself and the influence of its neighbors on the node [18, 45] Importantly, this is a more efficient way to quickly target groups of neurons, even individual neurons or synapses that have core controlling power. In this view, our theoretical approach contributes to exploring the importance and the evolution mechanism of the individual neurons, and proves how to calculate the control level with the dynamics. These observations are not limited to nematode, but are also applicable to connectome model organism with complete data and clear background. However, the real interactions of nematode connection network are much more complicated than our dynamical model [46, 47]. We use this model for theoretical exploration, there may be some difference between our parameter setting and the real system. We need to refer to more biological experiments for further parameter setting and to use more complex dynamical equations or generative models to capture the interactions in our next work.

 Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://doi.org/10.1016/j.cell.2012.12.013 and https://doi.org//10.1038/s41586-019-1352-7.

 Author Contributions

DD and WX performed the analysis. WX validated the analysis and drafted the manuscript. DD and SS reviewed the manuscript. WX designed this study. All authors have read and approved the content of the manuscript.

 Funding

This work was supported by the National Natural Science Foundation of China (projects no. 72071153, the Natural Science Foundation of Shaanxi Province (project n.2020JM-486), the China Postdoctoral Science Foundation (project n.2017M613336).

 Conflict of Interest

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.

 Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

 References Dirven L Petersen M Petersen MA Aaronson NK Chie W-C Conroy T Development and Psychometric Evaluation of an Item Bank for Computerized Adaptive Testing of the Eortc Insomnia Dimension in Cancer Patients (Eortc Cat-Sl). Appl Res Qual Life (2019) 16:827–44. 10.1007/s11482-019-09799-w Bargmann CI . Genetic and Cellular Analysis of Behavior in C. Elegans. Annu Rev Neurosci (1993) 16:47–71. 10.1146/annurev.ne.16.030193.000403 Cook SJ Jarrell TA Brittin CA Wang Y Bloniarz AE Yakovlev MA Whole-animal Connectomes of Both Caenorhabditis elegans Sexes. Nature (2019) 571:63–71. 10.1038/s41586-019-1352-7 Sarma GP Lee CW Portegys T Ghayoomie V Jacobs T Alicea B Openworm: Overview and Recent Advances in Integrative Biological Simulation of Caenorhabditis elegans . Philos Trans R Soc Lond B Biol Sci (2018) 373(1758):20170382. 10.1098/rstb.2017.0382 Witvliet D Mulcahy B Mitchell JK Meirovitch Y Berger DR Wu Y Connectomes Across Development Reveal Principles of Brain Maturation. Nature (2021) 596:257–61. 10.1038/s41586-021-03778-8 Morita S Oshio K-i. Osana Y Funabashi Y Oka K Kawamura K . Geometrical Structure of the Neuronal Network of Caenorhabditis elegans . Physica A: Stat Mech its Appl (2001) 298:553–61. 10.1016/s0378-4371(01)00266-7 Alicea B . Raising the Connectome: The Emergence of Neuronal Activity and Behavior in Caenorhabditis elegans . Front Cel Neurosci (2020) 14(524791):524791. 10.3389/fncel.2020.524791 Duan D Lv C Si S Wang Z Li D Gao J Universal Behavior of Cascading Failures in Interdependent Networks. Proc Natl Acad Sci USA (2019) 116:22452–7. 10.1073/pnas.1904421116 Chen BL Hall DH Chklovskii DB . Wiring Optimization Can Relate Neuronal Structure and Function. Proc Natl Acad Sci (2006) 103(12):4723–8. 10.1073/pnas.0506806103 Lav R Varshney LB Chen EP Hall DH Chklovskii DB . Structural Properties of the Caenorhabditis elegans Neuronal Network. PLoS Comput Biol (2011) 7(2):e1001066. 10.1371/journal.pcbi.1001066 Chalfie M Sulston J White J Southgate E Thomson J Brenner S . The Neural Circuit for Touch Sensitivity in Caenorhabditis elegans . J Neurosci (1985) 5(4):956–64. 10.1523/jneurosci.05-04-00956.1985 Costa AC Ahamed T Stephens GJ . Adaptive, Locally Linear Models of Complex Dynamics. Proc Natl Acad Sci U S A (2019) 116(5):1501–10. 10.1073/pnas.1813476116 Baltzley MJ Gaudry Q Kristan WB . Species-specific Behavioral Patterns Correlate with Differences in Synaptic Connections Between Homologous Mechanosensory Neurons. J Comp Physiol A (2010) 196(3):181–97. 10.1007/s00359-010-0503-y Liu Y Sanhedrai H Dong GG Shekhtman LM Wang F Buldyrev SV Efficient Network Immunization Under Limited Knowledge. Natl Sci Rev (2020) 10.1093/nsr/nwaa229 Dong G Wang F Shekhtman L Danziger M Fan J Du R Optimal Resilience of Modular Interacting Networks. Proc Natl Acad Sci 118(22):e1922831118. 10.1073/pnas.1922831118 Morone F Del Ferraro G Makse HA . The K-Core as a Predictor of Structural Collapse in Mutualistic Ecosystems. Nat Phys (2018) 15:95–102. 10.1038/s41567-018-0304-8 Morone F Makse HA . Symmetry Group Factorization Reveals the Structure-Function Relation in the Neural Connectome of Caenorhabditis elegans . Nat Commun (2019) 10(1):4961. 10.1038/s41467-019-12675-8 Yan G Vértes PE Towlson EK Chew YL Walker DS Schafer WR Network Control Principles Predict Neuron Function in the Caenorhabditis elegans Connectome. Nature (2017) 550(7677):519–23. 10.1038/nature24056 Blaxter ML De Ley P Garey JR Liu LX Scheldeman P Vierstraete A A Molecular Evolutionary Framework for the Phylum Nematoda. Nature (1998) 392:71–5. 10.1038/32160 Duan D-L Tao C Wu X-X Wu C Shu-bin S Gen-qing B . Identification of Unstable Individuals in Dynamic Networks. Chin Phys B (2021) 30(9):090501. 10.1088/1674-1056/abe92f Gao J Barzel B Al B . Universal Resilience Patterns in Complex Networks. Nature (2016) 530(7590):307–12. 10.1038/nature16948 Bumbarger DJ Riebesell M Rödelsperger C Sommer RJ Sommer RJ . System-wide Rewiring Underlies Behavioral Differences in Predatory and Bacterial-Feeding Nematodes. Cell (2013) 152(1):109–19. 10.1016/j.cell.2012.12.013 Franks CJ Holden-Dye L Bull K Luedtke S Walker RJ . Anatomy, Physiology and Pharmacology of Caenorhabditis elegans Pharynx: A Model to Define Gene Function in a Simple Neural System. Invert Neurosci (2006) 6(3):105–22. 10.1007/s10158-006-0023-1 Mörck C Axäng C Pilon M . A Genetic Analysis of Axon Guidance in the C Elegans Pharynx. Dev Biol (1200) 260:158–75. 10.1016/s0012-1606(03)00238-0 Pilon M . Developmental Genetics of the Caenorhabditis Eleganspharynx. Wires Dev Biol (2014) 3:263–80. 10.1002/wdev.139 Gruninger TR Gualberto DG LeBoeuf B Garcia L . Integration of Male Mating and Feeding Behaviors in Caenorhabditis elegans . J Neurosci (2006) 26:169–79. 10.1523/jneurosci.3364-05.2006 Avery L You YJ . C. elegans Feeding. WormBook (2012) 1–23:1–23. 10.1895/wormbook.1.150.1 Avery L Horvitzt HR . Pharyngeal Pumping Continues After Laser Killing of the Pharyngeal Nervous System of C. elegans . Neuron (1989) 3:473–85. 10.1016/0896-6273(89)90206-7 Pervez Rashid . Attraction of Allodorylaimus americanus (Nematoda: Dorylaimida) towards different prey nematodes and factors influencing this attraction. Archives of Phytopathology and Plant Protection 42:344–351. 10.1080/03235400601070421 (2009). Lichtman JW Livet J Sanes JR . A Technicolour Approach to the Connectome. Nat Rev Neurosci (2008) 9(6):417–22. 10.1038/nrn2391 Haspel G O'Donovan MJ Hart AC . Motoneurons Dedicated to Either Forward or Backward Locomotion in the Nematode Caenorhabditis elegans . J Neurosci (2010) 30(33):11151–6. 10.1523/jneurosci.2244-10.2010 Albertson D Thomson JN . The Pharynx of Caenorhabditis elegans . Phil Trans R Soc Lond B (1976) 275:299–325. 10.1098/rstb.1976.0085 Sherman D Harel D . Deciphering the Underlying Mechanisms of the Pharyngeal Motions in Caenorhabditis elegans . Ithaca, NewYork: Neurons and Cognition (2019). Mapes CJ . Structure and Function in the Nematode Pharynx. Parasitology (1965) 55:583–94. 10.1017/s0031182000069286 Mayer WE Herrmann M Sommer RJ . Phylogeny of the Nematode Genus Pristionchus and Implications for Biodiversity, Biogeography and the Evolution of Hermaphroditism. BMC Evol Biol (2007) 7:104. 10.1186/1471-2148-7-104 Dieterich C Clifton SW Schuster LN Chinwalla A Delehaunty K Dinkelacker I The Pristionchus Pacificus Genome Provides a Unique Perspective on Nematode Lifestyle and Parasitism. Nat Genet (2008) 40:1193–8. 10.1038/ng.227 Stern S Kirst C Bargmann CI . Neuromodulatory Control of Long-Term Behavioral Patterns and Individuality across Development. Cell (2017) 171:1649–62. 10.1016/j.cell.2017.10.041 Avery L Horvitz HR . A Cell that Dies During Wild-type C. elegans Development Can Function as a Neuron in a Ced-3 Mutant. Cell (1987) 51(6):1071–8. 10.1016/0092-8674(87)90593-9 Avery L Horvitzt HR . Pharyngeal Pumping Continues After Laser Killing of the Pharyngeal Nervous System of C. elegans . Neuron (1989) 3(4):473–85. 10.1016/0896-6273(89)90206-7 Wicks S Rankin C . Integration of Mechanosensory Stimuli in Caenorhabditis elegans . J Neurosci (1995) 15(3):2434–44. 10.1523/jneurosci.15-03-02434.1995 Trojanowski NF Padovan-Merhar O Raizen DM Fang-Yen C . Neural and Genetic Degeneracy Underlies Caenorhabditis elegans Feeding Behavior. J Neurophysiol (2014) 112(4):951–61. 10.1152/jn.00150.2014 Rudiger B Colomb J Pankratz B Schröck A Stocker RF Genetic Dissection of Neural Circuit Anatomy Underlying Feeding Behavior Indrosophila: Distinct Classes Ofhugin-Expressing Neurons. J Comp Neurol (2007) 502(5):848–56. 10.1002/cne.21342 Whalen AJ Brennan SN Sauer TD Schiff SJ . Observability and Controllability of Nonlinear Networks: The Role of Symmetry. Phys Rev X (2015) 5:011005. 10.1103/PhysRevX.5.011005 Tsalik EL Hobert O . Functional Mapping of Neurons that Control Locomotory Behavior in Caenorhabditis Elegans. J Neurobiol (2003) 56(2):178–97. 10.1002/neu.10245 Gao J Liu YY D'Souza RM Barabási AL . Target Control of Complex Networks. Nat Commun (2014) 5:5415. 10.1038/ncomms6415 Stephens GJ Johnson-Kerner B Bialek W Ryu WS . Dimensionality and Dynamics in the Behavior of C. elegans . Plos Comput Biol (2008) 4:e1000028. 10.1371/journal.pcbi.1000028 Zhen M Samuel AD Samuel T . C. elegans Locomotion: Small Circuits, Complex Functions. Curr Opin Neurobiol (2015) 33:117–26. 10.1016/j.conb.2015.03.009