Download Central Pattern Generators: Understanding the Behavior of Coupled Oscillators in Biology and more Exams Architecture in PDF only on Docsity!
Central Pattern Generators
SHARON CROOK and AVIS COHEN
8.1 Introduction
Many organisms exhibit repetitive or oscillatory patterns of muscle activity that produce rhythmic movements such as locomotion, breathing, chewing and scratching. Examples include the escape swimming of the mollusc Tritonia diomedia , the digestive rhythms of the lobster, the undulatory swimming movements of the fish or the lamprey, the stepping movements of the cockroach, the rapid wing motion of the locust during flight, and the more complicated locomotion of a quadruped mammal such as the domestic cat. The neuronal circuits that give rise to the patterns of muscle contractions which produce these movements are referred to as central pattern generators , or CPGs. Various experimental preparations in which the CPG is isolated from external influence demonstrate that these circuits require no external control for the generation of temporal sequences of rhythmic activity. However, these animals move through the world in an adaptive manner where the same motorneurons are involved in the production of a variety of rhythmic behaviors. Thus, many CPGs are capable of producing multiple patterns of activity in the intact behaving animal (Getting 1989). The ability to switch between different motor behaviors and blend different rhythms relies on feedback from proprioceptors and influence from higher centers of the nervous system; therefore, it is most appropriate to view every CPG as one piece of a distributed control system (Cohen 1992). One would like to understand how the neurons in a CPG interact and influence one another, how the underlying circuitry of the network produces the collective behavior of the cells, what mechanisms might allow the network to switch among various patterns of
131
132 C h apter 8. Central Pattern Generators
activity, and whether the oscillatory patterns are due primarily to the activity of individual intrinsically oscillatory neurons or to oscillations that are a product of the entire network. The number of cells composing a network that functions as a CPG often determines the manner in which the CPG is studied and the choice of a modeling strategy. Some CPG circuits are anatomically localized and contain a small number of neurons. This occurs most often in CPGs that produce rhythmic behaviors in invertebrates. In these small networks, neurons can be individually identified from animal to animal, permitting detailed circuit descriptions that include cellular and synaptic properties. In contrast to these invertebrate CPGs, there are possibly millions of neurons involved in the production of rhythmic patterns of motor activity in most vertebrates (Murray 1989). In this case, modelers often categorize the neurons into classes that share similar properties so that large networks can be simulated by relatively few cell types (Getting 1989). The small localized CPGs that occur in invertebrate preparations make it possible to study the relationship between the emergent collective behavior of the biological network and the network’s underlying circuitry (Getting 1989). The dynamical properties of many invertebrate CPGs have been analyzed using such techniques as experimental manipula- tions of cellular, synaptic, and connectivity properties, detailed simulations of the cell in- teractions within the network, and analytical studies of equations that might describe the network dynamics. For example, Getting created a network simulation of the escape swim- ming rhythm of the mollusc Tritonia diomedia (Getting 1989). This simulation relies on a compartmental model of the network cells with appropriate passive membrane properties, repetitive firing characteristics, and synaptic actions. In addition, the input to the model corresponds to the normal sensory activation of the actual CPG. Some of the properties of Getting’s model are demonstrated in the GENESIS simulation Tritonia. Another example of an invertebrate CPG that has been studied and modeled extensively is the lobster stomato- gastric ganglion. This region contains the neurons that are involved in the generation of the slow rhythm that fires the muscle contractions of the lobster gastric mill and also those that generate the rhythm that controls the muscles of the pyloric region of the lobster stomach (Shepherd 1994). Experimentation with this system has shown that even when a detailed study of the network circuitry provides a qualitative description accounting for the presence of a given motor pattern, there is often no precise explanation for the mechanisms that con- trol the frequency, duration, and phase relations of the motor pattern (Marder and Meyrand 1989). Studies of the invertebrate CPGs mentioned above show that the generation of these rhythms is a complicated process involving the influence of multiple neurotransmitters and modulators that modify the output of the circuit (Marder and Meyrand 1989). Due to the large number of neurons present in most vertebrate CPG circuits, mod- els of CPGs in vertebrates often involve simplified mathematical representations where a single oscillator may represent many neurons. For example, most models of mammalian locomotion attempt to create an oscillatory network that can account for the production of the alternating flexor-extensor activity responsible for limb coordination during locomo-
134 C h apter 8. Central Pattern Generators
8.2 Two-Neuron Oscillators
In this section, we consider models that mimic the behavior of a system of two coupled os- cillators in an attempt to understand some of the ways in which biological oscillators may influence one another. Each oscillator can represent a single neuron or a network of cells that collectively function as an oscillator. We first present phase equation models that do not depend on the oscillator structure so that the model results apply to both single cell and network oscillators. We show that these models are simple enough to be analyzed mathe- matically yet complex enough to capture some of the underlying principles that govern the behavior of a two-oscillator network. Next, we briefly mention a modeling option that uses higher-order systems of equations and incorporates more information about the oscillator structure.
8.2.1 Phase Equation Model of Coupled Oscillators
First consider a general mathematical model due to Rand, Cohen, and Holmes (1988) for a network of two oscillators where each is treated as a simple biological oscillator, ignoring the structure of the oscillation and the mechanisms that produce it. In actuality, the behavior of each oscillatory neuron or network oscillator is determined by a multitude of parameters that can be used to represent the state of the oscillator at any given time. Due to the cyclic behavior of each oscillator, if we draw an orbit in the parameter space that shows how the parameters change with time, the oscillator eventually returns to the same state. This type of orbit is known as a limit cycle. In addition, for small perturbations away from the orbit in the parameter space, the orbit of the oscillator returns to this cycle; that is, the orbit is locally asymptotically stable as shown in Fig. 8.1. These assumptions allow us to test assertions regarding the coordinating system while knowing little about the individual oscillators.
Figure 8.1 Locally asymptotically stable limit cycle in two-dimensional space.
8.2. Two-Neuron Oscillators 135
Since we assume that each of the two oscillators in this model can be represented by a structurally stable dynamical system that exhibits a locally asymptotically stable limit cycle, the behavior of each can be represented by a single variable, θ i
t^^ for^ i^^ ^1 ^2 ^ , which specifies the position of the oscillator around its limit cycle at time t ; that is, θ i
t^ specifies the phase of the limit cycle. In addition, we rescale θ i
t^^ so that it flows uniformly around the limit cycle taking values from 0 to 2π radians over one cycle. Thus, θ i
t^^ is proportional to the fraction of the period that has elapsed and the behavior of each oscillator is characterized by the differential equation
θ˙ i t ω i (8.1)
where ω i is the constant frequency of the oscillator, and 2π ω i is the period. We use mod- ular arithmetic so that θ i
t^^ always lies between 0 and 2π, and the solution to Eq. 8.1 is θ i
t^^
ω it θ i
0 ^
mod 2 π (8.2)
where θ i
0 is the initial value of^ θ i. When two such oscillators are coupled, we obtain a system of equations
θ˙ 1 t ω 1 h 12 θ 1 θ 2 (8.3) θ˙ 2 t ω 2 h 21 θ 2 θ 1 (8.4)
where hi j
θ i θ^ j represents the coupling effect of the^ j th oscillator on the^ i th oscillator. This coupling term must be 2π periodic since we would like the rate of change to depend only on the oscillator phase and not on the number of cycles that have already occurred or the amplitude of the oscillation. Since the behavior of each cell is described using only the phase, the coupling also depends only upon the phase. It is often convenient to define a quantity φ
t θ 1
t^^ ^ θ 2
t^^ (8.5)
which represents the difference between the phases or the phase lag of oscillator 2 relative to oscillator 1. Note that like θ i , φ
t 0 2 π. Combining Eqs. 8.3 and 8.4, we obtain
φ˙ t θ˙ 1 t θ˙ 2 t (8.6)
ω 1 ω 2
h 12
θ 1 θ 2 ^ h 21
θ 2 θ 1 (8.7)
Much of the analytical work to date on coupled oscillators deals with interactions that de- pend only on the difference between the phases. Consider the simplest case where hi j depends only on the phase lag, and in addition hi j
θ i θ^ j ^ 0 when the phase lag between the i th and j th oscillators is zero. This is known as diffusive coupling. For example, if we let hi j ai j sin
θ (^) j θ i , we obtain
φ˙ t ω 1 ω 2 a 12 a 21 sin φ t (8.8)
8.2. Two-Neuron Oscillators 137
R1 to 5 and making sure that the toggle for that connection is in the excitatory mode. Set the strength of the connection from R1 onto L1 to 5 in a similar manner. After running the simulation, the graphs reveal that the two oscillating cells are drifting with respect to one another. In other words, the coupling is too weak to cause phase-locked behavior. Now gradually increase the coupling strength between the two cells by increments of 5,
clicking on the and buttons after each change in order to observe the behavior
of the system. As the coupling strength is increased, eventually the system exhibits be- havior that converges to phase-locked oscillations where the cell with the greater intrinsic frequency, namely, R1, leads. This causes the frequency of oscillation of the system to be greater than the intrinsic frequency of L1. Now repeat the process outlined above for inhibitory connections, beginning with weak mutually inhibitory coupling and gradually increasing the coupling strength. In this case, stronger coupling strengths are required for phase-locked behavior due to the nature of the inhibitory connections in the simulation en-
vironment as explained in the text of the help selection under . During
the phase-locked oscillations, the cell with the slower intrinsic oscillation will lead, causing the frequency of the system to be smaller than the intrinsic frequency of R1.
8.2.3 Initial Conditions
In the qualitative analysis of any system of differential equations, it is important to consider the role of the initial conditions since the system behavior depends on these initial values. Consider a CPG simulation architecture that demonstrates the crucial role of the initial conditions for the two-oscillator model presented in Sec. 8.2.1. Begin by changing the current injections to L1 and R1 so that they are identical, say, 0 00015 μA. In this case the two oscillators have identical intrinsic frequencies; that is, ω 1 ω 2. Since the frequency difference ω 1 ω 2 is zero, Eq. 8.9 indicates that there are two solutions that provide phase- locked behavior, φ 0 and φ π. Thus, the two oscillators may be phase-locked with no phase lag so that they demonstrate synchronous behavior, or they may be phase-locked 180 degrees out of phase. We can obtain these two types of behavior by varying the initial conditions of the system. First, set the delays for both cells to zero and make sure that the durations are set to 200 msec. Begin with symmetric mutually excitatory coupling with the coupling strength of each connection set to 20. Run the simulation and observe that the two cells begin firing in phase due to their identical initial values and continue to fire in phase due to the symmetry of the network. Next change the delay of R1 to 15 msec. In this case, we are effectively changing the initial conditions of the system so that we begin with two cells that are firing out of phase at time t 15 msec. As predicted for the case where the intrinsic frequencies of the two oscillators are identical, the cells still exhibit phase-locked motion; however, they fire out of phase. This demonstrates how important initial conditions may be for the production of different temporal patterns of activity in these simplified models. However,
138 C h apter 8. Central Pattern Generators
the initial conditions do not play such a crucial role in biological system where noise and inherent differences prevent any two oscillators from having identical frequencies.
8.2.4 Synaptic Coupling
Naturally, the phase equation model described above is not sufficient for predicting or study- ing all possible types of behavior exhibited by two coupled neurons or network oscillators. In particular, the choice of diffusive coupling may be limiting. Although diffusive coupling probably models the behavior of electrical coupling between multiple cells accurately, in general we would not expect the synaptic interactions between two neurons to depend ex- clusively on the difference between the phases. For example, with diffusive coupling a phase difference of zero results in a coupling term with a value of zero so that two neurons that are exhibiting identical behavior have no influence on one another. Because this does not seem biologically plausible for synaptic connections between neurons, Ermentrout and Kopell (1990) consider a model that differs slightly from the previous one and contains coupling terms that behave more like synaptic coupling between cells. The assumptions for this coupling hold true for Hodgkin-Huxley-like neural models in which the coupling is due to voltage only (Ermentrout and Kopell 1990). In the following description of the model, we assume that the oscillators are identical for ease of analysis. The equations for this model are
θ˙ 1 ω 1 p θ 2 r θ 1 (8.10) θ˙ 2 ω 2 p θ 1 r θ 2 (8.11)
which are simply Eqs. 8.3 and 8.4 with hi j
θ i θ^ j p
θ (^) j r
θ i , where p is a periodic smooth pulse function and r plays a role that is analogous to that of a phase response curve. A phase response curve, or PRC, for a given oscillator and a given brief stimulus is determined experimentally by stimulating the oscillator and waiting until the system relaxes back to its oscillation with a shift in phase. The PRC can be represented by a function ˆ r
θ that gives the phase shift and depends upon the phase θ at which the stimulus is adminis- tered. Consider the situation where the relaxation time is short relative to the time between stimuli τ. Let θ k^ be the phase just before the k th stimulus so that we obtain a sequence of phases θ k^^ that are related by the difference equation
θ k ^1 ωτ θ k^ r ˆ
θ k^^ (8.12)
where ω is the natural frequency of the oscillator. We can rewrite Eq. 8.12 as the differential equation θ˙ ω δ t mod τ r ˆ θ (8.13)
where δ is the Dirac delta function that provides a pulse stimulus at intervals a time τ apart. According to Ermentrout, since real stimuli are not instantaneous, it is reasonable
140 C h apter 8. Central Pattern Generators
the simulation is more complex than the behavior of the phase equation model of coupled oscillators discussed above; however, some basic interactions can be understood in the con- text of the phase equation model.
8.3 Four-Neuron Oscillators
Now that we have discussed some of the basic principles that govern the behavior of two coupled oscillators, we would like to apply these ideas to slightly larger networks. As previously mentioned, larger networks are more difficult to analyze mathematically since they require more equations, and thus we will eventually rely on simulations in our attempt to understand their behavior. In this section we examine chains of four oscillators with nearest-neighbor coupling as well as more complex architectures with various types of cou- pling among four oscillatory cells. We attempt to explore some of the basic ideas used to formulate the vertebrate CPG models mentioned in Sec. 8.1. In particular, the discussion of chains of four coupled oscillators that follows provides insight into various models for lam- prey locomotion and the locomotion of some types of fish. These models are inspired by the fact that the waves of muscle contractions that produce locomotion in these organisms are induced by periodic bursts of activity in the ventral roots that emerge at each segment of the spinal cord. The phase lag in activity between segments is proportional to the dis- tance between the points, indicating a constant speed traveling wave of contractions (Rand et al. 1988, Kopell and Ermentrout 1986). Following the discussion of chains of coupled oscillators, we construct networks capable of producing patterns of activity that mimic var- ious gaits for tetrapod locomotion. In these gait simulations, the phase relations among the various oscillators in the CPG model are important since they represent the phase relations among the limbs during locomotion.
8.3.1 Chains of Coupled Oscillators
We begin by extending the concepts outlined in Sec. 8.2.1 in order to develop a phase equation model for a chain of four coupled oscillators. The following treatment appears in Rand et al. (1988) in a more general form for chains of n oscillators. We would like the chain to have only nearest-neighbor coupling with no long-range connections. Thus, Eqs. 8.3 and 8.4 are extended to a system of four differential equations of the form
θ˙ 1 t ω 1 h 12 θ 1 θ 2 (8.14) θ˙ 2 t ω 2 h 21 θ 2 θ 1 h 23 θ 2 θ 3 (8.15) θ˙ 3 t ^ ω 3 h 32 θ 3 θ 2 ^ h 34 θ 3 θ 4 (8.16) θ˙ 4 t ω 4 h 43 θ 4 θ 3 (8.17)
8.3. Four-Neuron Oscillators 141
As before, we assume diffusive coupling and let hi j ai j sin
θ j θ i . However, we use
a uniform coupling strength ai j a for ease of analysis. Again, we define the quantity
φ i
t^^ ^ θ i
t^^ ^ θ i 1
t^^ (8.18)
to represent the phase lag between oscillators for i 1 2 3 . So, our system of equations
can be represented by the vector equation
φ˙ t Ω AS t (8.19)
where
φ
t
φ 1
t^
φ 2
t^
φ 3
t^
(8.20)
(8.21)
A a
2 1 0
(8.22)
and
S
sin
φ 1
t^^
sin
φ 2
t^^
sin
φ 3
t^^
(8.23)
As in the previous sections, we are interested in 1 : 1 phase locked motion where the phase
lags remain constant; therefore, we would like to find solutions to the equation ˙φ 0. Note
that ˙φ 0 when S A 1 Ω, and we know that
A ^1
^1
4 a
3 2 1 2 4 2 1 2 3
(8.24)
If we assume that there is a smooth gradient in the intrinsic frequency of oscillation along
the cord, then the frequency difference along the chain is constant, and ω 1 ω 2 ω 2 ω 3
ω 3 ω 4 c. In this case,
Ω c
1 1 1
(8.25)
8.3. Four-Neuron Oscillators 143
L1 R
L2 R
Figure 8.2 Simulation architecture for a chain of four oscillators with nearest-neighbor coupling. All con- nections are excitatory with identical strengths.
end conditions. This is the behavior predicted by the model in the case when there is no difference in oscillator frequencies along the chain so that c 0.
Next change the injection currents to the cells causing the cells to have varying intrinsic frequencies with a constant difference in frequency along the chain. For example, you might use values of 0 00035, 0 0003, 0 00025, and 0 0002 μA as the injection inputs along the chain. In this case, the model predicts a traveling wave of activation provided that the value of the connection strength is large enough. Begin with very small values for the identical connection strengths and gradually increase the value until you see a wave of activation. It may be difficult to see the wave since each one begins before the previous wave has reached the end of the chain.
It is also possible to obtain a wave of activation in a chain of neurons simply by causing the first cell in the chain to fire. Each cell activates the next so that the activation moves down the length of the chain. To simulate this type of activity, set the injection input to the soma of the first cell in the chain to 0 0002 μA for a duration of several milliseconds, and set all other possible inputs to zero. This will cause a single wave of activation; however, this model requires outside control for continuous behavior since multiple waves require repeated inputs to the first oscillator at the desired temporal intervals. For this reason, the behavior is less robust than that of the model described above which uses intrinsic oscillators to create continuous waves of behavior in the absence of any external influence to the CPG.
144 C h apter 8. Central Pattern Generators
8.3.3 Modeling Gaits
During tetrapod locomotion, interlimb coordination causes different combinations of limbs to be on or off the ground at the same time as the gait of an animal varies. Each limb moves through a step cycle that consists of a swing phase and a stance phase. During the swing phase, the limb is lifted and brought forward, mostly due to the effort of the flexor muscles. During the stance phase, extensor activity is dominant and provides force to thrust the animal forward (Grillner 1981). The most common mode of coordination among limbs is strict alternation between the two limbs of the same girdle; that is, the movements of the two hindlimbs alternate with each other, and the movements of the two forelimbs alternate with each other. These are known as the alternating gaits which include the walk, pace and trot and are seen in a variety of organisms from the millipede to humans. The differences among these three gaits depend on the timing between the forelimbs and hindlimbs. We refer to all other gaits as non-alternating gaits which include the gallop, canter, half bound, and leaping gaits. Most animals prefer to use different gaits at different velocities, but the gait is not directly linked to the speed (Grillner 1981). We want to model the patterns of limb coordination required for some of the various gaits mentioned above. In these simulations we let one intrinsic oscillator represent the population of cells that controls the cyclic flexor-extensor activity of a single limb. To
begin, click on the button in the main control panel to return to the default
architecture provided with the CPG simulation. In this default network architecture, the cell labeled L1 represents the left forelimb, R1 represents the right forelimb, L2 represents the left hindlimb, and R2 represents the right hindlimb. Note that all connection strengths are
set to a value of 100. Hit and observe the pattern of activation that emerges from the
network circuitry. As shown in Fig. 8.3, this pattern of activity is the sequence observed in a walking gait. The lateral symmetric inhibitory connections between the two forelimbs and between the two hindlimbs encourage the alternating phase-locked behavior required for all alternating gaits. The excitatory connections help maintain a robust pattern of activation where each oscillator stimulates the next one in the sequence. The remaining inhibitory connections discourage the oscillators from becoming active out of turn. Also, note that the delays are set to values that initialize the simulation with the sequence of activation required for a walking gait. Experimenting with the delays reveals that a variety of initial conditions will lead to this same pattern of behavior although the amount of time required for convergence to phase-locked oscillations varies for different initial conditions. It is important to make sure that the delays are not all identical since beginning the oscillators in phase results in completely synchronous firing due to the symmetry of the network. Now, alter the architecture so that the vertical connections between L1 and L2 and between R1 and R2 are inhibitory, and make all of the cross-connections excitatory. In addition, change the delays so that L1 and R2 have delays of zero, and L2 and R1 have
delays of 15. Click on and then run the simulation. This pattern of activity mimics
146 C h apter 8. Central Pattern Generators
for changing gaits.
8.4 Summary
We are not suggesting that the models and simulations described in this chapter accurately describe how CPGs produce observed patterns of behavior. It is certain that the actual mechanisms are much more complicated and that feedback from proprioceptors and influ- ence from higher centers of the nervous system play an enormous role in the selection and maintenance of robust patterns of activity. However, we believe that this discussion demon- strates how mathematical models and simulations can help us learn to ask the correct types of questions when probing biological systems, and perhaps provide insight into the complex interactions involved in the production of repetitive motor activity.
8.5 Exercises
- Consider a simulation architecture in which L1 and R1 are identical where each has an input of 0 00015 μA , a delay of zero, and a duration of 200 msec. Implement sym- metric mutually excitatory coupling with connection strengths of 1500, and gradually increase the strengths by increments of 100. How might one explain the observed be- havior? What happens when the symmetry is broken; that is, what happens when the coupling strengths are unequal, the injection currents are unequal, or the initial conditions differ?
- It is interesting to note that a fish is capable of swimming backwards when placed in a corner (Grillner 1974). Consider the architecture and parameters for a chain of oscillators described in Sec. 8.3.2. Try experimenting with the parameters to cause the wave to travel in the opposite direction and consider the implications for possible mechanisms for reversing the direction of the wave propagation in a fish.
- Figure 8.3 shows that the pattern of activity for the idealized gallop is identical to that for the walking gait except that the speed is more rapid. Alter the default architecture provided for the walking gait to create a faster gait that is comparable to a gallop. Begin by simulating an increase in the frequency of the intrinsic oscillations via an increase in the current injection.
- In a true gallop, there is a lag between the planting of the left forelimb and the right forelimb. Alter the architecture for the idealized gallop that was developed in the preceding question to obtain this asymmetry. Begin by introducing asymmetries in the connection strengths.
8.5. Exercises 147
- In Sec. 8.3.3 we mention that the differences between the various alternating gaits depend on the timing between the forelimbs and hindlimbs. In fact, the transition from walking to trotting is continuous. As the locomotion speed increases, each forelimb begins to step before the opposite hindlimb touches the ground until the opposite legs step at the same time resulting in a trotting gait (Pearson 1976). Try to mimic the transition from a walk to a trot by making incremental alterations to the architecture and initial conditions.
- In various experiments, both picrotoxin and strychnine, which block inhibitory synapses, induce changes in the motor output pattern of several different organisms (Rand et al. 1988). These changes involve a transformation from an out-of-phase mode of oscillation to an in-phase mode. Delete the inhibitory synapses in the various gait simulation architectures discussed in this chapter and in the exercises above. In each case, do the effects of the blocked inhibition support the choice of architecture or suggest that an alternate choice might be more biologically plausible?
- Use the insight gained from the various models in this chapter and the gait simulation architectures discussed above to create architectures that simulate the non-alternating gaits such as the canter, leap, and half bound.
