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


     


J Neurophysiol 92: 96-110, 2004. First published March 17, 2004; doi:10.1152/jn.01146.2003
0022-3077/04 $5.00
This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow All Versions of this Article:
92/1/96    most recent
01146.2003v1
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 (6)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Abarbanel, H. D. I.
Right arrow Articles by Talathi, S.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Abarbanel, H. D. I.
Right arrow Articles by Talathi, S.

Mapping Neural Architectures Onto Acoustic Features of Birdsong

Henry D. I. Abarbanel1,2,3, Leif Gibb3, Gabriel B. Mindlin3,4 and Sachin Talathi1,3

1Department of Physics, 2Marine Physical Laboratory (Scripps Institution of Oceanography), and 3Institute for Nonlinear Science, University of California, San Diego, La Jolla, California 92093-0402; and 4Departamento de Física, Facultad de Ciencias Exactas y Naturales, Universidad de Buenos Aires, Ciudad Universitaria, Pabellon I (1428), Buenos Aires, Argentina

Submitted 1 December 2003; accepted in final form 7 March 2004


    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
The motor pathway responsible for the complex vocalizations of songbirds has been extensively characterized, both in terms of intrinsic and synaptic physiology in vitro and in terms of the spatiotemporal patterns of neural activity in vivo. However, the relationship between the neural architecture of the song motor pathway and the acoustic features of birdsong is not well understood. Using a computational model of the song motor pathway and the songbird vocal organ, we investigate the relationship between song production and the neural connectivity of nucleus HVc (used as a proper name) and the robust nucleus of the archistriatum (RA). Drawing on recent experimental observations, our neural model contains a population of sequentially bursting HVc neurons driving the activity of a population of RA neurons. An important focus of our investigations is the contribution of intrinsic circuitry within RA to the acoustic output of the model. We find that the inclusion of inhibitory interneurons in the model can substantially influence the features of song syllables, and we illustrate the potential for subharmonic behavior in RA in response to forcing by HVc neurons. Our results demonstrate the association of specific acoustic features with specific neural connectivities and support the view that intrinsic circuitry within RA may play a critical role in generating the features of birdsong.


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
Birdsong is a rich example of behavior emerging from the interplay between neural circuits, motor functions, and environment. In the process of creating song, neural instructions drive a highly nonlinear device—the syrinx—capable, along with respiratory control, of generating simple whistles as well as rich, complex sounds (Fee et al. 1998Go).

Song requires the coordination of vocal and respiratory muscles. The neural circuitry generating the patterned activities of these muscles includes well-defined sets of interconnected neurons, called nuclei, constituting a song motor pathway (Nottebohm et al. 1976Go). The properties of these nuclei have been investigated in some detail in the last few years, and there is now an emerging descriptive model of this and other song pathways. This article explores a simple mathematical model of the song motor pathway, from nucleus HVc to the robust nucleus of the archistriatum (RA) to muscles, which allows the quantitative exploration of properties of the descriptive model.

Interest in songbirds has been stimulated by the similarities between vocal learning in songbirds and humans (Doupe and Kuhl 1999Go). The neural processes involved in the production, learning, and control of birdsong have been carefully investigated because they constitute a paradigm for the learning of complex behavior and may have lessons for similar activity in humans (Perkel et al. 2002Go).

The purpose of the present work is to study the way in which the architectural properties of the neural circuits of the song motor pathway map onto acoustic features of the song. After reviewing the present knowledge on the structure and connectivity of the main nuclei in this pathway, we build a mathematical model of its principal neural circuits. In this model, each unit in the circuit displays a dynamics ruled by the Hodgkin–Huxley equations.

A previous effort (Laje and Mindlin 2002Go) attempted to build this link between neural circuits and acoustic features. In this previous work, the activity in the nuclei of the song motor pathway was simulated in terms of rate models. However, these simple models fail to capture recently described aspects of the dynamics found in the motor pathway, in particular that the activity of RA neurons seems to be strongly determined by sparsely spiking RA-projecting HVc neurons (Hahnloser et al. 2002Go). For this reason, in this work we explore the dynamics displayed by small arrays of spiking excitable units compatible with these constraints.

We use the activities of different sets of units of our model of RA as a forcing function for a mathematical model of the syrinx, which is capable of generating synthetic sounds (Gardner et al. 2001Go; Laje et al. 2002Go). We also explore the different acoustic features produced when the connectivity of our model song motor pathway is changed. In this way we are able to close the link between neurons and song.

The model for the syrinx not only permits the generation of synthetic songs but also suggests that many different syllables can be produced by controlling the relative timing between the oscillations that take place in lung pressure and syringeal muscle tension during the vocalization of syllables, as hypothesized (Gardner et al. 2001Go). According to this hypothesis, the mapping from central neural activity onto acoustic features is equivalent to an appropriate mapping from central neural activity onto the region of the song motor pathway where these timing differences originate. We focus on the analysis of the RA nucleus (Spiro et al. 1999Go) and explore how its internal architecture translates into acoustic features of the song. More specifically, we test the hypothesis that the intrinsic circuitry within RA, both excitatory and inhibitory, is important for the generation of appropriate sound output. The nonlinear nature of the RA model reveals itself in its dynamic behavior under particular conditions. This allows us to elaborate some conjectures concerning the structure of complex syllables.

The song control nuclei

The brain structure that controls the production of song in songbirds consists of discrete sets of neurons (nuclei) and the axons that interconnect them, forming a pathway (Catchpole and Slater 1995Go). This pathway runs from nucleus HVc to RA. In RA there are neurons projecting to the motor neurons that innervate the muscles of the syrinx as well as neurons projecting to areas that control respiration (Spiro et al. 1999Go; Vicario 1993Go; Wild 1993Go). Repetitive syllables, characteristic of many songs, are created when these projection populations of RA display patterns of activity. Syllables combined with quiet periods between them compose the motifs that are combined in sequence to produce a song.

It has been suggested that the patterns of activity in RA neurons are the result of interaction between local circuits and input instructions from HVc (Yu and Margoliash 1996Go). The intrinsic physiology and morphology of RA interneurons and projection neurons, and the synaptic connections between these 2 classes of neuron within RA, have recently been described in some detail (Spiro et al. 1999Go). In particular, it was reported that interneurons make long-range inhibitory connections onto projection neurons within RA. These data provide a picture of how the intrinsic circuitry of RA is organized.

RA receives excitatory glutamatergic input from HVc (Kubota and Saito 1991GoGo; Mooney 1992Go; Mooney and Konishi 1991Go; Spiro et al. 1999Go). Just as RA is closely tied to the peripheral motor system, HVc is believed to be closely tied to higher-level functions (Suthers and Margoliash 2002Go; Yu and Margoliash 1996Go). A song is organized in terms of minimal continuous utterances that we call syllables, combined with quiet periods to build motifs. It has been proposed that this hierarchy in the song has a neural correlate based on the observation that electrical stimulation of HVc causes a resetting of the motif being sung, whereas stimulation of RA causes a distortion of the vocalization (Vu et al. 1994Go). This suggests that the HVc activity is sequenced to orchestrate the production of syllables and motifs, so interrupting its progress along such a sequence causes a resetting or reinitiation of the activity. It also suggests that RA is a "junction box," where HVc signals are combined to produce song, but does not by itself initiate or command the song; song production can be influenced in RA but not started there. The connections in the "junction box" (RA) can be altered by modification of its input connections from HVc and the anterior forebrain pathway.

Work on zebra finches (Taeniopygia guttata) has analyzed in detail the time relation between the firing of some neurons in HVc and RA during song (Hahnloser et al. 2002Go). It shows that during the vocalization of a motif, lasting approximately 1,000 ms, a given RA-projecting neuron in HVc is active only during one window of temporal size 6.1 ± 2 ms and produces 4.5 ± 2 spikes in that temporal window. If a motif is repeated, any given HVc neuron will repeat its spiking in the same time window.

On the other hand, during singing RA neurons generate highly stereotyped sequences of action potential bursts, typically a few bursts of about 10-ms duration per motif, well time locked with parts of syllables (Chi and Margoliash 2001). The picture emerging from these experiments is that during each time window in the RA sequence, RA neurons are driven by a subpopulation of RA-projecting HVc neurons that are active only during that window of time (Hahnloser et al. 2002Go).

These experiments suggest that the premotor burst patterns in RA are driven by the activities of HVc neurons (McCasland et al. 1981Go). In that case, the architecture of the connectivity between the HVc nucleus and the RA nucleus will play an important role in determining the complex patterns of activity in RA neurons, which are responsible for driving the avian vocal organ. The different populations of RA excitatory projection neurons that control the respiratory rhythms and the muscular activity in the syrinx are connected through long-range inhibitory interneurons local to RA (Spiro et al. 1999Go).

The data of Herrmann and Arnold (1991)Go suggest that the majority of synapses onto RA neurons are intrinsic connections from other neurons within RA. Some of these intrinsic connections are from inhibitory interneurons (Spiro et al. 1999Go), whereas others are excitatory synapses made by local axon collaterals of RA projection neurons (Perkel 1995Go; Spiro et al. 1999Go). To unveil the relative weight of these mechanisms, we explore the effect of different architectures on the acoustic properties of synthetic song generated by a model of the syrinx driven by our model nuclei.

We now turn to a simplified model of the HVc and RA nuclei that accounts for the essential aspects of HVc command to the RA neurons. The synaptic connections from HVc to RA are one to many, so that the firing of an HVc neuron influences the firing of several RA neurons in the same time window. The RA projection neurons are split into 2 populations, one of which influences the muscles of the syrinx and is expressed as a change in the effective time-dependent "spring constant" T(t) of the syrinx apparatus. The other population of RA projection neurons influences respiratory dynamics, expressed through the time-dependent driving pressure P(t) from the respiratory system. In addition to excitatory neurons in both HVc and RA, we include in a simplified manner local interneurons, which provide inhibition. The role of the inhibition is explored within a class of models, and it appears to be important in producing the timing of T(t) relative to P(t) during song. Beyond inducing changes at the level of an individual syllable, inhibition plays an important role in the generation of complex sequences of syllables when RA is periodically driven by HVc in our model. Finally, we explore the dynamic and acoustic consequences of varying the strength of the local excitatory connections between the projection neurons in RA.


    METHODS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
We model the RA and HVc nuclei as sets of Hodgkin–Huxley (HH) neurons with synaptic interconnections. We initiate the sequence of bursting of the RA-projecting, excitatory HVc neurons through the induction of spiking activity in one member of the HVc population. Then through feedforward connections, the remainder of the HVc population bursts sequentially. We do not inquire here how the sequential firing is initiated by command signals from other higher centers in the nervous system, nor do we attempt to produce a comprehensive model of HVc dynamics. The inhibitory connections in HVc (Mooney 2000Go) are represented in a simplified manner by an HH interneuron receiving excitatory connections from, and sending inhibitory connections to, the excitatory HVc population. The RA projection neurons are also represented as HH neurons, which are driven by excitatory inputs from HVc, have excitatory interconnections among themselves, and receive inhibition by an HH interneuron.

Modeling the HVc nucleus

The excitatory neurons of HVc project either to Area X (Perkel et al. 2002Go) or to RA, but not to both. We consider only the HVc -> RA connections here. Our analysis of Area X and other nuclei of the anterior forebrain pathway, and of how that set of nuclei influences RA operation, will be presented in a forthcoming report. The different classes of neurons in HVc have been characterized physiologically, pharmacologically, and morphologically in the zebra finch (Dutar et al. 1998Go; Kubota and Saito 1991aGo; Kubota and Taniguchi 1998Go). These elements will eventually allow us to understand the activity observed in HVc in terms of the emergence of complex neural dynamics. At this moment, we are interested only in generating the simple activity of the projection neurons that drive RA during song production. To emulate such activity, we formulate a most simple model.

We simulate a set of NHVc = 10 HVc -> RA projection neurons, denoting the membrane potential for the jth HVc excitatory neuron as VHj(t) with j = 1, 2, ... , NHVc. It satisfies the HH equation

(1)
where we have included a leak, sodium, and potassium current in the standard fashion. gL is the conductance for the leak current, and gNa and gK are the maximal conductances for the voltage-dependent currents. At IDC = 0 each isolated HH neuron is silent and has a resting potential of {approx} –65 mV. The activation and inactivation variables mj(t), hj(t), and nj(t) satisfy the usual first-order kinetics

(2)
where Xj(t) stands for mj(t), hj(t), or nj(t). The functions appearing here are given in an APPENDIX.

The function S[Vpre(t)] represents the fraction of postsynaptic receptor channels that are open in response to the binding of excitatory or inhibitory neurotransmitter. It is taken to satisfy the first-order kinetics

(3)
where {alpha} and {beta} are rate constants. We use the values {alpha} = 0.15 ms–1, {beta} = 0.2275 ms–1, and Vp = 10 mV in our calculations. This equation for S[Vpre(t)] has the property that when the presynaptic voltage rises, as during an action potential, S[Vpre(t)] rises toward the value {alpha}/({alpha} + {beta}) with a time constant 1/({alpha} + {beta}). When the presynaptic voltage Vpre(t) is large and negative (i.e., near the resting potential), S[Vpre(t)] -> 0 on a time scale 1/{beta}.

The membrane voltage VHI(t) of the inhibitory neuron satisfies an HH equation

(4)
At IDC = 0, the inhibitory HVc neuron is also at a resting potential of about –65 mV. The last term in this equation represents the input to the inhibitory unit from all the excitatory HVc units in the model.

These equations describe a set of HVc neurons that are connected in a feedforward fashion so that activation of neuron j = 1 excites j = 2 and so forth to j = NHVc, with S(VH0) = 0. We use constants guaranteeing that an initial excitation of the first unit will excite the remaining 9 units sequentially. The parameters of the system are adjusted so that the activity of an RA-projecting HVc unit consists of 4 to 7 spikes (Hahnloser et al. 2002Go).

When GHVcI != 0 and GHI != 0, the inhibitory neuron receives excitatory input in a global fashion from the excitatory HVc neurons and each HVc neuron receives inhibition as well as sequential excitation. In our calculations we use GHE = 17.7 mS/cm2, GHI = 7.5 mS/cm2, GHVcI = 10.5 mS/cm2, ERevE = 0 mV, and ERevI = –80 mV. The DC of the interneuron was set as IDC,I 0.

The number of units (10) used to emulate the activity of RA-projecting HVc neurons allows us to represent sequential units sparsely spiking during the nearly 40 ms that a brief syllable lasts in our models. If we wish to generate more than one syllable, longer syllables, or bursting patterns in HVc with several neurons active simultaneously, the model will have to be extended. In our simulations, the most important parameter when adjusting the HVc network to generate sequential patterns of activity with the appropriate number of spikes was GHE. We obtained acceptable sequences for GHE in the interval [15 mS/cm2, 18 mS/cm2].

Modeling the RA nucleus

A set of NRA = 10 units is used to represent the projection neurons in RA. Five of these units are taken to be responsible for the control of the syringeal muscle that determines the frequency of the vocalizations (Goller and Suthers 1996Go; Laje and Mindlin 2002Go), and the other 5 units are assumed to force the expiratory rhythm. Each of these 10 RA units is innervated by a subset of the HVc units. The number of units used to represent RA projection neurons in our study (10) is the smallest such that the activity in each subpopulation ("syringeal" and "respiratory") is capable of generating motor instructions leading to noticeable and continuous changes in the acoustic features of the synthetic syllables generated by our syringeal model. In this work we focus on understanding the basic connections between neural architecture and song. Further work with larger models will allow us to reproduce more subtle acoustic features.

Developmental changes in the density and number of synaptic connections between HVc and RA neurons, and other changes in RA during song learning, have been investigated (Herrmann and Arnold 1991Go; Kittelberger and Mooney 1999Go). Changes in the synaptic connectivity within RA can also be induced by juvenile lesions of the anterior forebrain nucleus, the lateral magnocellular nucleus of the anterior neostriatum (LMAN; Kittelberger and Mooney 1999Go). The relative contributions of AMPA and N-methyl-D-aspartate (NMDA) receptors at these synapses have been described (Kubota and Saito 1991bGo; Mooney 1992Go; Mooney and Konishi 1991Go; Stark and Perkel 1999Go). Both AMPA and NMDA receptors contribute to synapticcurrents at HVc–RA synapses, although the contribution of NMDA-mediated synaptic currents decreases as maturation takes place (Stark and Perkel 1999Go). Here we use a simplified model of RA connectivity in the adult.

The ath projection neuron in RA satisfies an HH equation with the addition of synaptic currents from HVc projection neurons, inhibitory interneurons, and other projection neurons in RA. The membrane voltage of the ath RA neuron, VRA(t), with a = 1, 2, , NRA, satisfies

(5)
where, again, we have the usual leak, sodium, and potassium currents, as well as a DC current. The matrix {Gamma}aj dictates the strength of the synaptic connections between HVc and RA neurons.

The matrix {Gamma}ab determines the strength of the excitatory synaptic connections between projection neurons within RA. Based on the observation that axon collaterals of RA projection neurons extend <200 µm from the soma (Spiro et al. 1999Go), we consider only excitatory connections between close neighbors (within the same subpopulation of five units that controls either tension or respiratory activity). The matrix elements corresponding to other excitatory connections between RA units are zero.

In addition there is a global inhibitory coupling among the RA neurons. The membrane voltage of the inhibitory neuron in RA satisfies

(6)
These equations thus describe a set of RA excitatory neurons that are driven by connections from HVc neurons, excite each other locally, and are inhibited by the inhibitory RA neuron with membrane potential VRI(t). The inhibitory RA neuron, in turn, is driven by synaptic input from all the RA excitatory neurons, as well as by projection neurons in HVc. In our calculations reported below, we always take {Gamma}aj = {Gamma}Aaj with {Gamma} = 18.62 mS/cm2 and Aaj a matrix with entries of 0 or 1. Similarly, {Gamma}ab = {Gamma}'Bab, with {Gamma}' = 3.5 mS/cm2 and Bab a matrix with entries 0 or 1. Within each of the 2 populations of excitatory units in RA, the local nature of their interconnections is emulated by choosing Bab = 0 whenever ||ab|| > 1. Unless stated otherwise, the value {Gamma}' = 3.5 mS/cm2 is used. This value guarantees that the activity of an HVc-driven RA unit induces spikes only in its first neighbors. Finally, GHI = 5 mS/cm2 is used in our simulations to represent the forcing of the inhibitory neuron in RA by projection neurons in HVc. The value of IDC,R = 1.93 µA/cm2 is chosen so that the units representing excitatory neurons in RA display, before HVc activity starts, spontaneous spiking at approximately 20 Hz. In this work, we explore the different dynamics that emerge in RA as a function of the long-range inhibition. Numerical simulations explore the dynamics presented by our model as GRIa takes values in the range 0–75 mS/cm2.

In the calculations below we have 10 excitatory units representing RA projection neurons and one unit representing inhibitory interneurons. In Fig. 1, we display the basic elements used to build our model. To make the figure readable, we choose to illustrate one connection as representative of each of the types included in our model. The connection from the inhibitory unit in RA (dark gray) to excitatory unit 1 in RA represents long-range inhibition (in the mathematical implementation of the model, all excitatory units receive inhibition from this unit). The excitatory connections from RA unit 7 to RA units 6 and 8 represent the local interconnections between excitatory units (in the mathematical implementation of our model, every excitatory unit synapses on one or 2 other excitatory units within its subpopulation of 5 units). The excitatory connection from unit 5 in RA to the inhibitory neuron in RA is representative of the connections from each excitatory unit to the inhibitory unit in RA in the mathematical implementation of the model.



View larger version (21K):
[in this window]
[in a new window]
 
FIG. 1. Hodgkin–Huxley (HH) neuron units representing neurons in a songbird's song motor pathway. Ten units in HVc generate a sparse sequence of action potentials initiated by firing unit 1; these units are linked sequentially by feedforward synaptic connections (a representative connection is shown between units 9 and 10). These HVc units excite robust nucleus of archistriatum (RA) units by another set of synaptic connections. In this example, HVc unit 1 connects to RA units 6 and 7. Five units in RA drive the motor units controlling the tension T(t) of the labia in the syrinx, and the other 5 units control the bronchial pressure P(t). In both the RA and the HVc nuclei we represent the global inhibition by a single inhibitory neuron, which receives excitation from the other excitatory neurons in the nucleus and in turn inhibits them. A representative inhibitory connection in RA is shown between the inhibitory unit and unit 1. Excitation of the inhibitory unit by the excitatory unit is illustrated in the figure by the connection between unit 5 and the inhibitory unit. PN, projection neuron; I, interneuron.

 
HVc -> RA neural connections

According to the experiments on the time relation between firing of neurons in HVc and RA (Hahnloser et al. 2002Go), a sequence of RA-projecting HVc neurons fire briefly during a motif, each of them "overseeing" (Nottebohm 2002Go) the activity of a set of RA neurons. Therefore our first concern will be to investigate the way in which different connectivities between HVc and RA translate into acoustic features of song in this pure feedforward scheme.

The matrix A defines the connectivity between the 2 nuclei: the values of Aa,j can be either 0 or 1. If Aa,j = 1, there is a connection between the jth neuron in HVc and the ath neuron in RA; otherwise, these 2 neurons are not connected. The units representing neurons in HVc are ordered according to the time at which they produce a burst of spikes: the jth unit is the jth to spike in the temporal sequence described by Hahnloser et al. (2002)Go. The units representing RA neurons, on the other hand, are not ordered in any particular way in our model. Each HVc neuron fires once in a motif while inducing more than one RA neuron to fire.

In Fig. 2 we display an example of HVc-driven RA dynamics. In Fig. 2A, we show the membrane potential of a pair of HVc excitatory units (number 1 and number 10), and in Fig. 2B we show the response of an RA excitatory unit (number 1), for a connectivity matrix A1,1 = A1,10 = 1 (and Ai,j = 0, otherwise). Note that before the activity starts in HVc, the unit in RA is spiking at a low frequency. As the sequence of bursts is generated in HVc, the RA inhibitory unit is driven into a spiking regime (not shown), which suppresses the spontaneous activity of the RA excitatory unit. The RA excitatory unit spikes when HVc units 1 and 10 in HVc drive it beyond threshold. In this simulation GRIa = 0.



View larger version (23K):
[in this window]
[in a new window]
 
FIG. 2. Two bursts of activity in 2 different HVc units, unit 1 and unit 10 (A), which drive a single unit in RA, unit 1 (B). A1,1 = A1,10 = 1 (and Ai,j = 0, otherwise). Note that both before the bursting in HVc starts and a few milliseconds after it ends, the units in RA are spiking in a low-frequency dynamic regime.

 
Neural forcing of the syrinx model

A subpopulation of 5 of our 10 RA neurons is taken to project onto motor neurons in the tracheosyringeal portion of the hypoglossal nucleus (nXIIts), which drives the syringeal muscles, and another subpopulation of 5 is assumed to project onto neurons influencing the respiratory rhythm. The dynamics of the syringeal muscles is given by a tension T(t), whereas the respiratory dynamics is represented by a pressure P(t). We assume that all the motor neurons driving the syringeal muscles make synaptic contact with a similar number of muscle fibers. The tension of the muscles can be increased by increasing the number of firing motor neurons (recruitment). We further assume that to recruit a large number of motor neurons at a given time, a large number of RA neurons must be firing.

Under these hypotheses, at any given time the sum of the activities within each one of the 2 RA subpopulations (the one projecting onto motor neurons controlling the syringeal muscles, and the other locking the expiratory rhythm; Gardner et al. 2001Go) will constitute the driving input to the corresponding muscles. Let us denote by VaT,RA(t) the voltage of the ath RA unit representing a neuron in RA that controls syringeal tension, and by VaP,RA the voltage of the ath unit representing a neuron in RA that controls bronchial pressure. Then, the instructions controlling the syrinx model will be


(7)
with {tau}RA = 200 ms and Pos(x) = x if x > 0; Pos(x) = 0, otherwise. This means that when the syringeal or respiratory population is excited above a threshold {theta} (taken as {theta} = –200 in our simulations), T(t) or P(t) changes from its rest value of 0, and when the excitation is completed, relaxation proceeds to the rest state on a time scale of {tau}RA.

T(t) and P(t) are used to drive a model of the syrinx. In previous work, the same amount of filtering was applied to electromyographic data to generate synthetic signals that were later used to drive a mathematical model of the syrinx. This procedure led to the generation of synthetic sound that was acoustically very similar to the actual recordings of song performed while the electromyographic data were being taken (Mindlin et al. 2003Go).

Syrinx model

In the syrinx model the labia are assumed to support both lateral oscillations and an upward propagating surface wave (Titze 1988Go). In this way, the opposing labia have a convergent profile when they move away from each other and a divergent profile when they move toward each other. This results in a higher pressure on the labia during the opening phase, and an overall gain in energy during each cycle of oscillation (Gardner et al. 2001Go). The reason is that, while presenting a convergent profile, the average pressure between the labia is closer to the bronchial pressure, whereas interlabial pressure is closer to atmospheric pressure for a divergent profile. This results in a force in the same direction as the velocity of displacement of the labia, which might overcome the dissipation. The labia are also subjected to other forces: restitution attributed to the elasticity of the tissues and eventually a high dissipation, if the labia displace too much from their equilibrium position. Labial oscillation can be interrupted if an additional force pushes the labia against each other (or against the walls). Adding up all these effects, a simple set of equations can be written to describe the motion of the variable describing the departure of the midpoint of a labium from the prephonatory position, x(t) (Gardner et al. 2001Go; Laje et al. 2002Go)


(8)
where {alpha}1T(t) + {alpha}0 stands for the restitution constant of the labium, {beta}1P(t) + {beta}0 is the negative dissipation induced by the interlabial pressure as discussed above, and b is the linear dissipation. The parameters ({alpha}1, {alpha}0) and ({beta}1, {beta}0) are linear coefficients that convert the dimensionless activities T(t) and P(t) into a restitution constant and a negative dissipation, respectively, and are adjusted to produce vocalizations in a realistic frequency range. In this work we use {alpha}1 = 1.1 x 105, {alpha}0 = 0.9 x 108, {beta}1 = 8.75, and {beta}0b = –15 x 10–3. The parameter C is the nonlinear dissipation constant. For large departures from the equilibrium position [i.e., large absolute values of x(t)], the nonlinear dissipation term becomes very large. This models the effect of containing walls, as described in Laje et al. (2002)Go. Throughout this work, we take C = 2 x 108. The instructions drive the syrinx model to generate a syllable using an integration step dt = 1/44,100 s. The choice of parameters of our mathematical model of the syrinx is such that the typical frequency of oscillation during the vocalization is around 2 kHz.

Once the dynamics of the labia is given, it is possible to generate synthetic songs (Laje et al. 2002Go). The cyclic activities of P(t) and T(t) periodically drive the system within the oscillating region, producing syllables. The dynamic evolution of the fundamental frequency within a syllable is one of the most noticeable acoustic features of birdsong. This evolution can be exhibited by means of spectrograms.

The model used to generate synthetic song was originally derived to emulate the song of canaries (Gardner et al. 2001Go). Songs of other species were fit as well: chingolo sparrows (Laje et al. 2002Go), white-crowned sparrows (Laje and Mindlin 2002Go), and cardinals (Mindlin et al. 2003Go). In this last case, the integration of our model, driven by experimentally recorded parameters P(t) and T(t) from singing birds, gave rise to synthetic songs that were very similar to the recorded ones. It is interesting to point out that sparrows producing the kind of tonal sounds generated by our model produce up and down "slurs" as well as whistles. Yet, it is likely that for zebra finches, many of the uttered sounds could require additional modeling efforts. We are assuming throughout this work that the neurophysiological results used to build our model, obtained mostly for zebra finches, will be representative of the neural behavior in a wide set of songbird species.

Numerical methods

To integrate the equations describing the dynamics of the RA and HVc nuclei we used a Runge–Kutta 6(5) scheme with variable time step and a relative error of 10–10. The synthetic instructions generated by this algorithm were properly sampled to drive our equations for syringeal motion, which were integrated with a Runge–Kutta 4 scheme. The synthetic songs were analyzed with the software snd, running on a Linux platform. Epochs of 1,000 ms were integrated before the onset of syllable production to incorporate into our description the spontaneous spiking of RA projection neurons.


    RESULTS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
In the framework of a completely feedforward model in which the populations of RA projection units are not coupled through inhibition, it is easy to understand how an architecture of connectivity maps onto motor instructions. We will compare 2 cases. In the first case, the architecture of connectivity between HVc and RA will be such that the first RA-projecting neurons in HVc to burst during the production of a syllable recruit fewer syringeal-projecting units in RA than the last ones. Mathematically, this can be achieved in our description by matrices Aij such that {sum}i=15 Aij monotonically increases with j. In the second case, the first RA-projecting HVc neurons active within a syllable recruit more syringeal units in RA than the last ones. In this case, the matrices Aij satisfy that {sum}i=15 Aij monotonically decreases with j. The rationale behind these choices is that an important component of the relationship between the architecture of neural connections and acoustic output is the number of RA units in the different regions of RA that get recruited by the RA-projecting HVc neurons that sequentially fire during a syllable. Because in our mathematical implementation, the first 5 units in RA project to motor neurons controlling the syringeal muscles, the conditions described above will translate into the frequency of the vocalizations.

Let us analyze the example presented in Fig. 3, which displays the instructions P(t) and T(t) generated when the architecture of connectivity between the units in HVc and those in RA is determined by the matrix presented in Table 1. A sequence of brief bursts is generated in the HVc population. The indices of the columns in Table 1 indicate the temporal order of the bursting in HVc. The matrix displayed in Table 1 implies the following dynamics: in the first 3 time intervals of HVc activity, only the RA neurons involved in the control of expiratory pressure will be active. When the fourth neuron in HVc bursts, the subpopulation of RA neurons in charge of the syringeal tension will begin to activate. When the 7th HVc neuron bursts, both subpopulations will be fully active. Note that as a result of this architecture, the pressure builds up faster than the tension. In Table 2, we show a different connectivity matrix between HVc and RA, and the corresponding instructions T(t) and P(t) are displayed in Fig. 4. In this case, T(t) increases faster than P(t).



View larger version (22K):
[in this window]
[in a new window]
 
FIG. 3. Synthetic motor control instructions in a completely feedforward model. A: bronchial pressure P(t) and labial tension T(t), for the connectivity matrix displayed in Table 1. B: sound pressure generated when these instructions drive the mathematical model of the syrinx, and the sonogram. Frequency is lower at the beginning than at the end of each syllable. Units representing HVc neurons were forced 3 times, to generate 3 syllables.

 

View this table:
[in this window]
[in a new window]
 
TABLE 1. HVc–RA connectivity in simulations shown in Figs. 3, 5, 6, and 10

 

View this table:
[in this window]
[in a new window]
 
TABLE 2. HVc–RA connectivity in simulations shown in Fig. 4

 


View larger version (24K):
[in this window]
[in a new window]
 
FIG. 4. Synthetic motor control instructions in a completely feedforward model, with a different connectivity matrix than the one used in Fig. 2 (see Table 2). A: P(t) and T(t) stand for the bronchial pressure and labial tension, respectively. B: sound pressure of the syllables generated when the instructions displayed in a drive the model of the syrinx, and their sonograms. Note that the frequency is lower at the end than at the beginning of each syllable. The units representing HVc neurons were forced 3 times, to generate 3 syllables.

 
It is possible to understand qualitatively why these different motor instructions will lead to different sounds. To generate sound, the bronchial pressure P(t) has to overcome a certain threshold (Gardner et al. 2001Go). Also the larger the syringeal tension T(t), the higher will be the frequency of the vocalized sound. Therefore a configuration such that the pressure builds up faster than the tension will generate a syllable with an "upsweep" (Fig. 3B), whereas a configuration such that the tension builds up faster will generate a syllable with a "downsweep" (Fig. 4B). We denote as a "downsweep" a vocalization fragment in which the frequency of the sound at the beginning is higher than the frequency of the sound at the end. Based on the results obtained in the 2 cases analyzed so far, we can deduce that a sequence of RA-projecting HVc neurons that each excite approximately the same number of syringeal RA units during the production of a syllable will give rise to a sound of constant frequency.

The observation of HVc neurons driving the activity of RA neurons brings our understanding of the mechanisms behind the production of birdsong a step closer, but it has to be complemented by an understanding of the role of the local circuitry within the process. It has been suggested that the discharge patterns of neurons in RA are the result of interaction between local circuits and input instructions from HVc (Yu and Margoliash 1996Go). It was only recently that the interneurons in RA were characterized in detail (Spiro et al. 1999Go). These long-range inhibitory interneurons project widely throughout RA without exiting the nucleus. They provide inhibition to projection neurons and may play a role in the creation of the patterns of bursts and pauses characteristic of RA during singing (Spiro et al. 1999Go). In our model, these are the 2 main factors in the determination of the phase difference between the tension and pressure gestures driving the syrinx: extrinsic input to RA from HVc and local inhibitory circuitry within RA. To explore the dynamic consequences of the latter factor in the generation of motor signals driving the syrinx, we include inhibition in our simulations.

The role of inhibitory neurons in RA: the level of the syllable

The introduction of an inhibitory unit into the model of RA introduces important changes in the spiking activity of the network. To explore this issue, we describe the activities of projection neurons in RA and how these activities are modified by inhibition. Figure 5 displays spike raster plots of our model neurons for 2 different values of GRIa, which parameterizes the inhibition in RA. In Fig. 5A, the inhibition was set to GRIa = 50 mS/cm2, whereas in Fig. 5B, we used GRIa = 0, which corresponds to the complete feedforward model of interaction between HVc and RA. In both panels, the first 10 rows display the times at which each of the 10 RA-projecting HVc neurons spike. The next 10 rows in each panel display the spike times of the 10 RA projection neurons. Both numerical experiments were generated with the matrix of connections between HVc and RA defined in Table 1. When the inhibition is taken into account, the RA neurons do not display the same level of activity under the influence of the same pattern of activity in HVc.



View larger version (28K):
[in this window]
[in a new window]
 
FIG. 5. Raster plots of the units in our models, for GRIa = 50 mS/cm2 (A), and GRIa = 0 mS/cm2 (B). In both panels, the first 10 rows display the times at which each of the 10 HVc projection neurons spike. Next 10 rows in each panel display the spike times of the 10 RA projection neurons. Both numerical experiments were generated with the matrix of connections between HVc and RA defined in Table 1.

 
Beyond changes in the global activity, the inhibitory coupling of the 2 subpopulations of projection units in RA leads us to expect important changes in the phase difference between the gestures of pressure and tension when intrinsic inhibitory connections in RA are taken into account. Moreover, because the model of the syrinx produces qualitatively different syllables depending on the relative phase between the cyclic activity driving T(t) and P(t), the inhibitory coupling of the subpopulations in RA can translate into important differences in the produced syllables.

To explore this issue, we repeated the numerical simulations that gave rise to Fig. 3, but including global inhibition (GRIa = 75 mS/cm2). In the numerical experiments displayed in both Fig. 3 and Fig. 6, the connectivity between HVc and RA was determined by Table 1. The instructions for P(t) and T(t) that are obtained as inhibition is considered (Fig. 6A) are considerably different from the ones displayed in Fig. 3A. The syllables that are generated when these instructions drive our model of the syrinx are displayed in Fig. 6B. To generate this figure, the parameters T(t) and P(t) are swept cyclically 3 times each. Although the architecture of connectivity between HVc and RA is the same for the cases shown in Figs. 3 and 6, the effect of adding a long-range inhibitory unit in RA is important and translates into the acoustic properties of the syllable. The syllables produced are qualitatively different, as can be seen in the time evolution of the fundamental frequencies. A comparison between the syllables generated with and without inhibition is illustrated in Fig. 7 A.



View larger version (23K):
[in this window]
[in a new window]
 
FIG. 6. Effect of inhibition on the structure of a syllable. A: instructions P(t) and T(t) that our model generates when inhibition is included. B: sound pressure and sonogram generated by these instructions. Architecture of connectivity between HVc units and RA units is defined by the matrix in Table 1 (the same one used to generate Fig. 3). Parameter controlling the level of the inhibition within RA is GRIa = 75 mS/cm2. Units representing HVc neurons were forced 3 times, to generate 3 syllables.

 


View larger version (21K):
[in this window]
[in a new window]
 
FIG. 7. In A, we generate 2 syllables without inhibition, and 2 syllables with GRIa = 75 mS/cm2. Note that the average frequency of the vocalization decreases as the inhibition is increased. In B we display the average frequency of the vocalizations as a function of GRIa.

 
The morphological changes observed in the frequency representation of the song are not the only changes that can be observed. As inhibition is increased, the activity of both RA projection neuron populations is decreased. Because in our model the tension of the syringeal labia depends on the recruitment of units in RA, it can be expected that the average frequency of the vocalizations will diminish as the inhibition is increased. To test this hypothesis, we generated a set of songs for a given architecture of connectivity between HVc and RA, for different values of GRIa. The average frequency of the vocalization was computed as the main peak in the power spectrum of the sound pressure produced by the model. The result is displayed in Fig. 7B. As the inhibition in RA was increased from GRIa = 0 to GRIa = 75 mS/cm2, the average frequency of the vocalization decreased monotonically.

The role of inhibitory neurons in RA: subharmonics

The results described so far emphasized the effect of inhibition at the level of a single syllable. Yet, the complex architecture of nonlinear units used to model the activity of RA and the architecture of units used to emulate HVc activity anticipate the possibility of complex dynamic behavior under the forcing of HVc. To repeat a syllable in our model, it is necessary to repeat a sequence of bursts in HVc. For low enough forcing frequencies, the activity of the projecting units in both HVc and RA will recover from the forcing, and successive syllables will be identical. If the units used to model the RA nucleus were linear devices, their activities would lock their dynamics to the forcing received from HVc for any value of the forcing frequency. Since the units used to model RA, like the biological neurons, are nonlinear devices, it is possible to find more subtle dynamics.

We explored the behavior of our model for different values of the inhibition in RA and different frequencies used to force HVc. In Fig. 8, we display a nontrivial solution. In Fig. 8A we show the pressure P(t) generated by our RA model when the units used to emulate RA-projecting HVc neurons were repeatedly forced with an intersyllabic time of 60 ms. The level of the inhibition within RA was set as GRIa = 75 mS/cm2. Note that the values of the pressure do not repeat themselves after one period of the forcing. The pressure repeats itself after a time interval equal to 2 times the forcing period. This behavior is typical of forced nonlinear systems (Guckenheimer and Holmes 1983Go).



View larger version (23K):
[in this window]
[in a new window]
 
FIG. 8. For high values of the inhibition, the response of RA to the periodic forcing of HVc displays subharmonic responses. In A, we show the computed bronchial pressure P(t) as a function of time, when a sequence of bursts from HVc is applied periodically to RA. Values of the pressure repeat themselves after a time interval equal to 2{tau}, where {tau} is the period of the forcing from HVc. In B, we illustrate the resulting vocalizations, in which an alternation of notes can be observed.

 
When these synthetic motor instructions are used to drive the model of the syrinx, the result is an alternation of similar syllables spanning different frequency ranges, as illustrated in Fig. 8b. In fact, an alternation of syllables of this kind has been reported in recordings of the white-crowned sparrow (Laje and Mindlin 2002Go), although it is not possible so far to associate this feature to the particular mechanism that our model of RA suggests.

A critical ingredient in the generation of this subharmonic behavior in our RA model is the global inhibition. The inhibition leaves the system far away from the "nonexcited" state after a sequence of bursts in HVc. In this way, a new sequence of HVc bursts starting shortly after the previous syllable was generated will find the units in RA in a different dynamic state. This is a necessary condition for subharmonic behavior.

If the intersyllabic time is increased, this effect will not be observed. Now the units in the nuclei will relax to their "nonexcited" dynamic state after the driving stops, and each forcing of HVc will simply lead to a repetition of the same syllable. In Fig. 9 we display the results of repeating our numerical simulation, but with the larger intersyllabic time of 65 ms. Now the response of RA to forcing by HVc is simpler. The activity of RA is locked one to one to that of HVc. For each repetition of an HVc burst sequence, the syllable generated by the syrinx model is identical. Qualitative changes in the dynamics of a system as a parameter is varied are called bifurcations and are a signature of nonlinear behavior. The precise value at which the system changes its behavior is called a bifurcation point.



View larger version (26K):
[in this window]
[in a new window]
 
FIG. 9. Subharmonic nature of the response of RA under the periodic forcing from HVc disappears as the intersyllabic forcing frequency from HVc is decreased. Here the intersyllabic time interval was chosen as 65 ms. Now RA locks in a one-to-one fashion to the driving from HVc. Pressure P(t) and tension T(t) are displayed in A, and the resulting vocalization is shown in B. Note that all syllables are now identical.

 
The role of excitatory connections among projection neurons in RA

If the conductance {Gamma}' characterizing the excitatory interconnections among RA units is large, the activity of each unit is driven both by the direct activity of an HVc unit and by its neighbors in RA. We explored the effect of this coupling in our model. To choose a relevant range of {Gamma}', we designed a numerical experiment in which only one RA unit was driven directly by a unit in HVc. For 2 < {Gamma}' < 4.5 a burst in the HVc unit induces bursting activity in the RA unit that is directly connected to it, and the activity of this RA unit can induce activity only in its first neighbors. For {Gamma}' > 8 the burst in HVc induces activity in the whole population of RA units by the excitatory RA–RA connections. Therefore we explored the dynamics of our model when {Gamma}' was varied within the range 0–8 mS/cm2.

Once an interesting range of {Gamma}' had been chosen, we repeated the numerical simulations with the connectivity between HVc and RA prescribed by Table 1. The results of these numerical simulations are summarized in Fig. 10. For larger values of the strength of excitatory RA–RA interconnections, a larger number of units activate when the same initial RA units are driven by HVc. This results in larger values of both P(t) and T(t). Given that T(t) is translated into the frequency of the resulting vocalization, a direct consequence of increasing the strength of excitatory interconnections is an overall increase in the average frequency of the generated syllables. This effect is illustrated in Fig. 10A. However, we did not find any qualitative change in the morphology of the vocalizations. We show the qualitatively similar syllables generated with {Gamma}' = 2.5 (Fig. 10B) and {Gamma}' = 5.5 (Fig. 10C), given the architecture of connectivity of Table 1. This negative result was obtained under the hypothesis that every excitatory unit in RA projects onto a similar number of excitatory units. Under this hypothesis, connections between excitatory units in RA lead to an amplification of the acoustic features determined by the connections between HVc and RA.



View larger version (21K):
[in this window]
[in a new window]
 
FIG. 10. Increasing the strength of the excitatory interconnections among units in the model of RA changes the average frequency of the vocalizations (A), which results in an increase in the number of RA units being recruited under the command of each HVc unit within the population controlling the tension T(t). Morphology of the syllables did not show qualitative changes. In B, {Gamma}' = 2.5; in C, {Gamma}' = 5.5. HVc–RA connectivity matrix is shown in Table 1.

 

    DISCUSSION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
The recent descriptions of the spiking patterns in HVc and RA, the main nuclei involved in the production of birdsong, suggest that the dynamics at the time scale of a syllable originates in HVc (Yu and Margoliash 1996Go). It is natural then to assume that the connectivity between HVc and RA, where the specific instructions driving the syrinx are generated, plays a major role in determining the acoustic features of syllables. However, our results support the view that these features may emerge from the interaction of the connectivity between HVc and RA and the internal architecture of RA. We have explored the effect of both long-range inhibition and local excitatory connections among RA projection neurons in a model of the song motor pathway and syrinx. Particularly important in our model is the global inhibition, which constrains the activity of the different RA neurons projecting to muscles. Our model lacks details of the HVc and RA nuclei, but we suggest that the importance of long-range inhibition will persist in an augmented model that accounts for this internal structure.

Experimentally we would expect that a reversible blocking of the inhibition in RA by the application of a GABAA receptor antagonist will dramatically alter the acoustic features of the song. This is consistent with recent experiments in which microinjections of bicuculline in RA led to changes in the song (Vicario and Raksin 2000Go). This would not be the case in a complete feedforward model with no interneurons modulating the dynamics.

Spontaneous spiking of RA neurons has been recorded in vitro (~15 Hz; Mooney 1992Go; Spiro et al. 1999Go) and in vivo after the end of a song (Yu and Margoliash 1996Go). This nonsinging, spiking activity is qualitatively different from the higher-frequency bursting that is generated in RA neurons during singing (Hahnloser et al. 2002Go; Yu and Margoliash 1996Go). In our model, inhibitory RA interneurons suppress this spontaneous activity during singing. Spiking is then evoked in RA projection neurons only when the direct excitatory input from HVc is strong enough to overcome the inhibition.

The model of the syrinx that we use to generate synthetic song produces qualitatively different syllables for different values of the phase difference between the cyclic instructions controlling P(t) and T(t). Therefore a change in the inhibition, which in principle affects different RA units in different ways, results in important changes in the song. Because inhibition also affects the total activity of the RA nucleus, the average frequency of the vocalizations, proportional to the average value of the tension, changes with the level of global inhibition in RA. The intrinsic connections within RA are not the only factor determining the phase delays between respiratory and vocal muscle activity. Connections between the nucleus RAm (which contains premotor neurons projecting to expiratory motor neurons) and nXIIts have been recently described (Sturdy et al. 2003Go). These connections may play an important role in coordinating respiratory activity with the activity of motor neurons controlling the dorsal syringeal muscles (not included in our model; Laje et al. 2002Go).

Finally, our model predicts complex dynamic responses of RA under the forcing of HVc. We illustrated this behavior by showing that periodic forcing of HVc can produce an alternation of different syllables. This complex response emerges from the interaction between simple instructions and a complex connectivity of nonlinear units within RA.

Birdsong is a test bed for the learning of complex behavior, in that songbirds require a tutor for learning their vocalizations, just as humans do (Brainard and Doupe 2002Go; Doupe and Kuhl 1999Go). Our simple model suggests that in addition to the study of how the HVc -> RA projections may be modified by synaptic plasticity under control of the anterior forebrain pathway during song learning, it will be important to consider the internal architecture of RA.


    APPENDIX
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
The dynamics of the units representing neurons in our model is governed by the HH equations. The parameter values used in our simulations are CM = 1 µF/cm2 for the membrane capacitance, and gNa = 215 mS/cm2, gK = 43 mS/cm2, and gL = 0.813 mS/cm2 for the maximal conductances. The values of the reversal potentials are ENa = 50 mV, EK = –95 mV, EL = –64 mV.

The parameters {alpha}X and {beta}X, present in the equations ruling the dynamics of the activation and inactivation variables, are






(A1)
The other constants used in our calculations are {alpha} = 0.15 ms–1, {beta} = 0.2275 ms–1, Vp = 10 mV, GHE = 17.7 mS/cm2, GHI = 7.5 mS/cm2, GHVcI = 10.5 mS/cm2, ERevE = 0 mV, ERevI = –80 mV, {Gamma} = 18.62 mS/cm2, {Gamma}' = 3.5 mS/cm2, GRIa = 75 mS/cm2, GRAI = 7.5 mS/cm2, {alpha}1 = 3.75 x 108, {alpha}0 = 5.93 x 105, {beta}1 = 30, {beta}0 = 80, and C = 2.0 x 108.


    GRANTS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
This work was partially funded by the University of Buenos Aires, ANCyP PICT 03-08133, Fundacion Antorchas; by the U.S. Department of Energy, Office of Basic Energy Sciences, Division of Engineering and Geosciences, under Grants DE-FG03-90ER14138 and DE-FG03-96ER14592; by a grant from the National Science Foundation, NSF PHY0097134; by a grant from the Army Research Office, DA-AD19-01-1-0026; by a grant from the Office of Naval Research, N00014-00-1-0181; and by National Institutes of Health Grants R01 NS-40110-01A2 and R01 GM-65830-01.


    ACKNOWLEDGMENTS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 RESULTS
 DISCUSSION
 APPENDIX
 GRANTS
 ACKNOWLEDGMENTS
 REFERENCES
 
Constructive input from two anonymous referees is acknowledged.


    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.

Address for reprint requests and other correspondence: H.D.I. Abarbanel, Institute for Nonlinear Science, University of California, San Diego, 9500 Gilman Dr., La Jolla, CA 92093-0402 (E-mail: hdia{at}jacobi.ucsd.edu).


    REFERENCES