JN Fuel your research with LabChart
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
 QUICK SEARCH:   [advanced]


     


J Neurophysiol 93: 1523-1556, 2005. First published October 6, 2004; doi:10.1152/jn.00475.2004
0022-3077/05 $8.00
This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow All Versions of this Article:
93/3/1523    most recent
00475.2004v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Right arrow Citation Map
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via ISI Web of Science (17)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Brezina, V.
Right arrow Articles by Weiss, K. R.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Brezina, V.
Right arrow Articles by Weiss, K. R.

Modeling Neuromuscular Modulation in Aplysia. III. Interaction of Central Motor Commands and Peripheral Modulatory State for Optimal Behavior

Vladimir Brezina1, Charles C. Horn2 and Klaudiusz R. Weiss1

1Department of Physiology and Biophysics and Fishberg Research Center for Neurobiology, Mount Sinai School of Medicine, New York, New York; and 2Monell Chemical Senses Center, Philadelphia, Pennsylvania

Submitted 7 May 2004; accepted in final form 4 October 2004


    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
Recent work in computational neuroethology has emphasized that "the brain has a body": successful adaptive behavior is not simply commanded by the nervous system, but emerges from interactions of nervous system, body, and environment. Here we continue our study of these issues in the accessory radula closer (ARC) neuromuscular system of Aplysia. The ARC muscle participates in the animal's feeding behaviors, a set of cyclical, rhythmic behaviors driven by a central pattern generator (CPG). Patterned firing of the ARC muscle's two motor neurons, B15 and B16, releases not only ACh to elicit the muscle's contractions but also peptide neuromodulators that then shape the contractions through a complex network of actions on the muscle. These actions are dynamically complex: some are fast, but some are slow, so that they are temporally uncoupled from the motor neuron firing pattern in the current cycle. Under these circumstances, how can the nervous system, through just the narrow channel of the firing patterns of the motor neurons, control the contractions, movements, and behavior in the periphery? In two earlier papers, we developed a realistic mathematical model of the B15/B16-ARC neuromuscular system and its modulation. Here we use this model to study the functional performance of the system in a realistic behavioral task. We run the model with two kinds of inputs: a simple set of regular motor neuron firing patterns that allows us to examine the entire space of patterns, and the real firing patterns of B15 and B16 previously recorded in a 21/2-h-long meal of 749 cycles in an intact feeding animal. These real patterns are extremely irregular. Our main conclusions are the following. 1) The modulation in the periphery is necessary for superior functional performance. 2) The components of the modulatory network interact in nonlinear, context- and task-dependent combinations for best performance overall, although not necessarily in any particular cycle. 3) Both the fast and the slow dynamics of the modulatory state make important contributions. 4) The nervous system controls different components of the periphery to different degrees. To some extent the periphery operates semiautonomously. However, the structure of the peripheral modulatory network ensures robust performance under all circumstances, even with the irregular motor neuron firing patterns and even when the parameters of the functional task are randomly varied from cycle to cycle to simulate a variable feeding environment. In the variable environment, regular firing patterns, which are fine-tuned to one particular task, fail to provide robust performance. We propose that the CPG generates the irregular firing patterns, which nevertheless are guaranteed to give robust performance overall through the actions of the peripheral modulatory network, as part of a trial-and-error feeding strategy in a variable, uncertain environment.


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
The ultimate goal of nervous systems is to generate successful adaptive behavior, behavior that enhances the survival and reproduction of the animal. However, no disembodied nervous system can achieve this goal: "the brain has a body" (Chiel and Beer 1997Go). The motor commands of the nervous system must be such as to engage the properties of the animal's body, and those of its environment, so as to change the body, the environment, or their mutual relationship in the required adaptive way. Adaptive behavior is thus not simply "commanded" by the nervous system, but "emerges from interactions of nervous system, body, and environment" (Beer 1995Go; Beer et al. 1999Go; Chiel and Beer 1997Go).

The body, indeed, can have "a mind of its own" (Raibert and Hodgins 1993Go). In vertebrates as well as invertebrates, the structural and dynamical complexity of the periphery can be as large as that of the central nervous system, so that, seen more abstractly, the computational capability of the periphery rivals that of the nervous system that is attempting to control it. The center and periphery are connected through a channel, the neuromuscular transform (NMT) (Brezina et al. 2000bGo), which can be very narrow relative to the rich computational capabilities at either end (Brezina and Weiss 2000Go; Brezina et al. 2000cGo). The semiautonomous periphery may then facilitate and simplify, but can also considerably complicate and constrain, motor control (Beer 1990Go; Brezina and Weiss 2000Go; Chiel and Beer 1997Go; Dickinson et al. 2000Go). A fundamental question is how the two halves of the system are matched in their properties, and coordinate in their dynamics, so as to ensure the production of functional behavior.

Here we continue our long-standing study of this question in an advantageous model system, the accessory radula closer (ARC) neuromuscular system of Aplysia. The ARC muscle is controlled by just two motor neurons. In contrast to this central simplicity, there is great peripheral complexity, generated by a dynamically complex network of actions of numerous endogenous modulatory neurotransmitters and neuropeptides that shape the contractions of the muscle. For any central motor command, this modulatory network multiplies the range of possible peripheral outputs of the NMT. Extensive studies over the past 25 years (reviewed by Kupfermann et al. 1997Go; Weiss et al. 1993Go; see also Brezina et al. 2003aGo) have characterized many aspects of the ARC system, in particular the various modulatory actions in quantitative terms, to the point where integrative understanding has become possible. We now briefly summarize the details essential for the work in this paper.

The ARC muscle is one of the muscles of the buccal mass, a complex structure that produces rhythmic, cyclical movements of the animal's food-grasping organ, the radula, in several types of consummatory feeding behavior such as biting, swallowing, and rejection of unsuitable food (Kupfermann 1974Go). All of these behaviors are driven by a single, multitasking central pattern generator (CPG) (for recent review see Elliott and Susswein 2002Go) whose patterned activity is then distributed to the motor neurons of the various buccal muscles. The ARC muscle's two motor neurons, B15 and B16, release their classical transmitter, ACh, to depolarize and contract the muscle (Cohen et al. 1978Go; Kozak et al. 1996Go). In addition, however, both motor neurons release modulatory peptide cotransmitters that shape the basal ACh-induced contractions. B15 releases, in particular, the small cardioactive peptides (SCPs) (Cropper et al. 1987aGo, 1990bGo; Lloyd et al. 1984Go; Vilim et al. 1996bGo; Whim and Lloyd 1989Go), and B16 releases the myomodulins (MMs) (Brezina et al. 1995Go; Cropper et al. 1987bGo, 1991Go; Vilim et al. 2000Go). The SCPs and MMs shape contractions through three main actions on the muscle: 1) they potentiate the contractions by enhancing the muscle's depolarization-activated Ca current that supplies Ca2+ essential for contraction (Brezina and Weiss 1995Go; Brezina et al. 1994aGo, 1995Go; Cropper et al. 1987a, bGo, 1991Go; Lloyd et al. 1984Go; Whim and Lloyd 1990Go); 2) they depress the contractions by activating in the muscle a K current that opposes the ACh-induced depolarization, activation of the Ca current, and Ca2+ influx (Brezina et al. 1994bGo, 1995Go; Cropper et al. 1987bGo, 1991Go; Orekhova et al. 2003Go); and 3) they increase the relaxation rate of the contractions probably by modulating the muscle's contractile machinery (Brezina et al. 1995Go; Cropper et al. 1990bGo; Heierhorst et al. 1995Go; Probst et al. 1994Go; Whim and Lloyd 1990Go). The balance of the competing potentiating and depressing actions determines the net modulation of contraction size, which, together with the modulation of the relaxation rate, then produces the final shape of the contraction. [It should be mentioned that all this happens against the background of yet other, although probably less critical, actions of these and other modulators that are present in or act on the ARC system (see Brezina et al. 2003bGo).] The main elements of this modulatory network are graphically summarized within the gray box in Fig. 1.



View larger version (33K):
[in this window]
[in a new window]
 
FIG. 1. Schematic summary of our model of the modulated B15/B16-ARC (accessory radula closer) neuromuscular transform and our use of it in this paper. For details see INTRODUCTION, METHODS, and The model, functional task, and performance measure in RESULTS.

 
Our understanding of how the ARC system operates has over the years been embodied in a series of models, at first purely conceptual (e.g., Weiss et al. 1993Go) but increasingly, with more quantitative measurements made in the real system, more realistic models (Brezina and Weiss 2000Go; Brezina et al. 1996Go, 2000b, cGo). Most recently, in Parts I and II of the present study, we have constructed a fully quantitative dynamic model of the network of modulation in the ARC system and of its dependency on the firing of the two motor neurons B15 and B16—essentially, a model of the dynamic interdependence of all of the variables X(t) down the right-hand side of Fig. 1 (Brezina et al. 2003aGo). Joining this modeled "neuro-modulation transform" to a model of the basal, ACh-mediated neuromuscular transform—down the left-hand side of Fig. 1—gave a complete model of the modulated B15/B16-ARC NMT (Brezina et al. 2003bGo). For any input pattern of firing of the two motor neurons (top of Fig. 1), the model outputs the resulting modulated contraction shape (bottom right of Fig. 1). [This input–output relationship is, in essence, the definition of the NMT (Brezina et al. 2000b, cGo).] In particular, the model makes it possible to follow the modulated contraction waveform—and, of course, all of the internal modeled variables—as the motor neurons B15 and B16 fire in different patterns, as they do in the intact feeding animal, over hundreds of cycles that constitute a typical meal (Kupfermann 1974Go; Susswein et al. 1976Go, 1978Go; Weiss et al. 1982Go).

For adaptive behavior, however, what is important is not the contraction waveform per se, but rather how well it satisfies the functional requirements of the desired behavior (Brezina and Weiss 2000Go). As a final step, therefore, we will in our work here pass the contraction waveform through a realistically modeled behavioral task (schematized in Fig. 1, bottom right, and explained further in RESULTS) to evaluate its functional performance.

Running the model will immediately answer basic (although, given the complexity of the system, by no means trivial) questions such as: Given a firing pattern, what contractions and functional performance will it produce? Conversely, given a functional requirement, what firing patterns are optimal? Beyond this, however, the model now allows us to answer some general questions of motor control in this representative brain–body system. Some of these questions were raised in our previous work (see Brezina and Weiss 2000Go; Brezina et al. 2000b, cGo) but could then be addressed only conceptually or under unrealistically restricted conditions. Now, for the first time, we have the tools to answer them in a reasonably realistic fashion. Furthermore, for the first time, we can feed as input into the model an extended dataset of the real firing patterns of the motor neurons B15 and B16, recorded with chronic electrodes implanted in freely feeding animals by Horn et al. (2004)Go. With these tools, we will in this paper focus, in particular, on the following questions.

1) Why is the additional layer of peripheral modulation necessary? Why cannot the motor neurons fire in such patterns as to immediately give the required contraction shapes through the basal NMT alone? The conceptual answer given before (Brezina et al. 2000cGo) is that, with fixed properties, the channel of the NMT is too narrow to allow even the entire set of firing patterns to create the broad range of contraction shapes that is required by the range of feeding behaviors in which the neuromuscular plant participates. Is this still true in the realistic model, where, for example, the modulators are no longer exogenously superimposed but rather are released endogenously by the activity of the system? Note that an affirmative answer requires demonstrating not merely that modulation improves the performance of a specific firing pattern, in particular the real pattern recorded in the feeding animal—although that is certainly a prerequisite—but furthermore that no other pattern exists that, without modulation, could have given a better performance.

2) The previous work (e.g., Brezina and Weiss 1997Go; Brezina et al. 1996Go, 2000bGo, 2003a, bGo) has suggested that the various components of the system—the two motor neurons, their modulators, and the different modulator effects—act in interdependence and mutual complementarity, one or another playing a key role depending on the circumstances. Furthermore, there is strong nonlinear coupling between the components, such that the functional role of each component depends on the context of the others. To what extent is this important for performance in the realistic model, with real motor neuron firing patterns?

3) The dynamics of the three main modulator effects have very different timescales. The modulation of K current is fast enough to react to each cycle of behavior, whereas the effects on Ca current and relaxation rate are slowly integrated over many cycles (Brezina et al. 2003a, bGo). Accordingly, it was proposed that whereas the former provides local adjustments, the persistence of the latter effects provides a peripheral memory of past behavior, which then prepares the system for future behavior of the same kind. However, these proposed roles—and other previous ideas—must be reexamined in light of the finding of Horn et al. (2004)Go that the real firing patterns of motor neurons B15 and B16, along with practically all other parameters of Aplysia feeding behavior, are not highly stereotyped as was assumed, but show very large, quasi-random variability from cycle to cycle.

4) A fundamental question then is: With stochastic central motor commands and slow, history-dependent dynamics in the periphery, to what extent is the system actually under control of the nervous system? When the nervous system sends motor neuron firing patterns through the NMT, to what extent can it predict the contractions and performance that will result? To what extent can it, through those patterns, control the various desired behaviors?

To answer these questions, we proceed in this paper in two complementary ways. First, as in the previous work, we run the model with a set of simple, regular motor neuron firing patterns. Although idealized, these patterns allow us to examine the complete space of firing pattern parameters as we must do if we are to demonstrate a global optimization of performance. Then, we use a dataset of the real, irregular firing patterns recorded by Horn et al. (2004)Go. These patterns, conversely, are fully physiological but their parameters cannot easily be manipulated. Although each approach has its deficiencies, in combination they provide, as will be seen, considerable insight into motor control in the Aplysia B15/B16-ARC neuromuscular system and, we believe, other complex brain–body systems.


    METHODS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
Model

The equations and parameter values of the model were taken from Brezina et al. (2003a, bGo; also 1996Go, 2000aGo), and used without modification, except where explicitly stated. The major components of the model are shown in Fig. 1. Below, after an initial overview, each section gives first the model and its implementation as used with the ideal set of regular firing patterns in Figs. 2 9, then any variations used with the real dataset of irregular firing patterns in Figs. 1019. The meanings of certain symbols—principally those that are often used in the RESULTS—are repeated in Table 1.



View larger version (25K):
[in this window]
[in a new window]
 
FIG. 2. Release of endogenous small cardioactive peptide (SCP) by motor neuron B15 improves performance: steady-state simulations with regular firing of B15 only. A: summary maps of steady-state performance over the space of parameters of the regular B15 firing pattern. Each of the 10 maps (identified by the white numbers 1–10 in the top right corner) contains a 10 x 10 array of pixels, each pixel representing a separate simulation with the particular values of the burst duration dintra,B15, intraburst firing frequency fintra,B15, and cycle period P given by its location within the map. Across the maps from left to right, P ranges from 1 to 100 s. On the x-axis of each map, dintra,B15 ranges from 0 to P, except with P = 100 s, where dintra,B15 ranges only to P/10. On the y-axis of each map, fintra,B15 ranges from 0 to 30 Hz. Here and in subsequent figures, where an axis of a map is unlabeled, it is the same as the axis labeled at the beginning of that row or bottom of that column of maps. On the z-axis of each map, the color represents the magnitude of the steady-state performance measure, M{infty}, evaluated in the last 1–2 cycles of the contraction waveform that was computed long enough to be close to the steady state (see METHODS). M{infty} ranges from 0 (blue) to 0.121, the highest performance in any of the simulations in this figure (red). Each map thus shows M{infty} over the space dintra,B15 x fintra,B15. Top row of maps shows performance without modulation (see METHODS), the bottom row the corresponding performance when the B15 firing also released endogenous SCP. B: representative contraction waveforms, from the locations marked a–f in the maps in A. Shown are the last 2 cycles close to the steady state, all with P = 3 s and fintra,B15 = 21 Hz, and with dintra,B15 = 1.5, 2.1, and 2.7 s, unmodulated (waveforms a–c) and modulated (waveforms d–f). Four vertical lines superimposed over waveforms b and d mark the timepoints t1–t4 at which M{infty} > 0 was computed relative to the contraction threshold = 0.5 (see Fig. 1 and METHODS); waveforms a, c, e, and f had M{infty} = 0.

 


View larger version (38K):
[in this window]
[in a new window]
 
FIG. 9. Different dynamical time scales uncouple components of the modulated B15/B16-ARC neuromuscular transform: transient simulations with regular firing of both B15 and B16. A: representative example. Motor neurons B15 and B16 were fired with the pattern dintra,B15 = 3.5 s, fintra,B15 = 24 Hz, dintra,B16 = 2.0 s, fintra,B16 = 21 Hz, P = 5 s, repeating constantly except for one short pause, as shown by the records of fB15 and fB16 at the bottom. Records above show (bottom to top) the corresponding evolution of the modulatory effects: the enhancement of Ca current C, activation of K current K, modulation of contraction size S, and increase of relaxation rate R; the contraction amplitude c; and the performance M in each cycle. B: summary maps of performance. At the 3 timepoints in A marked by the vertical lines—at the very beginning, or at 2 times when the simulation was halted at the current values of all model variables—the performance M in the next cycle was evaluated over the space of B15 and B16 parameters dintra,B15 x fintra,B15 x dintra,B16 x fintra,B16, with P = 5 s, as in Figs. 68. Again, the 2 maps in each box are simultaneous 2-dimensional sections dintra,B15 x fintra,B15 and dintra,B16 x fintra,B16 through the complete 4-dimensional space at the point of maximal M. This point of maximal M was: left, (dintra,B15 = 0.17 s, fintra,B15 = 30 Hz, dintra,B16 = 0.17 s, fintra,B16 = 20.0 Hz), M = 0.0868; middle, no specific point because M was zero everywhere; right, (dintra,B15 = 0.17 s, fintra,B15 = 30 Hz, dintra,B16 = 2.58 s, fintra,B16 = 25.9 Hz), M = 0.177.

 


View larger version (72K):
[in this window]
[in a new window]
 
FIG. 10. Representative dataset of irregular firing of motor neurons B15 and B16 in a freely feeding animal. A: in vivo chronic recording of electrical activity in the ARC muscle during an entire, approximately 21/2-h-long meal (record 1). Records 2–6 successively expand the indicated portions of the recording. Records 3 and 6 show presumed rejection cycles at the end of the meal (see Real, irregular motor neuron firing patterns in RESULTS). The assumption that only B16 fired during these cycles provided an amplitude threshold (horizontal dashed line in records 5 and 6), which was then used throughout the recording to allocate spikes smaller than the threshold to B16 and spikes larger than the threshold (asterisks in record 5) to B15. Gray rectangle in record 5 marks a burst of firing identified by our criteria (see METHODS). B: instantaneous firing frequencies of B15 and B16 over the entire meal (i.e., analyzed from the entire record 1 in A).

 


View larger version (51K):
[in this window]
[in a new window]
 
FIG. 19. Performance by the irregular and regular firing patterns in regular and irregular tasks: details of individual cycles. A: irregular task, regular firing pattern and contraction. Probability density distribution (top) and the performance M (bottom) as a function of the ratio of the mean (period-averaged) contraction amplitude to the threshold amplitude (plotted on a log scale) over the 749 cycles of the "Optimal, modulated" meal in Fig. 18B. B: regular task, irregular firing pattern and contraction. As in A, except over the 749 cycles of the standard 21/2-h irregular meal with fixed = 0.132, as for example in the "Irregular" case in Fig. 18A. C: irregular task as well as irregular firing pattern and contraction. As in A and B, except over the 749 cycles of the "Irregular, modulated" meal in Fig. 18, B and C. D: one cycle (in the steady state) of the "optimal" regular contraction waveform plotted in A. Values of , <c>, and M are given for one representative location marked by the blue circle in A. E: 9 representative cycles of the irregular contraction waveform in the standard 21/2-h irregular meal. Values of and M are given for the irregular task in C; <c> is the same for both B and C. Locations of waveforms 1, 2, 4, 6, 8, and 9 are marked by the red, orange, pink, green, purple, and yellow circles, respectively, in B and C. In D and E, the vertical gray rectangle indicates the burst duration dintra, to which each waveform has been scaled horizontally; the actual cycle period P in that cycle is given under the right-hand end of the waveform. See further With an irregular task, irregular firing patterns are optimal in RESULTS.

 

View this table:
[in this window]
[in a new window]
 
TABLE 1. Frequently used symbols

 
OVERVIEW. The input to the model is the firing pattern—the spike patterns or, more precisely, the instantaneous firing frequency waveforms fB15(t) and fB16(t)—of the motor neurons B15 and B16 (Fig. 1, top). This input sets in motion two cascades of action that are modeled as initially independent of each other but ultimately joining.

First, the motor neuron firing pattern drives the "neuro-modulation" transform, down the right-hand side of Fig. 1. The firing patterns f(t) elicit release of the modulators SCP and MM, r(t). The released modulator amounts are converted to effective concentrations, C(t). The concentrations elicit the three main modulator effects: the enhancement of Ca current, C(t), the activation of K current, K(t), and the increase of the muscle's relaxation rate, R(t). The Ca-current and K-current effects then jointly determine the modulation of contraction size, S(t). All of these variables and their relationships are explicitly modeled as a function of the time t.

At the same time, the motor neuron firing pattern drives the basal NMT, down the left- hand side of Fig. 1. This is modeled without any internal details, simply as the direct conversion of the firing patterns f(t) to a basal muscle contraction waveform, cbasal(t) (Fig. 1, bottom left).

Finally, the outputs of the "neuro-modulation" transform, that is, S(t) and R(t), are applied to convert cbasal(t) to the modulated contraction waveform, c(t) (Fig. 1, bottom right). Altogether, for any input motor neuron firing pattern fB15(t) and fB16(t), the model outputs the corresponding modulated contraction c(t). It is thus a model of the complete modulated B15/B16-ARC NMT.

FIRING PATTERNS OF MOTOR NEURONS B15 AND B16. As in Brezina et al. (2003a, bGo), we use here simple cyclical bursting patterns definable by the burst duration dintra, the interburst interval dinter, and the intraburst firing frequency fintra, with cycle period P = dintra + dinter. We also define the duty cycle D {equiv} dintra/P. The instantaneous firing frequency f, as a function of time t, is then given by

(1a)
where t*, the start time of the current cycle, is the largest number ≤t generated by the sequence

(1b)
where N is the cycle number, and L is the length of the entire sequence of firing which began at t = 0.

With the regular firing patterns, the parameter values P, dintra,B15 (and so, by dinter = Pdintra, dinter,B15), fintra,B15, dintra,B16 (and so dinter,B16), and fintra,B16 were set and varied as described in RESULTS and figure legends. Specification of these parameters gave two particular functions of the form of Eq. 1Go, fB15(t) and fB16(t). With a common cycle period P and a common sequence of cycle start times t*, the bursts of B15 and B16 firing always ended together at the end of each cycle, although, with different durations dintra,B15 and dintra,B16, they could start at different times [see Fig. 3B of Brezina et al. (2003a)Go]. L was set as described below.



View larger version (57K):
[in this window]
[in a new window]
 
FIG. 3. Analysis of the improvement of performance by SCP. A: summary maps of steady-state performance M{infty} over the space of regular B15 parameters dintra,B15 x fintra,B15, as in Fig. 2A but showing a finer, 30 x 30 array of simulations, with no endogenous SCP release but instead with exogenous application of 0, 10 nM, 100 nM, 1 µM, and 10 µM SCP (left to right), all with P = 5 s. Map 1, with no modulation at all, is thus a finer version of map 5 in Fig. 2A. Scale of M{infty} (shown below B, and applying to maps 1–6) is the same as that in Fig. 2A for consistency. B: as in A, without exogenous SCP application but with endogenous SCP release. Map 6 shows M{infty}, map 7 the steady-state mean concentration of SCP, <CSCP>{infty}, attained. Log scale of <CSCP>{infty} saturates below 1 nM (black) and above 10 µM (white). Red and black lines connect 2 representative triplets of values (see SCP release from motor neuron B15 improves performance in RESULTS). C: plot of the steady-state mean contraction amplitude, <c>{infty}, against fintra,B15 with fixed dintra,B15 = 2.5 s, unmodulated and with endogenous SCP released, i.e., along the white vertical lines in maps 1 and 6. Points with M{infty} > 0 are marked by the squares color-coded on the same scale as the maps.

 
The dataset of irregular firing patterns was derived from a single in vivo chronic recording of electrical activity in the ARC muscle during a spontaneous, approximately 21/2—h—long meal (see RESULTS and Horn et al. 2004Go), shown in Fig. 10A. Allocation of each spike in this extracellular record to either B15 or B16 (see RESULTS), followed by conversion of the spike times to firing frequencies by assigning to each timepoint in each interspike interval the reciprocal of the duration of that interspike interval, gave the instantaneous firing frequency functions fB15(t) and fB16(t) shown in Fig. 10B. These replaced Eq. 1Go as the primary input to the model in Figs. 1219. For the purposes of analysis, we further partitioned the record of all spike times (before their allocation to B15 or B16) into individual cycles. As in Eq. 1Go, we assumed that there was one burst of spikes per cycle. Based on an examination of the distribution of interspike intervals in the record, we defined a burst as a continuous sequence of ≥10 spikes with interspike intervals of <1 s. Automatic implementation of this criterion identified 749 bursts in the 21/2-h-long record. Note that application of this procedure before the allocation of B15 and B16 spikes made the B15 and B16 bursts in each cycle nominally identical in their timing and duration, in contrast to Eq. 1Go for the regular firing patterns. Between bursts, the interburst firing frequency for B15 and B16, jointly as well as individually, was by definition <1 Hz, and in most simulations was in fact explicitly set to zero—thus altering fB15(t) and fB16(t) somewhat—to allow implementation of the modulation of the muscle's relaxation rate (Eqs. 6c and e below).



View larger version (28K):
[in this window]
[in a new window]
 
FIG. 12. Trajectories of key components of the modulated B15/B16-ARC neuromuscular transform during the 21/2-h irregular meal. At the top are the inputs to the transform, the instantaneous firing frequencies of B15 and B16, fB15 and fB16, from Fig. 10BGo. Below are the instantaneous concentrations of SCP and MM, CSCP and CMM, and the instantaneous effects on Ca current, K current, contraction size, and relaxation rate, C, K, S, and R, respectively. At the bottom is the output of the transform, the instantaneous contraction amplitude c.

 
Equation 1a implies that the interburst interval precedes the burst in each cycle, and we generally applied this rule—with the irregular as with the regular patterns—in computing, for example, cycle periods or duty cycles (Figs. 11, 14, 16). However, because the muscle contraction caused by a spike burst outlasted the burst considerably into the following interburst interval, it was more appropriate to take the burst together with the following interval as the cycle over which to evaluate quantities such as the contraction area, the functional performance of the contraction waveform (see below), and in some simulations (Fig. 17A) the corresponding period-averaged modulatory effects. This distinction between the preceding and the following interburst interval only became an issue with the irregular firing patterns in Figs. 1318; with the regular firing in the steady state throughout Figs. 28, the distinction was irrelevant.



View larger version (58K):
[in this window]
[in a new window]
 
FIG. 11. Statistics of the irregular firing of B15 and B16. Two records at top left repeat the records of the instantaneous firing frequencies of B15 and B16, fB15 and fB16, from Fig. 10B. An expanded view of the indicated portion of the records, illustrating their constituent discrete bursts of firing, is shown at top right. Aligned below the records of instantaneous firing frequencies are time series of 4 parameters analyzed from them in each cycle (n = 749): top to bottom, the mean firing frequency of B15 in each burst, corresponding approximately to fintra,B15 of the regular firing patterns; the mean firing frequency of B16 in each burst, corresponding approximately to fintra,B16; the burst duration, corresponding to dintra; and the cycle period P. Points along the top of each time series are offscale. To the immediate right of each time series is a marginal histogram summarizing the distribution of the values, scaled so as to approximate the probability density function. The 25th, 50th (median), and 75th percentiles of the distribution are marked by the horizontal lines. Last (top) bin of each histogram contains offscale values. Histogram of P has been magnified 2.5 times for clarity. To the far right of each time series is a histogram, scaled so as to approximate the probability density function, summarizing the distribution of normalized cycle-to-cycle differences (see Real, irregular motor neuron firing patterns in RESULTS). Vertical lines mark the 25th, 50th, and 75th percentiles, and the standard deviation of the distribution, {sigma}, is given. First and last bins of each histogram contain offscale values.

 


View larger version (49K):
[in this window]
[in a new window]
 
FIG. 14. Modulation improves performance during the 21/2-h irregular meal. A: time series of the performance M computed from the modulated contraction record in Fig. 13 with threshold amplitude = 0.132. Horizontal line shows the mean M, as then plotted in Fig. 15A. B: performance values in A, and corresponding values computed from the unmodulated record in Fig. 13, plotted as a function of the cycle period P. Each point is the mean ± SE over a range of P, with n = 24–212. Statistical significance was tested with 2-way ANOVA followed by pairwise multiple comparisons using the Holm–Sidak test: for the modulated vs. unmodulated group at each P < 30 s, P < 0.001 (***). C: summary maps of the unmodulated (left) and modulated (right) performance. In each case the entire set of values of M has been projected onto the plane spanned, on the x-axis, by the duty cycle D {equiv} dintra/P and, on the y-axis, by either fintra,B15 (top) or fintra,B16 (bottom) as defined in Fig. 11. Thus, the top and bottom maps in each column are different projections of the same set of values. Points along the top of each map are offscale. Color scale of M extends to M = 0.3; several higher values, up to M = 0.447, are shown in white. D: maps as in C but showing only the subset of values with M > 0.15.

 


View larger version (42K):
[in this window]
[in a new window]
 
FIG. 16. Analysis of differences in performance reveals circumstances in which components of the modulation make particular contributions. A–C: three examples. Performance M in each of the 749 cycles of the 21/2-h irregular meal simulated with one set of modulatory components was subtracted from the corresponding performance with a different set of modulatory components. Resulting differences in performance, {Delta}M, are then plotted here as a function of the duty cycle D. In all cases = 0.132. A: endogenous modulators released – no modulation. Points with {Delta}M > 0 are cycles in which performance was better modulated than unmodulated, points with {Delta}M < 0 the converse. Points along the top of the plot are offscale. B: endogenous modulators released – R only. C: C only – K only. Two red circles in each of A–C show the mean ± SE (in most cases invisible because smaller than the circle size), of D weighted by {Delta}M in the horizontal dimension and of {Delta}M in the vertical dimension, separately of the points with {Delta}M > 0 and {Delta}M < 0 (points with {Delta}M = 0 were ignored). Statistical significance of the differences between the means in the horizontal dimension was tested with the t-test: ** indicates P < 0.01. D: representative cycles, corresponding to the points a and b in B and c and d in C. Light- and mid-gray rectangles mark respectively the interburst interval and the burst of the representative cycle. Four vertical lines superimposed over the contraction waveforms indicate schematically the timepoints t1–t4 at which M > 0 was computed (in each case from one of the pair of waveforms only: the other had M = 0) relative to the contraction threshold = 0.132 (see METHODS).

 


View larger version (40K):
[in this window]
[in a new window]
 
FIG. 17. Consequences of the dynamics of the modulatory state for performance and control of the system. A: means ± SE of the performance M during the 21/2-h irregular meal (bottom) and the percentage of cycles with M > 0 (top) as in Fig. 15, with = 0.132. Simulations were run under the following circumstances (left to right): "No modulation": repeating the "No modulation" values from Fig. 15A. "Actual modulation": as "Endogenous modulators" in Fig. 15A, except fixing, for consistency with the following simulations in this set, the values of C, K (and so S), and R throughout each cycle to the mean, period-averaged values <C>, <K> (<S>), and <R> computed over that cycle. "Scrambled modulation": having randomly mixed the set of modulatory states (<C>, <K>, <R>) with respect to the set of cycles. "Scrambled K only": having mixed only <K>. "Scrambled C + R": having mixed (<C>, <R>). "Slow K": with the mean <C>, <K>, and <R>, after dividing the rate constants kK+ and kK of K in Eq. 4b (METHODS) by 375. "Fast C": with the mean <C>, <K>, and <R>, after multiplying the rate constants kC+ and kC of C by 250. "Fast R": with the mean <C>, <K>, and <R>, after multiplying the rate constants kR+ and kR of R by 500. Statistical significance was tested with one-way ANOVA followed by pairwise multiple comparisons using the Holm–Sidak test. For the differences from "Actual modulation" in each case, *** indicates P < 0.001; **, P < 0.01; and *, P < 0.05. B: trajectory of the modulatory state during the 21/2-h irregular meal in the (S, R) plane. Records of S and R from Fig. 12 are here plotted against each other. Gray region is the region traced out by the steady-state effects (S{infty}, R{infty}) of all possible combinations of the SCP and MM concentrations CSCP and CMM [see, e.g., Eq. 8 in METHODS, and Eqs. 5d–i of Brezina et al. (2003a)Go]. C: correlation strengths between selected components of the modulated B15/B16-ARC neuromuscular transform during the 21/2-h irregular meal. Inset at right illustrates a typical correlation plot, from the composite parameter fintra,B16D to K, showing all 749 cycles of the meal (points) and the best linear fit (line) with a coefficient of determination R2 = 0.50. Main diagram summarizes 25 such plots; the thickness of each arrow is proportional to R2 and the value of R2 is given. See further Dynamics of the modulatory state and control of the system in RESULTS.

 


View larger version (56K):
[in this window]
[in a new window]
 
FIG. 13. Statistics of the irregular contractions produced by the neuromuscular transform. Counterpart to Fig. 11. The 21/2-h meal of B15 and B16 firing frequencies was run through the transform without modulation (blue) and with the endogenous modulators released (red). At top left are the resulting records of instantaneous contraction amplitude c (the red, modulated record is repeated from Fig. 12, bottom). Black horizontal lines mark the threshold amplitudes = 0.132 and = 0.5 that were used to compute performance in Figs. 1418. Below the contraction records are time series showing the peak contraction amplitude and the area under the contraction waveform in each cycle. To the right of the time series are histograms summarizing the distributions of the individual-cycle values (the contraction-area histogram has been magnified 3 times for clarity) and the cycle-to-cycle differences, as in Fig. 11.

 


View larger version (43K):
[in this window]
[in a new window]
 
FIG. 18. Performance by the irregular and comparable regular firing patterns in regular and irregular tasks: overall comparison. A: regular task, with fixed = 0.132. Means ± SE of the performance M (bottom) and the percentage of cycles with M > 0 (top), as in Figs. 15 and 17A. Simulations were as follows (left to right): "Irregular": the 21/2-h irregular meal with full modulation, repeating "Endogenous modulators" from Fig. 15A. "Median of irregular, regular": a regular pattern constructed with the median parameter values of the 749 cycles of the irregular meal from Fig. 11, that is, fintra,B15 = 2.27 Hz, fintraB16 = 8.46 Hz, dintra = 3.28 s, and P = 7.36 s. "Best of irregular, regular": a regular pattern constructed with the (mean, period-averaged) parameter values of the cycle that has the best performance, M = 0.447, in the irregular meal, that is, fintra,B15 = 0.53 Hz, fintraB16 = 10.29 Hz, dintra = 1.84 s, and P = 4.55 s. "Optimal, regular": the best-performing regular pattern that was found in a search through the whole parameter space dintra,B15 x fintra,B15 x dintra,B16 x fintra,B16 at P = 5 s, as in Figs. 69. This pattern was (dintra,B15 = 2.18 s, fintra,B15 = 7.7 Hz, dintra,B16 = 2.06 s, fintra,B16 = 8.0 Hz). Except for "Irregular", the steady-state performance is shown. All simulations were run with the full modulation. B: similar to A, but with an irregular, variable task. In each cycle, was randomly selected from the positive part of a Gaussian distribution with mean and SD both equal to 0.132. Simulations were as follows (left to right): "Irregular, unmodulated": the 21/2-h irregular meal without modulation. "Irregular, modulated": the 21/2-h irregular meal with the full modulation. "Median of irregular, modulated": 749 cycles of the regular median pattern from A, with full modulation. "Optimal, modulated": 749 cycles, with full modulation, of the regular pattern that gave the best mean performance, over 749 cycles of the same irregular task, in a search through the whole parameter space dintra,B15 x fintra,B15 x dintra,B16 x fintra,B16 at P = 5 s. This pattern was (dintra,B15 = 2.24 s, fintra,B15 = 8.0 Hz, dintra,B16 = 2.06 s, fintra,B16 = 8.0 Hz). Statistical significance was tested, separately in A and B, with one-way ANOVA followed by pairwise multiple comparisons using the Holm–Sidak test. For all differences from the "Irregular" case in A and the "Irregular, modulated" case in B, P < 0.001 (***). C: summary maps of the performance M in the "Irregular, modulated" meal in B, projected onto the plane spanned by the duty cycle D and fintra,B15 (top) or fintraB16 (bottom) as in Fig. 14, C and D.

 


View larger version (28K):
[in this window]
[in a new window]
 
FIG. 8. Simultaneous release of SCP and MM optimizes performance over a range of cycle periods: steady-state simulations with regular firing of both B15 and B16. Summary maps of steady-state performance M{infty} over the space of regular B15 parameters dintra,B15 x fintra,B15 (top row) or regular B16 parameters dintra,B16 x fintra,B16 (bottom row), as in maps 7 and 8 of Fig. 6 (repeated here as maps 5 and 6), with P ranging from 1 to 100 s (left to right). Again, the maps of each pair enclosed by the box are simultaneous 2-dimensional sections through the 4-dimensional dintra,B15 x fintra,B15 x dintra,B16 x fintra,B16 space at each P, with the 2 parameters not varied in the map fixed at values that, in a search through the space, gave the maximal M{infty}. This point of maximal M{infty} was: with P = 1 s, (dintra,B15 = 0.41 s, fintra,B15 = 30 Hz, dintra,B16 = 0.31 s, fintra,B16 = 20.2 Hz), M{infty} = 0.245; with P = 3 s, (dintra,B15 = 1.35 s, fintra,B15 = 30 Hz, dintra,B16 = 1.35 s, fintra,B16 = 15.0 Hz), M{infty} = 0.156; with P = 5 s, (dintra,B15 = 2.5 s, fintra,B15 = 30 Hz, dintra,B16 = 2.5 s, fintra,B16 = 15.0 Hz), M{infty} = 0.151; with P = 10 s, (dintra,B15 = 4.3 s, fintra,B15 = 30 Hz, dintra,B16 = 5.0 s, fintra,B16 = 14.2 Hz), M{infty} = 0.0987; and with P = 100 s, (dintra,B15 = 1.4 s, fintra,B15 = 30 Hz, dintra,B16 = 1.6 s, fintra,B16 = 9.8 Hz), M{infty} = 0.0461. Color scale of M{infty} extends to M{infty} = 0.151 for consistency with Figs. 6 and 7; several higher values, up to M{infty} = 0.245, are shown in white.

 
MODULATOR RELEASE. Two modulators are considered: SCP released from B15 and MM released from B16. With either the regular or the irregular f(t) as input, the release r is described by

(2a)

(2b)

(2c)
where S is the size of the releasable pool of modulator, p is the probability of release of the modulator in S, kp+ and kp– are the "ON" and "OFF" rate constants of p, y is an exponent providing nonlinearity of release, and Q10 provides for dependency on the temperature T. Here T = 15°C and so Q10(T–15)/10 = 1. The other parameter values are (Brezina et al. 2000aGo, 2003aGo) kp+,B15 = 4.0 x 10–10, kp+,B16 = 7.4 x 10–11, kp–,B15 = 3.4 x 10–3 s–1, kp–,B16 = 6.6 x 10–3 s–1, and y = 3. Equation 2c describes depletion of S from an initial pool size S0,SCP {equiv} SSCP(t = 0) = 541 fmol and S0,MM = 2690 fmol. For simplicity, in all simulations here we prevented depletion by setting dS(t)/dt = 0 in Eq. 2c and so S(t) = S0.

MODULATOR CONCENTRATIONS. With r(t) given by Eq. 2GoGo as input, the effective concentrations, C, of SCP and MM are described by

(3)
where {nu} is the effective volume of the released modulator in the muscle and kC is the rate constant of removal of the modulator from {nu}. The parameter values are (Brezina et al. 2003aGo) {nu}SCP = 10 µl, {nu}MM = 35 µl, kC,SCP = 0.1 s–1, and kC,MM = 0.35 s–1.

MODULATOR EFFECTS. Three independent effects are considered: enhancement of Ca current, C; activation of K current, K; and increase of the muscle's relaxation rate, R. With C(t) given by Eq. 3 as input, each effect X {C, K, R} is described by

(4a)

(4b)
where C'i {equiv} CihX,i is the concentration of the modulator Ci modified by a Hill coefficient hX,i, Xmax,i is the maximal effect of the modulator, kX+,i and kX–,i are "ON" and "OFF" rate constants, and Fji specifies the fraction of the effect of modulator j that occludes the effect of modulator i. The parameter values are (Brezina et al. 2003aGo) Cmax = 80%, kC+ = 750 s–1 MhC, kC = 2.25 x 10–3 s–1, hC = 0.75, FC = 1; Kmax,SCP = 27%, Kmax,MM = 100%, kK+ = 1.9 x 105 s–1 MhK, kK–,SCP = 3.5 s–1, kK–,MM = 0.35 s–1, h K = 0.85, FK,SCP,MM = 1, FK,MM,SCP = 0.27; Rmax = 100%, kR+ = 305 s–1 MhR, kR = 1.11 x 10–3 s –1, hR = 0.7, FR = 1.

In one set of simulations (Fig. 17A) we altered the rate constants of the effects, either dividing kK+ and kK by 375 to make K approximately as slow as C and R, or multiplying kC+ and kC by 250, or kR+ and kR by 500, to make C or R as fast as K.

MODULATION OF CONTRACTION SIZE. C(t) and K(t) given by Eq. 4Go combine to determine the modulation of contraction size S, according to the equation

(5)
with the parameter values Smax = 400%, KS = 8%, {kappa} = 9%, and a = 8/3% (see Brezina et al. 1996Go, 2003aGo).

THE MODULATED CONTRACTION WAVEFORM. The instantaneous modulated contraction amplitude c is modeled as a basal contraction amplitude cbasal, then modulated by S(t) and R(t), according to the equations

(6a)

(6b)

(6c)

(6d)

(6e)
where cbasal,{infty}(f), rcontr(f), and rrelax(c) are empirical functions plotted in Fig. 1, B–D of Brezina et al. (2003b)Go, t* is the start time of the current cycle given by Eq. 1b with the regular f(t) or the end of the previous burst with the irregular f(t), {delta}R = 0.03, and

(6f)

When fNMT(t) > 32 Hz or c(t) > 1, the functions cbasal,{infty}(f); rcontr(f); and rrelax(c) were extrapolated. For further details see Brezina et al. (2003b)Go.

We routinely ran parallel simulations with modulation, as above, and without modulation, setting S(t) = 100 and R(t) = 0 in Eqs. 6d and 6e. In some simulations (Figs. 1517) we disabled one or more of the effects X {C, K, R} by setting X(t) = 0 in Eq. 4a, or alternatively S(t) = 100 or R(t) = 0 in Eq. 6d or 6e, as appropriate.



View larger version (23K):
[in this window]
[in a new window]
 
FIG. 15. Performance during the 21/2-h irregular meal with and without different components of the modulated B15/B16-ARC neuromuscular transform. With the threshold amplitude set at either = 0.132 (A) or = 0.5 (B), the bottom plot shows (black circles and error bars, in some cases invisible because smaller than the circle size) the mean ± SE of the performance M from all 749 cycles of the meal. Top plot shows the percentage of cycles with M > 0. Simulation was run (left to right) without modulation; then with the endogenous modulators released but with the Ca-current effect C only; with the K-current effect K only; with both C and K (i.e., with the complete modulation of contraction size S); with the relaxation-rate effect R only; and finally with both S and R (i.e., with the full modulation of the transform). (The unmodulated and fully modulated cases are the same as presented in Fig. 14.) In all of these cases the inputs to the simulation were the 21/2-h records of the instantaneous firing frequencies of B15 and B16 as in Figs. r10B14. In addition, we ran the simulation with the full modulation but allocating all of the spikes detected in Fig. 10A either entirely to B15 ("All B15") or entirely to B16 ("All B16"). Two horizontal dashed lines in each plot mark the unmodulated and fully modulated levels for reference. Two gray circles in each of A and B show the predicted M if the contributions, relative to the unmodulated case, of C and K to give S, and of S and R to give the full modulation, summed linearly (gray "+" signs and arrows; see Components of the NMT contribute to performance in nonlinear, context- and task-dependent combinations in RESULTS). Statistical significance of the actual M (i.e., the black circles) was tested, separately in A and B, with one-way ANOVA followed by pairwise multiple comparisons using the Holm–Sidak test. In A, all comparisons were significantly different with P < 0.001 except "C + K = S" vs. "C only" and "All B15" vs. "C only", P < 0.01, "All B15" vs. "No modulation", P < 0.05, and "C + K = S" vs. "No modulation", "All B15" vs. "C + K = S", and "K only" vs. "C only", P > 0.05. In B, all comparisons were significantly different with P < 0.001 except "C + K = S" vs. "C only", P < 0.01, "R only" vs. "C + K = S", P < 0.05, and "R only" vs. "No modulation", "C only" vs. "No modulation", "R only" vs. "C only", and "All B16" vs. "C + K = S", P > 0.05.

 
Performance measure

To evaluate the functional significance of the contraction waveform c(t), we use the performance measure m, from Task VI of Brezina and Weiss (2000)Go (see their Appendix F), except modified for the single ARC muscle. Here we give just the essential implementation details; for the biological motivation of the task and performance measure, see RESULTS.

The contraction waveform in any cycle n, [c(t)]n, is required to alternately rise above and fall below a threshold amplitude, c = , during particular phases of the cycle (see Fig. 1). We define a template of 4 timepoints t1 = t* + {phi}, t2 = t* + {phi} + {delta}, t3 = t* + {phi} + {lambda}, and t4 = t* + {phi} + {delta} + {lambda}, where t* is the start time of the current cycle and {phi}, {delta}, and {lambda} are time intervals ≤P. With the regular firing patterns, we set {delta} = 0.15P and {lambda} = 0.5P so as to create in each cycle 2 evenly spaced windows of width 0.15P (Figs. 1, 2B). We then require that c(t1) < , c(t2) > , c(t3) > , and c(t4) < , and if these basic requirements are satisfied, that c(t1) and c(t4) should be as far below and c(t2) and c(t3) as far above as possible. A set of equations that incorporates all of these requirements is

(7a)

(7b)

(7c)

(7d)

(7e)

(7f)
[See the Appendices of Brezina and Weiss (2000)Go for the development of these equations. They are simply computationally compact forms of the required conditions. Equation 7a, for example, is equivalent to A1 = { c(t1) if c(t1) < , 0 if c(t1) ≥ .] In Eq. 7f, we vary {phi} from 0 to P to slide the template of the timepoints t1–t4 along [c(t)]n so as to maximize m'n({phi}). This yields the basic performance measure mn. In this paper, however, we work throughout with the modified performance measure Mn {equiv} mn/P <c>n, where <c