HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) All Versions of this Article: 96/4/2089    most recent 00025.2006v1 Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in Web of Science Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Web of Science (4) Citing Articles via Google Scholar Google Scholar Articles by Tobin, A.-E. Articles by Calabrese, R. L. Search for Related Content PubMed PubMed Citation Articles by Tobin, A.-E. Articles by Calabrese, R. L.

Endogenous and Half-Center Bursting in Morphologically Inspired Models of Leech Heart Interneurons

Department of Biology, Emory University, Atlanta, Georgia

Submitted 10 January 2006; accepted in final form 1 June 2006

	ABSTRACT

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Based on a detailed morphology "Full Model" of a leech heart interneuron, we previously developed a computationally efficient, morphologically inspired "Reduced Model" to expedite tuning the model to produce endogenous bursting and alternating bursting when configured as a half-center oscillator (paired with reciprocally inhibitory synapses). To find conductance density distributions that produce endogenous bursting, we implemented a genetic algorithm automated parameter search. With multiple searches, we found eight parameter sets that produced endogenous bursting in the Reduced Model. When these parameter sets were applied to the Full Model, all produced endogenous bursting, although when the simulation time was extended from 80 to 300 s, only four parameter sets produced sustained bursting in the Reduced Models. All parameter sets produced alternating half-center bursting in the Reduced and Full Models throughout the entire 300 s. When conductance amplitudes were systematically varied for each of the four sustained burster sets, the effects on bursting activity differed, both for the same parameter set in the Reduced and Full Models and for different parameter sets with the same level of morphological detail. This implies that morphological detail can affect burst activity and that these parameter sets may represent different mechanisms for burst generation and/or regulation. We also tested the models with parameter variations that correspond to experimental manipulations. We conclude that, whereas similar output can be achieved with multiple different parameter sets, perturbations such as conductance variations can highlight differences. Additionally, this work demonstrates both the utility and limitations of using simplified models to represent more morphologically accurate models.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
To understand network function in nervous systems, the neurons composing the networks must be identified and characterized at the cellular and subcellular levels. To this end, the most progress has been made in invertebrate motor pattern generating networks. Their few and individually identifiable neurons enable identification of the component neurons in these networks and characterization of their synaptic connections and membrane currents (for reviews, see De Schutter et al. 2005; Harris-Warrick 2002; Marder and Calabrese 1996; Marder et al. 2005). The plethora of data on intrinsic and synaptic conductances and their kinetics from detailed electrophysiological studies has necessitated modeling to elucidate how intrinsic and synaptic conductances interact to shape neuronal and network electrical activity (e.g., Hill et al. 2001; Jezzini et al. 2004; Marder et al. 2005; Nadim and Manor 2000; Nadim et al. 1998). Although such models have been very successful in promoting analysis and understanding of neuron and network activity, realistic morphological details of the component neurons have not been included. This omission is explained partly by the computational costs of detailed compartmental models of neurons and partly by a lack of detailed information about how conductances are distributed within a neuron. Yet structure does matter for neuronal and thus network function. For example, Traub and colleagues (1991) produced an eight-conductance, 19-compartment model of a CA3 hippocampal neuron to analyze several neural activities such as generation of dendritic calcium spikes, depolarizing afterpotentials following spikes, and rhythmic bursting with injected current that simply were not observable in single-compartment models containing the same conductances. This model was subsequently reduced to two compartments, maintaining most of the model activity by segregating fast and slow conductances in two electrically coupled compartments for computational efficiency and dynamical analysis (Pinsky and Rinzel 1994). These studies demonstrate how conductances and their distributions may influence neural activity and provide an inspiration for our exploration of effects of neuronal structure on activity in the heartbeat pattern generating network of the leech.

The network driving the rhythmic heartbeat of the leech is a pattern generating network where not only the neurons and their connections (Thompson and Stent 1976b,c) but also the rhythmic motor output (Thompson and Stent 1976a; Wenning et al. 2004) have been well characterized. Within this central pattern generating network, consisting of seven bilateral pairs of heart interneurons, the interneuron pairs in ganglia 3 and 4 are each linked with reciprocally inhibitory synaptic connections, forming two half-center oscillators that produce alternating bursts of action potentials and pace activity in the network (Peterson 1983). The alternating bursting of each of these oscillators is maintained when their ganglia are isolated (Masino and Calabrese 2002b; Peterson 1983). Extensive electrophysiological characterization of the ionic conductances in the oscillator interneurons has led to the development of a conductance-based single-compartment model of the oscillator interneurons (Hill et al. 2001). An earlier generation of this model (Nadim et al. 1995) guided experiments investigating which conductances are targeted by the endogenous leech peptide FMRFamide to accelerate the rhythmic bursting of the oscillator heart interneurons (Nadim and Calabrese 1997).

Building on our single-compartment model (Hill et al. 2001), we developed a morphologically detailed "Full Model" (871 compartments) of an oscillator heart interneuron, described in the companion paper (Tobin et al. 2006), that permits localized membrane potential changes (such as synaptic inputs) and nonhomogeneous conductance distributions. This hand-tuned model successfully captures the tonic spiking activity of pharmacologically isolated heart interneurons recorded with sharp microelectrodes, but does not express endogenous bursting activity of which the living neurons are capable (Cymbalyuk et al. 2002).

For computational efficiency in network simulations and in parameter searches we reduced the compartment number (nine compartments), preserving passive and active properties and the functional regions. This Reduced Model maps in a systematic way onto the Full Model so that conductance densities chosen for the Reduced Model can be applied to the Full Model (Tobin et al. 2006). Here we use this model for criterion-directed retuning that can be directly applicable to the Full Model. Because of large parameter spaces of active-membrane multicompartmental models, parameter optimization algorithms can be extremely useful for finding parameter sets producing desirable model activity. Vanier and Bower (1999) demonstrated that genetic algorithms are among the most effective parameter-search methods for large parameter spaces [23 parameters in Vanier and Bower (1999)]. This algorithm has been used by others to parameterize compartmental neural models (Hobbs and Hooper 2003; Keren et al. 2005).

Our primary goal in the research described here was to recover parameter sets (conductance density distributions) from our parameter searches that capture endogenous bursting and half-center bursting in our Reduced and Full Model neurons, so that we might study the effects of neuronal morphology and conductance distribution on endogenous and network bursting activity. We implemented an automated genetic algorithm parameter optimization routine to find multiple parameter sets that produced endogenous bursting in a single-neuron Reduced Model. We recovered four sets that produced desired activity and requisite stability, and that were remarkably different, suggesting that different conductance distributions can lead to similar neuronal and network activity, as has been found in crustacean model neurons and networks (Prinz et al. 2003, 2004). The Full Models, both single-cell and half-center, generated from these parameter sets, were similar in activity to their corresponding Reduced Models but more stable to parameter variation, suggesting that neuronal morphology can stabilize endogenous and network activity. We then sought to determine whether, for each of these parameter sets, the models responded appropriately to parameter changes associated with known neuromodulatory effects. Our findings indicate that challenging models with such modulatory influences can be useful in distinguishing different parameterizations of morphological models.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Model parameters

The morphology, compartment functional/transitional groups, and passive parameter settings of the Full Model and Actively Reduced Model, heretofore called the "Reduced Model," are described in the companion paper Tobin et al. (2006). A specific axial resistivity of 250 cm and specific membrane capacitance of 2.2 µF/cm², chosen as described in Tobin et al. (2006), were used for all simulations described here. Leak reversal potential E_Leak and specific membrane resistance R_m were varied with the parameter optimization routine, as described below. We modeled nine active membrane conductances, as previously identified in leech heart interneurons using voltage-clamp procedures (I_Na was identified by current clamp). There are five inward currents: a fast Na current (I_Na), a persistent Na current (I_P), a rapidly inactivating low-threshold Ca current (I_CaF), a slowly inactivating low-threshold Ca current (I_CaS), and a hyperpolarization-activated cation current (I_h). There are three outward currents: a fast, transient K current (I_KA), a delayed-rectifier–like K current (I_K1), and a persistent K current (I_K2). Additionally, we included a spike-mediated synaptic current (I_Syn), with reversal potential –62.5 mV, as described in Hill et al. (2001). In brief, this postsynaptic conductance is triggered by presynaptic spikes and is modeled as the sum of a rising (2 ms) and a falling (11 ms) exponential. The peak of the synaptic conductance is modulated by presynaptic voltage, in accordance with experimental findings. The kinetics and reversal potentials of all the currents are as described in Hill et al. (2001), except for modifications noted in the companion paper (Tobin et al. 2006).

Conductance distributions were assigned to best approximate physiological data and to limit the parameter space for tuning the model. We restrict the Soma to only I_K1, I_K2, and I_KA and the Axon to only I_Na, I_K1, I_K2, and I_KA. In contrast to the Secondary Neurite conductance distributions implemented in the models described in Tobin et al. (2006), the Secondary Neurite compartments of these models have I_h and I_P in addition to I_K1, I_K2, I_KA, I_CaF, and I_CaS. Densities of the Synaptic Compartments are fixed to be identical to those in the Secondary Neurite, except that the synaptic conductance is present only in the Synaptic Compartments and was fixed to the conductance density of 100 mS/cm². The transitional compartments, Neurite–soma, Neurite–axon1, and Neurite–axon2, were designed to smooth transitions between compartments of adjacent functional regions that have different conductance densities. They roughly approximate linear changes in conductance densities that are often distributed in a linear fashion (Bekkers 2000; Hoffman et al. 1997; Korngreen and Sakmann 2000; Magee 1999; Stuart and Hausser 1994) and roughly match impedances to prevent reflection of propagating signals. The Neurite–soma conductance densities are an average of those in the Soma and Neurite. The Neurite–axon1 and Neurite–axon2 have the same density of K conductances as that of the Neurite. For I_h and I_P, the Neurite–axon1 conductance densities are the sum of two thirds of those in the Neurite; the Neurite–axon2 densities are one third of those in the Neurite (the Axon has neither I_h nor I_P). To set the spike initiation zone, fast Na density is set as one half of Axon density in the Neurite, three quarters in Neurite–axon1, and one in Neurite–axon2 to ensure action potentials initiate in the Neurite–axon2 compartments. Experimentally, spontaneous action potentials in heart interneurons are not elicited when the area near the proximal axon/distal neurite is removed by photoablation, thus indicating the region near the Neurite–axon2 compartments is likely the action potential initiation zone in heart interneurons, as it is in the model (Ivanov and Calabrese 1998).

Parameter optimization

We used a genetic algorithm parameter optimization routine implemented in MATLAB (The MathWorks, Natick, MA) by Houck and colleagues (1995). In brief, the genetic algorithm works by creating an initial population of parameter sets where the value of each parameter in a set is randomly chosen (within defined bounds). Each parameter set is rated (assigned a fitness value) by comparing the activity of the model with the parameter set against the desired activity (fitness criteria). In an evolution-inspired fashion, new parameters sets are created by "breeding" the highest-rated parameter sets. The new offspring are then rated, and the best-rated sets are bred to produce new offspring. The breeding continues throughout either a specified number of generations or until the fitness value passes a threshold. The specifics of the options used in calling the genetic algorithm MATLAB command are given in the Supplemental Methods, "Genetic Algorithm Command Options." Because of our large number of parameters, we chose a large initial population size of 10,000 unique parameter sets. Instead of specifying a convergence criterion based on goodness of fit, we specified the optimization routine to produce the best parameter set (evaluated by the fitness criteria described below) found after 200 generations.

When the genetic algorithm converged on a parameter set, to test whether the parameter set could be improved within local variations, we used a simplex search to find the local error minimum (described in LeMasson and Maex 2001). We chose the Nelder–Mead simplex (direct search) method, implemented in MATLAB as the standard function "fmins." Beginning with the parameter set optimized by the genetic algorithm, the simplex function adjusted parameters and ran the GENESIS simulations with the adjusted parameter set to produce new output from the model. This activity was evaluated with the same fitness criteria, described below, as that used for the genetic algorithm parameter optimization, although the final fitness value was multiplied by –1 for the simplex search. This simplex parameter optimization is an iterative process of adjusting parameters and evaluating how the parameter adjustments improve the fitness (bring the fitness value close to zero). We used all default options, as implemented in MATLAB.

Choice of boundaries for conductance densities

The upper and lower bounds for each conductance density were determined based on the conductance densities chosen while previously tuning the Full Model to synaptically isolated tonic spiking activity (Tobin et al. 2006). As described in the companion paper (Tobin et al. 2006), initial conductance densities were chosen to match voltage-clamp amplitudes from published voltage-clamp data, from simulated voltage-clamp in the single-compartment model, or values that represented a compromise between the two. After extensive conductance tuning to match synaptically isolated activity, the resulting parameter set produced conductance amplitudes, as measured by simulated voltage clamp, that were within 14% of initially chosen values, except I_CaS, which was increased by 42%. Because these conductance densities were close to those initially chosen to match voltage-clamp data, except for I_CaS, we used these initial conductance densities as the "default" values to determine the genetic algorithm boundaries. However, because tuning to synaptically isolated activity increased I_CaS densities in the Neurite compartments, we chose for I_CaS a default value that was a compromise between the initial and tuned values. Although these conductance densities were tuned for the Full Model, when they were applied to the Reduced Model the resulting conductance amplitudes, as measured by simulated voltage clamp, did not appreciably change (amplitudes were reduced between 0 and 15%; data not shown). Thus we used these default values to set the boundaries from which the genetic algorithm could choose conductance density values for tuning the Reduced Model. The conductance boundaries for the genetic algorithm were chosen as between about one third of the default value and about threefold the default value (Supplemental Table 1). However, we decreased the lower bound of I_h and I_P in the Secondary Neurite to 0 to allow for no conductance there, in accordance with the Full Model tuned parameter set (Tobin et al. 2006). To constrain action potential height and initiation in the Neurite–axon2 compartment, and to limit parameter space, we fixed I_Na to the default value (3.7 mS/cm²), as tuned in the Full Model (Tobin et al. 2006).

Evaluation of the fitness of a parameter set was based on whether it could produce bursting activity, as defined below, and how closely the burst characteristics matched ideal values, as defined below. Low fitness values were assigned to those parameter sets that most closely matched the desired characteristics. For each parameter set, model activity was simulated with an initial 1-s, –0.15-nA somatic current injection followed by 79 s of run time. Soma voltage was written to an ASCII file and the portion between 20 and 80 s was analyzed as follows. Spike detection and burst discrimination were performed as described in Masino and Calabrese (2002a). In summary, spikes were detected when the voltage crossed a threshold of –30 mV. To prevent double counts of a single spike, a 10-ms refractory period was used (Masino and Calabrese 2002a). Bursts were defined as groups of three or more spikes with 700 ms between each group of spikes. Burst period was calculated as the time between the median spikes of consecutive bursts. Burst duration is defined as the time from the first to the last spike of each burst, and duty cycle is the ratio of burst duration and period for each burst. All parameter sets that did not produce bursting activity were classified as "tonic spiking" or "no spiking," and were assigned the highest (worst) fitness value of 400.

For each parameter set for which bursting was detected, the mean period, SD of period, mean duty cycle, and mean maximum intraburst spike frequency were measured. The criteria for assigning fitness values to parameter sets were developed through successive trials of genetic algorithm runs to find criteria that successfully specified desirable bursting activity. Several types of undesirable bursting activity, described below, were assigned intermediate fitness values of 100, so that they were more likely to contribute to successive generations than the nonbursting parameter sets, but less than those with more desirable bursting activity. Irregular bursting was defined as bursting for which the SD of bursting was >10% of the mean period and was assigned the intermediate fitness value of 100. Nonsustainable bursting was defined as bursting for which there was more than 1.5-fold the duration of the mean period between the last spike of the last burst and the end of the trace (seen as bursting that falls into steady-state membrane potential before the end of the trace). Nonsustainable bursting was assigned a fitness value of 100. If the duty cycle was >80 or <40%, or if the time between bursts was <2.0 s, the parameter set was assigned a fitness value of 100.

Regular bursting that did not fit the above criteria was evaluated based on how closely the period and spike frequency matched the desired values of 8.0 s for period and 13.0 Hz for maximum spike frequency (Supplemental Methods, "Evaluation Criteria for Genetic Algorithm and Simplex Search"). These desired values were chosen to approximate the measured burst characteristics of endogenous bursters. As measured by Cymbalyuk et al. (2002), the period of endogenously bursting heart interneurons ranged from 4.8 to 8.6 s (mean: 6.5 ± 1.0 s) and mean spike frequency ranged from 6.6 to 14.7 Hz (mean: 9.3 ± 2.2 Hz). Our measurements of endogenous bursting in 0.5 mM bicuculline indicated burst periods from 5.0 to 8.4 s, with a mean of 6.6 ± 1.2 s; mean spike frequencies from 5.1 to 11.2 Hz, with a mean of 7.2 ± 2.4 Hz; and maximum spike frequencies from 7.7 to 16.6 Hz, with a mean of 11.3 ± 3.6 Hz (Tobin and Calabrese 2005). This variability in the rhythmic output of the heart interneurons indicates that no one parameter set could describe the range of oscillator heart interneurons, and thus we are seeking multiple parameter sets rather than only one ideal set.

The fitness value was calculated as a sum of the "period fitness" (a weighted difference between the desired period and actual period) and the "spike frequency fitness" (a weighted difference between the desired spike frequency and model spike frequency). For periods <3.0 or >14.0 s, "period fitness" was calculated as 16-fold the absolute value of the difference between the desired period and the period; if the period was between 3.0 and 14.0 s, "period fitness" was eightfold the absolute value of the difference between the desired period and the period. This weighting system was designed to more heavily penalize bursting with periods outside the range normally measured in heart interneurons. Likewise, if the maximum spike frequency was <4.0 Hz (much lower than experimentally measured), an overall fitness value of 100 was assigned. If the maximum spike frequency was >20.0 or <8.0 Hz (outside experimentally measured ranges), "spike frequency fitness" was twice the absolute value of the difference between the desired maximum spike frequency and the maximum spike frequency; if the spike frequency was between 8.0 and 20.0 Hz, the "spike frequency fitness" was the absolute value of the difference between the desired maximum spike frequency and the maximum spike frequency (MATLAB programming code provided in Supplemental Methods¹ ). In addition to penalizing values of burst characteristics that were not within the experimentally observed range, this weighting system also weighed differences in period more heavily than differences of spike frequency, as it was found in trial runs that spike frequency was consistently more closely matched than period.

Sustained endogenous and half-center bursting

We defined sustained bursting as bursting that continued throughout 300 s of simulation time. After parameter sets were chosen for their ability to sustain bursting for 300 s, all further simulations discussed in this text were run to 300 s and burst characteristics were averaged from 240 to 300 s of simulation time. When models were configured as half-centers, one cell received an initial 1-s, –0.15-nA somatic current injection, identical to that in the single-cell models. This asymmetry in starting conditions can cause asymmetric bursting between the two model neurons of the half-center. Unless otherwise noted, activity of models in a half-center was assumed to be symmetric (burst characteristics of noninjected cell were within 10% of those of the injected cell) and the values for the injected cell were quoted for comparison to the single-cell simulations with identical simulated current injection.

Mimicking microelectrode penetration

We modeled the electrode-induced soma leak conductance in a way similar to that described in Cymbalyuk et al. (2002), by introducing into the soma compartment a non-voltage-sensitive leak conductance with reversal potential . The original soma leak conductance and reversal potential values were replaced by the new soma leak conductance g_Leak and reversal potential E_Leak, calculated as (Cymbalyuk et al. 2002; their Eq. 2)

We increased the electrode-induced soma leak conductance until the endogenous bursters spiked tonically, and thus defined a critical value for .

Modeling Na/K pump inhibition

In previous studies exploring the effect of the neuropeptide myomodulin on membrane properties, we measured an inhibition of the Na/K pump, indicated by a myomodulin-induced average shift of 0.058 nA in the voltage-clamp holding current during hyperpolarizing voltage ramps in 0 Ca²⁺, 1.8 mM Mn²⁺ saline with 2 mM Cs⁺ (average holding current in myomodulin: 0.0423 ± 0.007 nA; control: 0.10 ± 0.008 nA; mean ± SE) (Tobin and Calabrese 2005). To simulate the effect of inhibiting the Na/K pump, we injected a constant current ("pump current") into each compartment, scaled by the compartment membrane area. The amplitude of simulated pump current was set for each model by simulating hyperpolarizing voltage-clamp ramps without and with pump current and choosing pump current amplitudes to shift the holding current by 0.058 nA, as empirically observed during myomodulin application. The effect of this pump current on burst characteristics was evaluated by adding the pump current to all compartments during simulated current clamp, with an initial 1-s, –0.15-nA somatic current injection, as in all other current-clamp simulations. The typical mechanisms by which voltage clamp and current clamp are simulated in GENESIS overwrite any additional current injections into the clamped compartment. Thus all simulations using pump current do not include pump current in the soma. Because nonsomatic pump current must be increased to compensate for the missing current in the soma, pump current amplitudes are likely overestimated.

	RESULTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Parameterization of the Reduced Model using a genetic algorithm

The genetic algorithm converged on a parameter set called "Parameter Set 1" (Supplemental Table 2) that produced an endogenously bursting Single-Cell Reduced Model with a period of 4.6 ± 0.1 s, a mean spike frequency of 9.2 ± 0.5 Hz, and a maximum spike frequency of 12.9 ± 1.2 Hz (Table 1). When we configured the Reduced Model into a half-center oscillator with reciprocally inhibitory synapses, the two model cells burst in alternation, with a period of 6.5 ± 0.0 s, mean spike frequency of 17.4 ± 0.3 Hz, and maximum spike frequency of 26.5 ± 0.4 Hz. This result corresponds to an increase in period and spike frequency compared with that in the Single-Cell Model. In heart interneurons bursting in a half-center configuration, we measured periods of 6.1 to 13.2 s (mean 10.2 ± 2.4 s) and mean spike frequencies from 6.8 to 15.3 Hz (mean 10.6 ± 1.5 Hz) (Tobin and Calabrese 2005), and from this same data set we here report maximum spike frequencies of 11.3 to 22.1 Hz (mean 15.8 ± 3.4 Hz). These models' burst characteristics for endogenous bursting are within experimental observations and, for half-center bursting, the period is similar to experimental observations and the spike frequencies are slightly higher.

To test whether the Reduced Model and Full Model are sufficiently similar that the former can be exploited with parameter optimization methods to parameterize the latter, we applied Parameter Set 1 to the Full Model and found that the Full Model did indeed burst endogenously. The burst period in the Full Model, 3.2 ± 0.2 s, was much less than that in the Reduced Model, whereas the mean (10.4 ± 3.2 Hz) and maximum (13.6 ± 6.0 Hz) spike frequencies were greater. When the Full Model was configured into a half-center oscillator, the burst period increased to 6.5 ± 0.03 s, the mean spike frequency to 20.4 ± 0.3 Hz, and the maximum spike frequency to 32.1 ± 0.4 Hz. The differences in period and spike frequency between the Full and Reduced Models demonstrate that the models are not identical but apparently similar enough to share/exchange parameter sets and achieve qualitatively similar activity patterns: endogenous and alternating bursting.

Multiple-parameter sets produce endogenous bursting

Modeling studies have shown that different conductance sets can produce similar bursting activity in models of rhythmically active neurons (Prinz et al. 2003). Moreover, experimental studies indicate that neurons may use a variety of conductance sets to produce the same activity (Golowasch et al. 1999, 2002). Therefore we asked whether a different parameter set could produce endogenous bursting in our Reduced Model. We ran the parameter optimization routine seven more times and assessed the resulting model parameter sets (Supplemental Table 2). From the eight total genetic algorithm runs, the resulting parameter sets all produced endogenous bursting with periods between 4.1 and 6.1 s, mean spike frequencies between 9.2 and 14.0 Hz, and maximum spike frequencies between 12.9 and 19.7 Hz (Table 1). As measured by Cymbalyuk et al. (2002), periods of endogenously bursting heart interneurons ranged from 4.8 to 8.6 s (mean: 6.5 ± 1.0 s) and mean spike frequencies ranged from 6.6 to 14.7 Hz (mean: 9.3 ± 2.2 Hz). Our measurements of endogenous bursting in 0.5 mM bicuculline indicate burst periods from 5.0 to 8.4 (mean: 6.6 ± 1.2 s), mean spike frequencies from 5.1 to 11.2 Hz (mean: 7.2 ± 2.4 Hz), and maximum spike frequencies from 7.7 to 16.6 Hz (mean: 11.3 ± 3.6 Hz) (Tobin and Calabrese 2005). Although a few model parameter sets produced maximum spike frequencies higher than those observed in the living system, overall the burst characteristics produced by these parameter sets are biologically realistic.

All parameter sets that produced endogenous bursting in the Reduced Model produced alternating bursting when the Reduced Model was configured as a half-center oscillator. When applied to the Full Model, all parameter sets produced endogenous and half-center bursting. In both the Reduced and Full Models, with all parameter sets, the period, mean, and maximum spike frequencies of alternating bursting were higher than those of endogenous bursting, in accordance with experimental data (Tobin and Calabrese 2005). In all models, half-center coupling increased the period by an average of 67.4 ± 21.1% and mean spike frequency by 74.5 ± 61.3% (as calculated from activity 20 to 80 s after the start of the simulation).

Sustaining bursting

Multiple parameter sets produced endogenous and half-center bursting in our models. However, to enable the genetic algorithm to converge within a reasonable time, we had limited the model simulation time to 80 s. To single out which of these parameter sets could produce sustained bursting, we ran the Reduced and the Full Models in single and half-center configurations for 300 s with each parameter set. Of the Single-Cell Reduced Models, only those with Parameter Sets 2, 3, 4, and 8 maintained endogenous bursting (Fig. 1), whereas those with Parameter Sets 1, 5, 6, and 7 developed a nonspiking steady resting potential within 300 s. All Half-Center Reduced Models maintained bursting (for examples, see Fig. 1) as did all Single-Cell Full Models and Half-Center Full Models (for examples, see Fig. 2). Because Parameter Sets 2, 3, 4, and 8 were the only ones to maintain endogenous and half-center bursting in both the Reduced and Full Models, we chose these sets as the most promising for realistically representing bursting in oscillator heart interneurons. That only these four parameter sets sustained endogenous bursting in the Full and Reduced Models, whereas all models sustained half-center bursting, supports findings from a single-compartment model of heart interneurons that the half-center configuration stabilizes bursting (Cymbalyuk et al. 2002). In all subsequent studies here only models with these four "robust" parameter sets are considered.

View larger version (26K): [in this window] [in a new window]   FIG. 1. Soma membrane potential is recorded for the Reduced Model with the parameter sets that sustained bursting in all models: Reduced 2 (A), Reduced 3 (B), Reduced 4 (C), and Reduced 8 (D) Models. For A–D: 1, endogenous bursting in Single-Cell Models, shown from 280 to 300 s from start of simulation; 2, alternating bursting in Half-Center Models, shown from 280 to 300 s from start of simulation. For 1 and 2: period and mean spike frequency are calculated for activity from 240 to 300 s after start of simulation.

View larger version (26K): [in this window] [in a new window]   FIG. 2. Soma membrane potential is recorded for the Full Model with the parameter sets that sustained bursting in all models: Full 2 (A), Full 3 (B), Full 4 (C), and Full 8 (D) Models. For A–D: 1, endogenous bursting in Single-Cell Models, shown from 280 to 300 s from start of simulation; 2, alternating bursting in Half-Center Models, from 280 to 300 s from start of simulation. For 1 and 2: period and mean spike frequency are calculated for activity from 240 to 300 s after start of simulation.

Comparing the conductance densities for all parameter sets chosen by the genetic algorithm (including parameter sets where bursting was sustained as well as those where it was not sustained), some trends are noticeable (Fig. 3). For example, although most conductance densities varied throughout the allowable range, neuritic I_KA values were clustered within the top half of the allowable range, suggesting that low neuritic I_KA may be unfavorable for bursting with our desired characteristics. I_K1 densities tended toward a bimodal distribution in the primary and secondary neurites, most noticeably in the axon. Of the inward currents, I_P values clustered in the lower half of the allowable range in both the primary and secondary neurites and, in the secondary neurite, I_h was relatively low, whereas I_CaS clustered in the upper half of the allowable range. Because the four parameter sets that produce sustained bursting likely represent only a few of the potentially large variety of bursters, we can highlight trends that appear, but cannot determine whether these trends represent rules of which distributions can or cannot produce bursting. Whereas a variety of axonal I_K2 conductance densities were observed among all bursters, the sustained bursters clustered near the top of the range. Similarly whereas I_KA and I_K1 spanned the range of allowable values in the secondary neurite, sustained bursters clustered at low values. The average I_h of sustained bursters was lower than that of the nonsustained bursters, as was the membrane resistance.

View larger version (51K): [in this window] [in a new window]   FIG. 3. For all parameter sets chosen by the parameter optimization routine for producing endogenous bursting in the Reduced Model, conductance densities are plotted, with black crosses, by compartment type. ELeak and Rm, which had the same value across all compartments, are plotted under the heading "homogeneous." Those parameter sets that exhibited sustained bursting in all models (Single-Cell, Half-Center, Reduced, and Full) are shown as black crosses highlighted by colored circles; legend displays parameter set number. Gray bars indicate the boundaries used by the parameter optimization routine (see Supplemental Table 1).

View larger version (42K): [in this window] [in a new window]   FIG. 4. Conductance density distributions of parameter sets 2, 3, 4, and 8 are compared. Icon represents compartment types as labeled. For each conductance, densities are normalized to the maximum density occurring in any of the 8 parameter sets from the genetic algorithm (regardless of compartment type). Thus densities are scaled with respect to only those values that are associated with endogenous bursting to enable comparison of densities within a model and across models.

To determine whether the bursting observed in both the single-cell and half-center configurations of the Reduced and Full Models was sensitive to small parameter variations, which would suggest "small hypervolume" in parameter space pinpointed by the genetic algorithm, we did a sensitivity analysis (Calabrese et al. 2001). Sequentially, we varied by ±5% the values of R_m and E_Leak and the density of each of the seven conductances, in the compartments where they were present, and we measured whether the Reduced and Full Models with Parameter Sets 2, 3, 4, and 8 could still produce endogenous and alternating bursting for 300 s.

When the conductance of I_K1 was varied, endogenous bursting was sustained for ±5% changes in Single-Cell Reduced 2, 3, and 8 and was not sustained in Single-Cell Reduced 4. For all other conductances, the Single-Cell Reduced 2, 3, 4, and 8 models were unable to sustain endogenous bursting when outward currents (I_K2, I_KA) were increased, when I_Leak was increased (modeled as a decrease in R_m), or when E_Leak was hyperpolarized. Similarly, decreasing inward currents (I_h, I_P, I_CaF, and I_CaS) resulted in nonspiking traces. The Single-Cell Reduced 8 differed from this trend only in that it was able to sustain bursting for both ±5% changes in I_CaF. All Half-Center Reduced Models were able to sustain bursting for all variations except that when E_Leak was hyperpolarized, the Half-Center Reduced 2, 3, and 4 stopped spiking. When the conductance variations were applied to the Full Model, all variations elicited endogenous bursting except for changes in E_Leak, for which hyperpolarizing caused no spiking in Single-Cell Full 2, 3, and 8, and depolarizing caused tonic spiking in Single-Cell Full 8. However, all Half-Center Full Models sustained bursting for all variations.

In general, the Single-Cell Reduced Models could not sustain endogenous bursting for hyperpolarizing changes: when outward currents were increased, E_Leak hyperpolarized or inward currents decreased. Half-Center Reduced Models and Single-Cell Full Models sustained bursting for all variations except E_Leak. Half-center bursting in the Full Models was robust to all changes. In experiments, endogenous bursting in heart interneurons is not always regular and simulations with the single-compartment model indicate it is very sensitive to leak parameters, whereas half-center bursting, in experiments and simulations, is far more regular (Cymbalyuk et al. 2002). Both the Reduced and Full Models corroborate these experimental and model results. Furthermore, these results indicate that bursting is more robust in the Full Model than in the Reduced Model, even though these parameter sets were tuned to the Reduced Model.

Assessing changes in burst characteristics with parameter variations

To assess changes in burst characteristics with parameter variations, we varied parameters by ±20% (to accentuate the small changes seen in 5% variations) in both Reduced and Full Models and calculated the percentage changes in average period and mean spike frequency. Full data sets are presented in Supplemental Tables 3–6. Because these burst characteristics have some variability within a simulated run, only those changes that are larger than the percentage SD of the burst characteristics measured from a simulation with the original parameter set are considered relevant. Because this criterion is less stringent than a statistical analysis (ANOVA with post hoc testing) for each variation, some statistically insignificant changes may be considered here.

Increasing the amplitude of the parameter variation to 20% led to a cessation of bursting for some models and some parameters, especially increases in outward currents (I_K1, I_K2, I_KA) and decreases of burst-supporting inward currents (I_P, I_CaF, I_CaS) or I_h. Single-Cell Reduced 4, on the other hand, unable to burst with a 5% decrease in I_K1, exhibited bursting with a 20% decrease in I_K1 (Fig. 5 B2). Otherwise, the responses of the models to 20% variations followed patterns similar to those with 5% variations (Figs. 5–8). In describing the responses of the models to the parameter variations here, we limit ourselves to models where bursting was maintained and focus on the general trends, mentioning exceptions only where they appear pertinent.

View larger version (30K): [in this window] [in a new window]   FIG. 5. Changes in burst characteristics arising from 20% conductance variations are measured in Single-Cell Reduced Models. Each conductance density, including Rm and ELeak, was varied by +20% (A1 and A2) and –20% (B1 and B2) in all compartments where the conductance was present, and changes in period (A1 and B1) and spike frequency (A2 and B2) were measured. Open shapes correspond to changes that were deemed nonrelevant, defined as the changes being within the SD of burst characteristics from the original parameter set. Where conductance variations stopped bursting, activity is indicated as "No spiking" or "Tonic spiking."

View larger version (31K): [in this window] [in a new window]   FIG. 8. Changes in burst characteristics arising from 20% conductance variations are measured in Half-Center Full Models. Each conductance density, including Rm and ELeak, was varied by +20% (A1 and A2) and –20% (B1 and B2) in all compartments where the conductance was present, and changes in period (A1 and B1) and spike frequency (A2 and B2) were measured. Open shapes correspond to changes that were deemed nonrelevant, defined as the changes being within the SD of burst characteristics from the original parameter set. Where conductance variations stopped bursting, activity is indicated as "No spiking" or "Tonic spiking."

The Half-Center Reduced Models generally exhibited qualitatively similar changes to one another and to their corresponding Single-Cell Models for all parameter variations (Fig. 6). A notable exception was for variation in I_h, where now increasing I_h decreased period and decreasing I_h increased period according to expectation from the single-compartment model (Hill et al. 2001).

View larger version (31K): [in this window] [in a new window]   FIG. 6. Changes in burst characteristics arising from 20% conductance variations are measured in Half-Center Reduced Models. Each conductance density, including Rm and ELeak, was varied by +20% (A1 and A2) and –20% (B1 and B2) in all compartments where the conductance was present, and changes in period (A1 and B1) and spike frequency (A2 and B2) were measured. Open shapes correspond to changes that were deemed nonrelevant, defined as the changes being within the SD of burst characteristics from the original parameter set. Where conductance variations stopped bursting, activity is indicated as "No spiking" or "Tonic spiking."

View larger version (32K): [in this window] [in a new window]   FIG. 7. Changes in burst characteristics arising from 20% conductance variations are measured in Single-Cell Full Models. Each conductance density, including Rm and ELeak, was varied by +20% (A1 and A2) and –20% (B1 and B2) in all compartments where the conductance was present, and changes in period (A1 and B1) and spike frequency (A2 and B2) were measured. Open shapes correspond to changes that were deemed nonrelevant, defined as the changes being within the SD of burst characteristics from the original parameter set. Where conductance variations stopped bursting, activity is indicated as "No spiking" or "Tonic spiking."

The effect of microelectrode penetration in Reduced and Full Models: endogenous bursting versus tonic spiking

Synaptically isolated oscillator heart interneurons can burst endogenously when recorded extracellularly, but when impaled by a sharp microelectrode, they spike tonically (Cymbalyuk et al. 2002). We tested whether the Reduced and Full Models would reproduce this phenomenon. We introduced a nonspecific (E_Leak = 0 mV) "electrode-induced" leak conductance , calculated as in Cymbalyuk et al. (2002) (see METHODS), into the soma compartment of each neuron model and increased the electrode-induced leak until the endogenous bursters spiked tonically. All Reduced and Full Models had a critical value of this electrode-induced leak for which the endogenous bursting became tonic spiking [Reduced 3 example: Fig. 9A; for comparison to biological data, see Fig. 2B in the companion paper Tobin et al. (2006)]. Normally connected oscillator heart interneurons (i.e., in half-center configuration) burst in alternation, even when both cells are penetrated with microelectrodes. Likewise, all Half-Center Models, with either Reduced or Full morphology, also continued bursting in alternation with the critical value of the somatic electrode-induced leak conductance [Reduced 3 example: Fig. 9B; for comparison to biological data, see Fig. 2A in the companion paper Tobin et al. (2006)].

View larger version (27K): [in this window] [in a new window]   FIG. 9. Soma voltage is recorded for the Reduced 3 Model in the (A) Single-Cell and (B) Half-Center configurations while mimicking a microelectrode-induced soma leak.

The effects of blocking I_h with Cs⁺ in the Half-Center Reduced and Full Models

Only a few types of conductance variations are useful for comparing the burst regulation of the models to experimental data; for example, blocking Na⁺ currents stops bursting in interneurons. We chose to test the effect of Cs⁺, which blocks I_h in heart interneurons (Angstadt and Calabrese 1989). Adding 2 mM Cs⁺ to the saline increases period (Masino and Calabrese 2002b; Tobin and Calabrese 2005) and decreases spike frequency of extracellularly recorded heart interneurons (Tobin and Calabrese 2005). To test whether our Half-Center Reduced and Full Models could replicate these effects, we mimicked the biophysical effects of Cs⁺ by setting h conductances to 0. As indicated by the parameter variations, decreasing I_h by 20% increased period and spike frequency in most Half-Center Reduced Models (Fig. 6 and Supplemental Table 4) and Half-Center Full Models (Fig. 8 and Supplemental Table 6). Completely blocking I_h, however, stopped bursting and spiking in all Half-Center Reduced and Full Models. Either our models do not fully represent the response of heart interneurons to the Cs⁺-induced I_h block or there are other membrane properties that Cs⁺ affects, for which we have not accounted. For example, Cs⁺ may enhance the Na/K pump in leech neurons (Catarsi and Brunelli 1991; Skou 1960, 1965). Potential nonspecific effects of Cs⁺ are indicated by the nonmonotonic nature of the dose–response relation of Cs⁺ on period; period reduction peaks at 2 mM and increasing concentrations decrease the effect (Masino and Calabrese 2002b). Moreover,intracellularly recorded heart interneurons fire tonically in Cs⁺-containing saline (Angstadt and Calabrese 1989), indicating that leak and I_h interact in a way that the models might not yet capture. The drug ZD7288, shown to block I_h selectively in vertebrate preparations (Gasparini and DiFrancesco 1997), does not block I_h in leech heart interneurons (Tobin and Calabrese 2005), precluding a more specific block at present.

The effects of myomodulin in the Single-Cell and Half-Center Reduced and Full Models

Myomodulin decreases period and increases spike frequency when bath applied at 1 µM to normally connected and synaptically isolated heart interneurons (Tobin and Calabrese 2005). We demonstrated that myomodulin enhances the maximum conductance of I_h by 30% and likely inhibits the Na/K pump (Tobin and Calabrese 2005). In the Full and Reduced Single-Cell and Half-Center Models, we simulated the effects of myomodulin on membrane properties to test whether they could account for the observed changes in burst period and spike frequency of endogenous and half-center bursting in heart interneurons.

To mimic the 30% increase of I_h caused by myomodulin, we increased I_h conductance uniformly until we measured a 30% increase in simulated somatic voltage clamp (for Parameter Sets 2, 3, 4, and 8, I_h conductance densities were increased by 30, 35, 50, and 30%, respectively, in Reduced Models and by 38, 33, 48, and 30%, respectively, in the Full Models). The Na/K pump inhibition, measured experimentally as a inward current shift in the linear portion of voltage ramp currents (Tobin and Calabrese 2005), was mimicked by current injections of 45, 54, 50, and 40 A/cm² in the Reduced 2, 3, 4, and 8 Models, respectively, and by 46, 55, 52, and 41 A/cm² in the Full 2, 3, 4, and 8 Models, respectively. In experiments, values for the inhibited pump current were measured only in Cs⁺ saline in which the pump may already be enhanced by the Cs⁺ (Catarsi and Brunelli 1991; Skou 1960, 1965). Thus we used these values as only an estimate of the effect myomodulin may have in normal saline, and for Single-Cell Models where the full amplitude of pump current caused tonic spiking (Reduced 8, Full 3, Full 8), we use half the amplitude of the inhibited pump current for both the Single-Cell and Half-Center Models. For Reduced 3, half-amplitude pump current was used only when both I_h and pump current were combined. For each of the modulatory effects, separately and combined, we measured burst characteristics for the Single-Cell and Half-Center Full and Reduced Models (Fig. 10 and Supplemental Table 7).

View larger version (23K): [in this window] [in a new window]   FIG. 10. Effects of the neuromodulator myomodulin are mimicked in the models by increasing Ih conductance by 30% and inhibiting the Na/K pump (simulated by current injection; see RESULTS). Each effect was applied separately and combined to qualitatively assess the individual and combined effects of the biophysical changes on period and spike frequency for (A) Single-Cell Reduced Models, (B) Single-Cell Full Models, (C) Half-Center Reduced Models, and (D) Half-Center Full Models. Where the full amplitude of pump current caused tonic spiking in Single-Cell Models (Reduced 8, Full 3, and Full 8), we use half the amplitude of the pump current for both the Single-Cell and Half-Center Models. Half-amplitude pump current was used for Reduced 3 only when both Ih and pump current were combined. B: Single-Cell Full 8 exhibited tonic spiking even when the pump current was reduced to half-amplitude. D: *Half-Center Full 4 exhibited asymmetric bursting such that the spike frequency of one model neuron decreased by 27.1%, whereas the spike frequency of the other model neuron increased by 17.4%. Periods of each model neuron decreased by 17.4 and 17.5%, respectively. Heart interneuron response to myomodulin is plotted with open triangles with SD bars. Response of synaptically isolated (bicuculline) heart interneurons to myomodulin is plotted with Single-Cell Models under the heading "both." Response of synaptically intact heart interneurons to myomodulin is plotted with Half-Center Models: under the heading "Ipump" when Ih is blocked by Cs+ and under the heading "both" when Ih and Ipump are intact (normal saline).

In the Single-Cell Full Models (Fig. 10B), increasing I_h, adding pump current, and both effects combined decreased period and spike frequency (except for Full 8, where even half the pump current caused tonic spiking). In each case, the effects were sublinearly additive. Whereas the decrease in period matches that measured experimentally, the decrease in spike frequency does not (Tobin and Calabrese 2005) (Fig. 10B). In the Half-Center Full Models (Fig. 10D), increasing I_h, adding pump current, and both effects combined decreased period and increased spike frequency, as seen experimentally with myomodulin (Fig. 10B) (Tobin and Calabrese 2005) (except for Half-Center Full 2, where increasing I_h slightly decreased spike frequency, and Half-Center Full 4, where the I_h and pump effects combined increased spike frequency in one model neuron and decreased it in the other). These effects were sublinearly additive as well. In general, the Single-Cell and Half-Center Full Models indicate that each of the biophysical effects of myomodulin, separately and combined, should decrease spike frequency in endogenously bursting cells, increase spike frequency of normally connected cells, and decrease period in either configuration.

	DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Using an automated parameter optimization routine, we recovered four parameter sets (conductance density distributions) that capture stable endogenous bursting and half-center bursting in our Reduced and corresponding Full Model of a leech oscillator heart interneuron. We then subjected the models with these parameter sets to parameter variations and tests corresponding to experimental measures: microelectrode-induced soma leak (Cymbalyuk et al. 2002), Cs⁺ application, and application of the neuromodulator myomodulin (Tobin and Calabrese 2005).

Choice of optimization algorithm and convergence criteria for model tuning

Our Full Model of a heart interneuron is too computationally intensive for effective automated parameter optimization, but the design of a corresponding Reduced Model with single compartments that map onto functional regions of the Full Model (Tobin et al. 2006) allowed us to use the Reduced Model as a proxy for the Full Model in parameter searches. Because of the high dimensionality of our Reduced Model (22 parameters), we chose a genetic algorithm for automated parameter optimization; this method has been shown to be among the most effective parameter search methods for large parameter spaces (Vanier and Bower 1999). We designed the fitness criteria to select for endogenous bursting with burst characteristics (period, spike frequency, and duty cycle) to match extracellularly recorded endogenous bursting in pharmacologically isolated heart interneurons (Cymbalyuk et al. 2002; Tobin and Calabrese 2005). We chose these characteristics as appropriate to the neurons/network studied. First, these measures are easily repeatable and represent native activity that is as unbiased as possible by the recording method. Second, because we sought to use these models in network simulations of the heartbeat pattern generator, period and spike frequency are important output parameters in the living system (Norris et al. 2006).

The specific target values chosen for the genetic algorithm were arrived at empirically. Trial runs indicated that constraining maximum spike frequency, rather than mean spike frequency, would more likely yield realistic intraburst spike frequency profiles. Such trial runs also indicated that spike frequency was more easily matched than period, leading us to weight period matching more heavily. Moreover, our experience with hand-tuning the single-compartment model on which the present morphological models are based (Hill et al. 2001) suggested that obtaining long periods was more difficult than short ones, so we picked a period target near the end of the observed range in living pharmacologically isolated extracellularly recorded heart interneurons (4.8–8.6 s; Cymbalyuk at al. 2002; Tobin and Calabrese 2005). Duty cycle was only loosely constrained, with equal fitness assigned for a broad range of duty cycle 40–80% because we reasoned that in half-center bursting strong mutual inhibition would constrain duty cycle to the nearly 50% observed in the living system. Because we used three criteria (period, maximum spike frequency, and duty cycle) for assigning fitness values, we did not expect to match all values exactly and we assigned weights for period and spike frequency in a two-tiered fashion, such that those values outside the biological range were penalized more heavily in their deviation from the desired value than those within the range of biological variability. Not unexpectedly we obtained model activity that did not closely match any of the target values, but exhibited biologically realistic activity. During the selection process, bursting activity was assessed over 80 s of simulation time for computational efficiency. This decision necessitated a further assay of stability; each of the eight sets independently chosen by the genetic algorithm was run for 300 s of simulation time. Only four sets produced sustained bursting and were thus chosen for further analysis. For these four sets the period was uniformly below and the maximum spike frequency uniformly above the target values and duty cycle was uniformly at the low end of the target range (Table 1). Judging from this limited number of searches, it appears that the searches converge on a competitive compromise of the target values.

Our experience suggests that in parameterizing a neuronal model, such optimization criteria and weighting should be chosen with regard to the system of interest and the desired goal. Characteristics other than those chosen here, such as amplitude of slow waveform or number of spikes per burst, in rhythmically active neurons, or spike timing, spike shape, or response to stimulation in input driven neurons may be considered important characteristics defining the activity and function of a neuron. Additionally, where possible, recordings from multiple sites in a neuron may be particularly useful in constraining multicompartmental models with nonhomogeneous conductance distributions (Keren et al. 2005).

Use of a Reduced Model for parameter optimization of a Full Model

The Reduced Model of an oscillator heart interneuron was created as a computationally efficient representation of the general morphology and the passive and active properties of the detailed morphology Full Model. The geometry of the Reduced Model was constrained to match both passive and active properties of the Full Model. Because these properties are determined by the conductance densities and distributions of the Reduced and Full Model, the geometry of the Reduced Model was dependent, to some degree, on the criteria used during the reduction (Tobin et al. 2006). Here, we generated eight parameter sets (conductance density distributions) in the Reduced Model that produce similar activity when implemented in the Full Model (Figs. 1 and 2). This observation suggests that the Reduced Model captures essential features of the Full Model and that our assertion that it could be used to tune the Full Model is supported. The Reduced and Full Models are not completely congruent, however, as evidenced by the four parameter sets that did not produce sustained endogenous bursting in the Reduced Model (and were thus eliminated from further analysis) but did in the Full Model. Although parameter sets produce bursting with somewhat different burst characteristics in the Reduced Model than in the Full Model, these sets can provide a useful starting point for further tuning in the Full Model.

General implications of multiple parameter sets producing similar activity

When characterizing K currents in an identified neuron in the rhythmic pyloric network of the crab stomatogastric nervous system, Golowasch et al. (1999) showed that homologous neurons from different animals differ greatly in the amplitudes of the conductances they possess yet produce similar activity. Moreover, modeling studies demonstrated that a model neuron with conductance amplitudes set to experimentally measured averages may not replicate the electrical activity of that neuron type (Golowasch et al. 2002). Recent studies in this system demonstrate that variation in specific conductance amplitudes is correlated with mRNA expression of the corresponding specific channel genes (Schulz et al. 2006). These studies suggest that the experimentally measured differences in conductance amplitudes reflect true differences and not simply measurement error. Therefore our finding of multiple parameter sets producing similar activity may indicate that similar variation occurs in living heart interneurons. By varying conductances in a three-cell model of the pyloric network, Prinz et al. (2004) demonstrated that multiple parameter sets could produce virtually indistinguishable model network output. Similar to the work of Prinz and colleagues, we found multiple parameter sets that produced similar activity, in this case endogenous and half-center bursting with specified characteristics. As emphasized by Hooper (2004), such findings demonstrate similar apparent steady-state output, but do not indicate how models would respond to perturbations. Indeed the parameter variations we imposed on our models have uncovered multiple qualitative and quantitative differences between the different parameter sets and between levels of morphological detail in our models (Figs. 5–8).

Significance of multiple parameter sets producing similar activity for understanding heart interneuron function

The four parameter sets that capture sustained endogenous bursting and half-center bursting in our Reduced and corresponding Full Models of a leech oscillator heart interneuron can seem bewilderingly different in both their conductance density distributions and responses to parameter variation and perturbation (Figs. 3–10). Yet some simple patterns emerge that provide some insight into heart interneuron function. As described in detail in Comparing conductance densities and distributions under RESULTS, there is a rough tendency for slow outward currents and inward currents that support burst formation to balance. This tendency is most easily seen in the secondary neurite, where high I_CaF and/or I_CaS correspond to high I_K2 (Fig. 4). This simple expectation on which the theory of activity-dependent homeostasis in neurons and networks is built (LeMasson et al. 1993; Turrigiano et al. 1994) is thus met in our models. I_h and I_P are low in all four parameter sets (Fig. 3), suggesting that these inward currents that are active during the interburst interval must be manageably small to sustain endogenous bursting; some I_h is nonetheless required or endogenous bursting ceases (see The effects of blocking I_h with Cs⁺ in the Half-Center Reduced and Full Models under RESULTS). Only detailed analysis of current flows within the models with these parameter sets can reveal how bursting arises with each, but the differences observed in responses to perturbation indicate different mechanisms.

Significance of parameter sets producing dissimilar changes to conductance variation

Although overall, the four stable parameter sets for the Single-Cell Reduced Model exhibited similar changes to conductance variations, each set differed from the others for at least one conductance variation. We expect that, if the parameter sets produce bursting by the same mechanism, they would respond in qualitatively similar manners. For example, previous analysis of the single-compartment single-cell and half-center heart interneuron model shows that I_h provides a "predictable" potential mechanism for regulating burst period by providing a depolarizing drive during the interburst interval (Cymbalyuk et al. 2002; Hill et al. 2001; Sorensen et al. 2004). A different role in bursting is needed to explain how increasing I_h would increase period (e.g., Parameter Sets 2 and 4, Fig. 5A1). Such responses are observed in the rhythmic respiratory network of the rat, where applying Cs⁺ or ZD7288 to block I_h decreases the period between inspiratory bursts (Thoby-Brisson et al. 2000), albeit the underlying mechanism is unknown.

When each dissimilar parameter set was applied to the Half-Center Model, for most conductance variations the Half-Center Models reacted similarly (Figs. 6 and 8), indicating that strong inhibitory synaptic interactions impose a bursting mechanism that causes pairs composed of dissimilar endogenous bursters to react similarly. Such a mechanism would be advantageous in the living system, where cells could have different mechanisms for endogenous bursting, but network interactions ensure the network responds predictably to modulator influences. Further experiments are necessary to determine how half-center and endogenously bursting heart interneurons respond to conductance variations, and whether qualitative differences arise when comparing the response of synaptically isolated oscillator interneurons from different animals or when comparing synaptically isolated oscillator interneurons to those in half-center configuration.

Implications of our results for the parameter space of heart interneuron bursting

Using the genetic algorithm multiple times allowed us to uncover multiple parameter sets that produced similar activity. The genetic algorithm samples many points in parameter space, initially evaluating a population of parameter sets randomly chosen from throughout parameter space. However, it is designed to narrow in on the single parameter set that produces the most desirable activity, and thus it is an inappropriate tool for mapping parameter space to uncover the multiplicity of solutions that can produce similar behavior. Therefore we do not purport that our four sets represent the range of potential "solutions" or even the most realistic solutions possible. In fact, given the dimensionality of the parameter space (22) we surmise that we have scarcely begun to plumb its depths. A more appropriate approach to delimiting the space would be to systematically vary parameters to determine the region(s) of parameter space containing the desired activity, ultimately producing a searchable database (Goldman et al. 2001; Prinz et al. 2003, 2004). Such previous studies explored at maximum 10 parameters, with sparse sampling (maximum of six values per parameter), our model with 22 parameters is likely too unwieldy for such analysis. However, by developing a group of multiple similar models from these and further genetic algorithm searches, we may begin to uncover parameters that appear to be most or least variable and thus constrain which parameters are best for systematic variations. The present results hint at potentially useful patterns (Fig. 3). For example, in the secondary neurite, I_CaS appears to be constrained to high values (with respect to our bounds), I_h and I_P to low values, whereas in Neurite, I_KA appears high, I_P low, and Axon I_K1 bimodal. Values that are more tightly constrained may be more important in determining the burst activity and characteristics of heart interneurons and thus may represent the most favorable objects of systematic parameter variation and database construction.

The effects of myomodulin in the models

Challenging a model with a parameter variation corresponding to a natural neuromodulation is an appealing test of model validity. Myomodulin increases I_h and decreases the Na/K pump, with associated decreases in period and increases in spike frequency (Tobin and Calabrese 2005). Mimicking these changes in all the Single-Cell Full Models decreased period, as expected, but also decreased spike frequency (Fig. 10B). In the Single-Cell Reduced Models, these changes increased spike frequency in nearly all models, but had mixed effects on period (Fig. 10A). Nearly all half-center oscillator models show that the membrane properties changed by myomodulin can account for the experimentally measured effects of myomodulin on burst properties both when I_h is present and when it is blocked by Cs⁺ (Fig. 10, C and D). When I_h is blocked, we hypothesize that the decrease in period and increase in spike frequency elicited by myomodulin are caused principally by myomodulin's inhibition of the Na/K pump. All but one Half-Center Model showed that inhibiting the pump can cause the observed change (qualitatively) in burst properties (Fig. 10, C and D and Supplemental Table 7). When I_h was increased (as is caused by myomodulin), nearly all half-center models show a decrease in period and a slight increase in spike frequency of oscillator interneurons bursting, corroborating the idea that increasing I_h in the living neurons would have these effects (Fig. 10, C and D and Supplemental Table 7). Furthermore, in most models, the combination of increasing I_h and inhibiting the pump in the half-center oscillator models decreased period and increased spike frequency more than either effect alone, supporting the idea that these two mechanisms work together to produce the myomodulin-induced change in bursting (Fig. 10 and Supplemental Table 7).

Although none of the sets perfectly captures myomodulin's effects, this challenge to the model identifies Parameter Set 3 as responding most in line with the observed effects of myomodulin, failing only to increase its spike frequency as a Single-Cell Full Model. Such challenges may be generally useful in distinguishing multiple parameter sets identified in modeling studies.

Influence of morphology on model activity

Electrical activity in a neuron with extended structure is a local balance of membrane and axial currents that are in turn affected by local conductance densities, membrane and axial resistances, and membrane capacitance. Each compartment in the Reduced Model approximately matches these properties, by preserving surface area and soma current load over large functional regions (e.g., neurite, axon, etc.) (Tobin et al. 2006), but does not explicitly match the ratio of the total membrane and axial current of the multiple compartments in a functional region. Stated alternatively, the total membrane current flowing into and axial current flowing to and from a function region is not the same as that of its representative reduced compartment. Thus the primary differences between a compartment in the Reduced Model and the functional region it represents in the Full Model is one of effective input impedance that determines the balance of membrane and axial current flow for the compartment/region. Because different parameter sets produced similar bursting activity in the Single-Cell Full and Reduced Models, it is possible that these considerations may not have an appreciable effect on the slow-wave endogenous bursting activity we seek, where spike frequency, but not precise spike timing, is important. However, endogenous bursting appears to be quite sensitive to the ratio of conductances in a compartment because even 5% variations in some conductances bring some parameter sets out of the bursting regime. Moreover, endogenous bursting is more robust to large parameter variations in the Full Model (over the different parameter sets) than in the corresponding Reduced Model (Figs. 5 and 7).

Interestingly, we find that in four cases, parameter sets that produced only transient bursting in the Single-Cell Reduced Model (and were thus eliminated from further analysis) could produce sustained endogenous bursting in the Full Model, suggesting that extended neuronal morphology as represented in the Full Model stabilized endogenous bursting activity. For some conductance variations, endogenous burst characteristics changed differently in Full and corresponding Reduced Models (e.g., compare Figs. 5 and 7). This puzzling result suggests that there are some conductance ratios that produce a form of bursting that is more sensitive to the effects of input impedance. Such results deserve more detailed analysis to understand how burst generation, and neuronal activity in general, is affected by detailed versus simplified model morphology.

For most parameter variations, however, the Reduced, Full (Figs. 6, 7, and 8), and even single-compartment (Hill et al. 2001) Half-Center Models reacted similarly to conductance variations. For example, I_CaS is similarly effective in regulating period and spike frequency in a predictable manner. These similarities again indicate an overriding influence of strong mutual inhibition and potentially general rules of how conductances influence such half-center bursting in heart interneurons—rules that are not dependent on specific parameter sets in the Full and Reduced Models or even the detail of morphological representation. However, the qualitative differences in reactions to certain parameter variations (seen especially for Single-Cell Models) and the quantitative difference in sensitivity to the parameter variation even when reactions are similar between Reduced and Full Models (Figs. 5–8) does indicate an influence of morphological representation and spatial distribution of conductances. For example, if the slow-wave oscillations are generated in a different area (perhaps the Secondary Neurite compartments) than the action potentials (initiated in Neurite–axon2), nonhomogeneous conductance distribution could cause the period to be differently sensitive to a particular conductance change than the spike frequency. This effect is not possible in an isopotential, single-compartment model.

In conclusion, we have successfully used a systematically constructed Reduced Model of a leech heart interneuron to find multiple parameter sets that produce endogenous and half-center bursting both in the Reduced Model and in the morphologically detailed Full Model from which it is derived. These parameter sets demonstrate that more than one mechanism of bursting may be possible in living oscillator heart interneurons. Furthermore, our findings indicate that the level of morphological detail represented in the model changes the sensitivity of bursting to parameter variations and may change the mechanisms by which bursting is generated and/or regulated. The responses to parameter variation observed in these models, including mimicking myomodulin-induced modulations of I_h and the Na/K pump, differentiate between them and, although indicating that each of our Full Models has flaws, show that we have made progress developing a morphological model that has experimental predictability. Additional parameter searching and possibly model modifications based on this analysis will certainly further this progress. The Reduced Model will continue to be an invaluable tool in future searches and potentially in network simulations.

Despite their drawbacks, these morphologically inspired models can now be used in future studies to evaluate the effect of morphology and conductance distribution on cell behavior and somatic electrophysiology recording. Such studies may include evaluation of how space-clamp error could affect measurements of conductance kinetics and whether morphologically distributed conductances lead to greater stability in the face of random variation.

	GRANTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
This work was supported by the National Institute of Neurological Disorders and Stroke National Research Service Award Predoctoral Fellowship NS-42983 to A.-E. Tobin and by NS-24072 to R. L. Calabrese.

	ACKNOWLEDGMENTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
We thank Dr. A. Wenning for critically reading the manuscript.

Present address of A.-E. Tobin: Volen Center for Complex Systems, MS 013, Brandeis University, 415 South Street, Waltham, MA 02454.

	FOOTNOTES

 
The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked "advertisement" in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

¹ The online version of this article contains supplemental data.

Address for reprint requests and other correspondence: R. L. Calabrese, Dept. of Biology, Emory University, 1510 Clifton Road N.E., Atlanta, GA 30322 (E-mail: rcalabre{at}biology.emory.edu)

	REFERENCES

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION GRANTS ACKNOWLEDGMENTS REFERENCES

 
Angstadt JD and Calabrese RL. A hyperpolarization-activated inward current in heart interneurons of the medicinal leech. J Neurosci 9: 2846–2857, 1989.[Abstract]

Bekkers JM. Distribution of slow AHP channels on hippocampal CA1 pyramidal neurons. J Neurophysiol 83: 1756–1759, 2000.[Abstract/Free Full Text]

Calabrese RL, Hill AA, and van Hooser SD. Realistic modeling of small neuronal circuits. In: Computational Neuroscience: Realistic Modeling for Experimentalists, edited by De Schutter E. Boca Raton, FL: CRC Press, 2001, p. 259–288.

Catarsi S and Brunelli M. Serotonin depresses the after-hyperpolarization through the inhibition of the Na⁺/K⁺ electrogenic pump in T sensory neurones of the leech. J Exp Biol 155: 261–273, 1991.[Abstract/Free Full Text]

Cymbalyuk GS, Gaudry Q, Masino MA, and Calabrese RL. Bursting in leech heart interneurons: cell-autonomous and network-based mechanisms. J Neurosci 22: 10580–10592, 2002.[Abstract/Free Full Text]

De Schutter E, Ekeberg Ö, Kotaleski JH, Achard P, and Lansner A. Biophysically detailed modelling of microcircuits and beyond. Trends Neurosci 28: 562–569, 2005.[CrossRef][Web of Science][Medline]

Gasparini S and DiFrancesco D. Action of the hyperpolarization-activated current (Ih) blocker ZD 7288 in hippocampal CA1 neurons. Pfluegers Arch 435: 99–106, 1997.[CrossRef][Web of Science][Medline]

Goldman MS, Golowasch J, Marder E, and Abbott LF. Global Structure, Robustness, and Modulation of Neuronal Models. J Neurosci 21: 5229–5238, 2001.[Abstract/Free Full Text]

Golowasch J, Abbott LF, and Marder E. Activity-dependent regulation of potassium currents in an identified neuron of the stomatogastric ganglion of the crab Cancer borealis. J Neurosci 19: RC33(1–5), 1999.[Abstract/Free Full Text]

Golowasch J, Buchholtz F, Epstein IR, and Marder E. Contribution of individual ionic currents to activity of a model stomatogastric ganglion neuron. J Neurophysiol 67: 341–349, 1992.[Abstract/Free Full Text]

Golowasch J, Goldman MS, Abbott LF, and Marder E. Failure of averaging in the construction of a conductance-based neuron model. J Neurophysiol 87: 1129–1131, 2002.[Abstract/Free Full Text]

Harris-Warrick RM. Voltage-sensitive ion channels in rhythmic motor systems. Curr Opin Neurobiol 12: 646–651, 2002.[CrossRef][Web of Science][Medline]

Hill AAV, Lu J, Masino MA, and Calabrese RL. A model of a segmental oscillator in the leech heartbeat neuronal network. J Comput Neurosci 10: 281–302, 2001.[CrossRef][Web of Science][Medline]

Hobbs KH and Hooper SL. Using a genetic algorithm to tune a lobster (panulirus) model pyloric neuron. Soc Neurosci Abstr 605.2, 2003.

Hoffman DA, Magee JC, Colbert CM, and Johnston D. K⁺ channel regulation of signal propagation in dendrites of hippocampal pyramidal neurons. Nature 387: 869–875, 1997.[CrossRef][Medline]

Hooper SL. Multiple routes to similar network output. Nat Neurosci 7: 1287–1288, 2004.[CrossRef][Web of Science][Medline]

Hooper SL and Marder E. Modulation of the lobster pyloric rhythm by the peptide proctolin. J Neurosci 7: 2097–2112, 1987.[Abstract]

Houck CR, Joines JA, and Kay MG. A Genetic Algorithm for Function Optimization: A Matlab Implementation. Raleigh, NC: Industrial Engineering Department, North Carolina State Univ., 1995, p. 1–14.

Ivanov AI and Calabrese RL. Spatial distribution and some properties of the synaptic connections between leech heart (HN) interneurons. Soc Neurosci Abstr 24: 751.715, 1998.

Jezzini SH, Hill AA, Kuzyk P, and Calabrese RL. Detailed model of intersegmental coordination in the timing network of the leech heartbeat central pattern generator. J Neurophysiol 91: 958–977, 2004.[Abstract/Free Full Text]

Keren N, Peled N, and Korngreen A. Constraining compartmental models using multiple voltage recordings and genetic algorithms. J Neurophysiol 94: 3730–3742, 2005.[Abstract/Free Full Text]

Korngreen A and Sakmann B. Voltage-gated K⁺ channels in layer 5 neocortical pyramidal neurones from young rats: subtypes and gradients. J Physiol 525: 621–639, 2000.[Abstract/Free Full Text]

LeMasson G and Maex R. Introduction to equation solving and parameter fitting. In: Computational Neuroscience: Realistic Modeling for Experimentalists, edited by De Schutter E. Boca Raton, FL: CRC Press, 2001, p. 1–24.

LeMasson G, Marder E, and Abbott LF. Activity-dependent regulation of conductances in model neurons. Science 259: 1915–1917, 1993.[Abstract/Free Full Text]

Magee JC. Dendritic Ih normalizes temporal summation in hippocampal CA1 neurons. Nat Neurosci 2: 508–514, 1999.[CrossRef][Web of Science][Medline]

Marder E, Bucher D, Schulz DJ, and Taylor AL. Invertebrate central pattern generation moves along. Curr Biol 15: R685–R699, 2005.[CrossRef][Web of Science][Medline]

Marder E and Calabrese RL. Principles of rhythmic motor pattern generation. Physiol Rev 76: 687–717, 1996.[Abstract/Free Full Text]

Marder E, Hooper SL, and Siwicki KK. Modulatory action and distribution of the neuropeptide proctolin in the crustacean stomatogastric nervous system. J Comp Neurol 243: 454–467, 1986.[CrossRef][Web of Science][Medline]

Masino MA and Calabrese RL. Phase relationships between segmentally organized oscillators in the leech heartbeat pattern generating network. J Neurophysiol 87: 1572–1585, 2002a.[Abstract/Free Full Text]

Masino MA and Calabrese RL. Period differences between segmental oscillators produce intersegmental phase differences in the leech heartbeat timing network. J Neurophysiol 87: 1603–1615, 2002b.[Abstract/Free Full Text]

Nadim F and Calabrese RL. A slow outward current activated by FMRFamide in heart interneurons of the medicinal leech. J Neurosci 17: 4461–4472, 1997.[Abstract/Free Full Text]

Nadim F and Manor Y. The role of short-term synaptic dynamics in motor control. Curr Opin Neurobiol 10: 683–690, 2000.[CrossRef][Web of Science][Medline]

Nadim F, Manor Y, Nusbaum MP, and Marder E. Frequency regulation of a slow rhythm by a fast periodic input. J Neurosci 18: 5053–5067, 1998.[Abstract/Free Full Text]

Nadim F, Olsen OH, De Schutter E, and Calabrese RL. Modeling the leech heartbeat elemental oscillator. I. Interactions of intrinsic and synaptic currents. J Comput Neurosci 2: 215–235, 1995.[CrossRef][Web of Science][Medline]

Norris BJ, Weaver AL, Morris LG, Wenning A, Garcia PA, and Calabrese RL. A central pattern generator producing alternative outputs: temporal pattern of premotor activity. J Neurophysiol 96: 309–326, 2006.[Abstract/Free Full Text]

Peterson EL. Generation and coordination of heartbeat timing oscillation in the medicinal leech. I. Oscillation in isolated ganglia. J Neurophysiol 49: 611–626, 1983.[Abstract/Free Full Text]

Pinsky PF and Rinzel J. Intrinsic and network rhythmogenesis in a reduced Traub model for CA3 neurons. J Comput Neurosci 1: 39–60, 1994.[CrossRef][Medline]

Prinz AA, Billimoria CP, and Marder E. Alternative to hand-tuning conductance-based models: construction and analysis of databases of model neurons. J Neurophysiol 90: 3998–4015, 2003.[Abstract/Free Full Text]

Prinz AA, Bucher D, and Marder E. Similar network activity from disparate circuit parameters. Nat Neurosci 7: 1345–1352, 2004.[CrossRef][Web of Science][Medline]

Schmidt J and Calabrese RL. Evidence that acetylcholine is an inhibitory transmitter of heart interneurons in the leech. J Exp Biol 171: 329–347, 1992.[Abstract/Free Full Text]

Schulz DJ, Goaillard JM, and Marder E. Variable channel expression in identified single and electrically coupled neurons in different animals. Nat Neurosci 9: 356–362, 2006.[CrossRef][Web of Science][Medline]

Skinner FK, Turrigiano GG, and Marder E. Frequency and burst duration in oscillating neurons and two-cell networks. Biol Cybern 69: 375–383, 1993.[Web of Science][Medline]

Skou JC. Further investigations on a (Mg²⁺+Na⁺)-activated adenosine-triphosphatase possible related to the active, linked transport of Na⁺ and K⁺ across the nerve membrane. Biochim Biophys Acta 42: 6–23, 1960.[CrossRef]

Skou JC. Enzymatic basis for active transport of Na⁺ and K⁺ across cell membrane. Physiol Rev 45: 596–617, 1965.[Free Full Text]

Sorensen M, DeWeerth S, Cymbalyuk G, and Calabrese RL. Using a hybrid neural system to reveal regulation of neuronal network activity by an intrinsic current. J Neurosci 24: 5427–5438, 2004.[Abstract/Free Full Text]

Stuart G and Hausser M. Initiation and spread of sodium action potentials in cerebellar Purkinje cells. Neuron 13: 703–712, 1994.[CrossRef][Web of Science][Medline]

Thoby-Brisson M, Telgkamp P, and Ramirez JM. The role of the hyperpolarization-activated current in modulating rhythmic activity in the isolated respiratory network of mice. J Neurosci 20: 2994–3005, 2000.[Abstract/Free Full Text]

Thompson WJ and Stent GS. Neuronal control of heartbeat in the medicinal leech. I. Generation of the vascular constriction rhythm by heart motor neurons. J Comp Physiol 111: 261–279, 1976a.[CrossRef]

Thompson WJ and Stent GS. Neuronal control of heartbeat in the medicinal leech. II. Intersegmental coordination of heart motor neuron activity by heart interneurons. J Comp Physiol 111: 281–307, 1976b.[CrossRef]

Thompson WJ and Stent GS. Neuronal control of heartbeat in the medicinal leech. III. Synaptic relations of the heart interneurons. J Comp Physiol 111: 309–333, 1976c.[CrossRef]

Tobin A-E and Calabrese RL. Myomodulin increases I_h and inhibits the Na/K pump to modulate bursting in leech heart interneurons. J Neurophysiol 94: 3938–3950, 2005.[Abstract/Free Full Text]

Tobin A-E, Van Hooser SD, and Calabrese RL. Creation and reduction of a morphologically detailed model of a leech heart interneuron. J Neurophysiol 96: 2107–2120, 2006.[Abstract/Free Full Text]

Traub RD, Wong RK, Miles R, and Michelson H. A model of a CA3 hippocampal pyramidal neuron incorporating voltage-clamp data on intrinsic conductances. J Neurophysiol 66: 635–650, 1991.[Abstract/Free Full Text]

Turrigiano G, Abbott LF, and Marder E. Activity-dependent changes in the intrinsic properties of cultured neurons. Science 264: 974–977, 1994.[Abstract/Free Full Text]

Vanier MC and Bower JM. A comparative survey of automated parameter-search methods for compartmental neural models. J Comput Neurosci 7: 149–171, 1999.[CrossRef][Web of Science][Medline]

Wenning A, Cymbalyuk GS, and Calabrese RL. Heartbeat control in leeches. I. Constriction pattern and neural modulation of blood pressure in intact animals. J Neurophysiol 91: 382–396, 2004.[Abstract/Free Full Text]

This article has been cited by other articles:

				R. Grashow, T. Brookings, and E. Marder Reliable neuromodulation from circuits with variable underlying structure PNAS, July 14, 2009; 106(28): 11742 - 11746. [Abstract] [Full Text] [PDF]

				D. J. Schulz, J.-M. Goaillard, and E. E. Marder Quantitative expression profiling of identified neurons reveals cell-specific constraints on highly variable levels of gene expression PNAS, August 7, 2007; 104(32): 13187 - 13191. [Abstract] [Full Text] [PDF]

This Article Abstract Full Text (PDF) All Versions of this Article: 96/4/2089    most recent 00025.2006v1 Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in Web of Science Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Web of Science (4) Citing Articles via Google Scholar Google Scholar Articles by Tobin, A.-E. Articles by Calabrese, R. L. Search for Related Content PubMed PubMed Citation Articles by Tobin, A.-E. Articles by Calabrese, R. L.