|Home | About | Journals | Submit | Contact Us | Français|
A biologically detailed model of the binaural avian nucleus laminaris is constructed, as a two-dimensional array of multicompartment, conductance-based neurons, along tonotopic and interaural time delay (ITD) axes. The model is based primarily on data from chick nucleus laminaris. Typical chick-like parameters perform ITD discrimination up to 2 kHz, and enhancements for barn owl perform ITD discrimination up to 6 kHz. The dendritic length gradient of NL is explained concisely. The response to binaural out-of-phase input is suppressed well below the response to monaural input (without any spontaneous activity on the opposite side), implicating active potassium channels as crucial to good ITD discrimination.
Computational neuroscience is a discipline that aims to understand how information is processed in the nervous system by developing formal models at many different structural scales. Using a multidisciplinary approach, the biological systems are modeled incorporating the constraints imposed by psychophysical, neurobiological, and formal computational analyses. To fully comprehend the very diverse and highly specialized tasks of the nervous systems, multiple levels of analysis have to be addressed ranging from detailing the biophysical mechanisms responsible for these computations to evaluating the overall system performance (see Hawkins and McMullen 1996). The end product of a computational analysis should be a sufficiently specified model, internally consistent and complete enough to enable formal mathematical characterization or computer simulation. Such properties ensure that the computed outcome accurately reproduces the functions of the modeled system. It is often the case that specific features of the real world are not fully understood by the modeler, who therefore must guess about them when designing and testing the model. Conjectures are to be made, nonetheless, following psychophysical and neurobiological constraints, as the theoretical outcomes of the model would otherwise fail to capture the nature of biological computation. A certain degree of ignorance about the real phenomenon is implicit in the modeling process because it is indeed the aim of such a process to gain a deeper insight into the mechanisms governing the object of study.
Characterized by processing information up to two orders of magnitude faster than the rest of the brain, the auditory system is responsible for providing a representation of the pattern of incoming sound signals and for specifying the locations of the corresponding sound sources relative to the listener. Specifically, sound localization can be computed from interaural temporal and level differences over the entire spectrum of the signal from the sound coming to both ears. The picture of what is accomplished in the brain stem auditory nuclei as acoustic information ascends the auditory pathway is still fragmentary, but there is no doubt that auditory neurons are specialized for precisely timed signaling.
The present study focuses on avian sound localization. There are strong similarities in the organization of the brain stem auditory nuclei among diverse animal classes, but the avian system is believed to have a simpler structure, and the function of each nucleus involved is better understood, than the mammalian analog medial superior olive (Rubel and Parks 1988; Oertel 1999). Timing information in sound is captured in the auditory sensory epithelium and conveyed to the auditory nerve, which projects to the cochlear nuclei. One branch of the auditory nerve enervates nucleus angularis (NA) and another enervates nucleus magnocellularis (NM), whose afferents project to nucleus laminaris (NL). NL is dominated by cells that extend their dendrites dorsally and ventrally in a bipolar fashion. Each cluster of dendrites receives inputs from one ear: the ipsilateral NM projects onto the dorsal side, whereas the contralateral NM innervates the ventral side. Additionally, neurons within both nuclei are arranged tonotopically: their neuronal location depends on the frequency band to which they are tuned (Rubel and Parks 1975). Evidence is strong that NM and NL are adapted for analysis of interaural time differences (ITDs) (Rubel and Parks 1988), whereas all other information is processed through NA. Detection and computation of ITDs is crucial to azimuthal sound localization by the brain.
The most prominent model of binaural interaction was put forth by Jeffress (1948) as a mechanism for sensitivity to interaural time delays. The model postulates the existence of neuronal fibers that contain temporal information about the waveform of the acoustic stimulus in their discharge pattern. These fibers project with ladder-like branching patterns to cells in a binaural nucleus that only discharge when receiving coincident spikes from their monaural afferents. The varying lengths or thicknesses of these fibers introduce internal delays that compensate for the relative acoustic delay of incoming signals between both ears. The model therefore computes ITDs by means of delay lines and coincidence detectors and transforms a time code into a space code. Anatomical and physiological evidence shows that the cochlear nucleus is the likely site of Jeffress’ model: NM axons are hypothesized to work as delay lines, and NL cells are thought to play the role of coincidence detectors (Carr and Konishi 1990; Overholt et al. 1992; Joseph and Hyson 1993).
In NL, strong functional hypotheses can be formulated and biological data can be feasibly obtained from experiments in vitro and in vivo. As a result, NL, primarily in chick and barn owl, has already been the objective of many anatomical (Rubel and Parks 1975; Smith and Rubel 1979; Smith 1981; Deitch and Rubel 1984, 1989a, b) and electrophysiological experiments (Carr and Konishi 1990; Warchol and Dallos 1990; Overholt et al. 1992; Pena et al. 1996; Reyes et al. 1996; Viete et al. 1997; Bruckner and Hyson 1998; Funabiki et al. 1998; Yang et al. 1999; Monsivais et al. 2000; Pena et al. 2001; Kuba et al. 2002; Cook et al. 2003). This wealth of information has also made it the target of several modeling efforts (Grün et al. 1990; Agmon-Snir et al. 1998; Dasika et al. 1999; Simon et al. 1999a, b; Dasika et al. 2001; Simon et al. 2001a, b; Cook et al. 2003; Carr et al. in press). Models have yielded significant results in explaining some particular features of NL: the role of dendrites in auditory coincidence detection (Agmon-Snir et al. 1998), phase-locking capabilities and synaptic location (Simon et al. 1999a, b), and short-term synaptic plasticity as an adaptive mechanism for preserving ITD information (Cook et al. 2003).
The present model incorporates the varied biological data available, and the dendritic length gradient in particular, to construct a biophysically realistic model of neuronal coincidence detection that is valid across all NL. This neuron-level approach aims to faithfully characterize the morphology, electrical properties, and channel dynamics that underlie the characteristic behavior of NL cells throughout the nucleus.
In addition to those models of NL referred to above, there are related models of coincidence detection in mammalian medial superior olive (MSO), or models general enough to apply equally well to both MSO and NL (Colburn et al. 1990; Han and Colburn 1993; Brughera et al. 1995; Brand et al. 2002; Colburn and Zhou 2002; Svirskis et al. 2002; Svirskis and Rinzel 2003; Zhou and Colburn 2003). For a dated but valuable review, see Colburn (1996). Additionally, the model discussed here has also been generalized to other animals: crocodilians (Soares et al. 1999), and emu (Soares et al. 2001, 2002).
The model presented here presumes the two-part model of Jeffress (1948) with delay lines and coincidence detectors. In the model, NL neurons are studied as coincidence detectors, and it is assumed that delays are external to the NL neurons (e.g., adjusted by fine-tuning axonal delays in the incoming NM neurons). An entire two-dimensional lamina is modeled, with one dimension corresponding to varying ITDs (the ITD axis) and the other corresponding to varying best frequency (BF) (the tonotopic axis). The model's computer code is written in the neuronal modeling language NEURON (Hines and Carnevale 1997, 2000, 2001).
The model is primarily of NL in chick since the biological data used to parameterize the model are taken primarily from chick studies. The large quantity of data from chick NL is not enough, however, to describe all crucial parameters, and data must be taken from other sources as well. The greater part of the remaining parameters come from studies of chick NM (which is similar in some electrophysiological properties to NL) and from barn owl NL (which is similar in some electrophysiological and anatomical properties to chick), but some from further afield. It is hoped, of course, that those parameters are not absolutely crucial to the model's gross behavior. Additionally, to the extent that barn owl NL and chick NL are sufficiently similar, generalizing the model to owl NL from chick is straightforward.
Each individual neuron in the two-dimensional array is modeled with multiple sections: a soma; an axon with an axon hillock, a myelinated segment, and a node of Ranvier; and an adjustable number of dendrites. Each section has an adjustable number of equipotential compartments. Each dendrite has an adjustable distribution and number of excitatory synapses.
The excitatory synapses are modeled as time-varying conductances, with adjustable time constant, peak conductance, reversal potential, and refractory period. The excitatory synapses are driven by simulated NM neurons whose firing is modeled as a constant frequency (periodic) inhomogeneous Poisson process with adjustable average firing rate, vector strength (Goldberg and Brown 1968), frequency, and phase. Typically the ipsilateral and contralateral dendrites receive inputs with the same average firing rate, vector strength, and frequency but different phases; it is straightforward, though, to modify the simulated stimuli, e.g., for ipsilateral monaural stimulation with adjustable contralateral spontaneous activity, or stimulation with nonzero interaural level difference.
The soma also has inhibitory synapses, parameterized similarly to the excitatory synapses. Unlike the fast, phase-locked inhibitory input in mammalian MSO (Brand et al. 2002), inhibitory input to NL is slow and not phase-locked (Funabiki et al. 1998; Yang et al. 1999; Monsivais et al. 2000). The input to the inhibitory synapses is from a modeled superior olivary nucleus (SON) (Takahashi and Konishi 1988; Carr et al. 1989; Lachica et al. 1994). The modeled SON can be either a simple integrate-and-fire neuron whose inputs come from a constant BF array of NL neurons (i.e., of varying ITD) or as a group of multicompartment neurons receiving input from neurons of constant BF in NL, NM, and a simulated nucleus angularis (Rodriguez-Esteban 2002). The inhibitory inputs project back equally to all neurons of the same BF in NL.
Channel types include active potassium channels in both high- and low-voltage-activated forms (and Hodgkin-Huxley delayed rectifier), active sodium channels, and passive leak channels. Voltage-dependent channels are specified by Hodgkin-Huxley-like parameters (Hodgkin and Huxley 1952).
The soma is modeled as a truncated ellipsoid of length and diameter 15 μm (Smith and Rubel 1979; Smith 1981; Parks et al. 1983; Deitch and Rubel 1989b). Dendritic diameter is 4 μm, and dendritic length is a function of position along the tonotopic axis (Rubel and Parks 1975; Smith and Rubel 1979),
which is simply a linear function of logarithmic frequency whose smallest length (20 μm) occurs at and above ~ 2500 Hz, and whose longest length (400 μm) occurs at and below ~ 300 Hz. The axon begins with a 30-μm-long, 8-μm-diameter initial segment (Carr and Boudreau 1993). It has a myelinated section with an attenuation of conductance and capacitance by a factor of 80 over those used throughout the rest of the model (standard values of Ri = 200 Ω-cm and Cm = 1 μF=cm2 (Koch 1999)).
Simulated NM neurons driving the NL excitatory synapses spike at an average firing rate independent of tonotopy (Carr and Konishi 1990;Warchol and Dallos 1990). The vector strength (phase-locking index) of the firing, modeled after Warchol and Dallos (1990) for chick and Köppl (1997) for barn owl, is a function of stimulus frequency
which is simply the linear function of logarithmic frequency whose weakest vector strength VSH occurs for frequencies at and above fVSH and whose strongest vector strength VSL occurs for frequencies at and below fVSL. Specific values are given in Tables 14 and and15.15. The firing pattern is modeled as an inhomogeneous Poisson process with an instantaneous firing probability proportional to an exponentiated sinusoid, also known as the von Mises distribution (Fisher 1993), but with a refractory period. The probability of firing at any time step after the refractory period is given by the product of the time-step interval and
where N is a normalization constant and γ is a monotonic function of vector strength.
The NM firing pattern is computed on the fly (Hill 1976, 1977a, 1977b,1981a, 1981b) so that the simulated run time does not need to be specified in advance. The NM spike probability rate gives the average NM firing rate approximately but reduced because of the minimum reexcitation time. The values used here produce a typical average NM firing rate of 350 spikes/s.
The conductance of all synapses has the time profile of an alpha function
but the time constants and strengths differ substantially between the excitatory and inhibitory cases (cf. Tables 2 and and8).8). The inhibitory parameters for these simulations use a reversal potential hyperpolarized with respect to the resting potential (–80 mV), but a depolarized reversal potential of –60 mV (accompanied by a doubling of the conductance to give the same effective driving force) produced almost identical behavior (see Bruckner and Hyson 1998; Funabiki et al. 1998; Yang et al. 1999; Monsivais et al. 2000).
The passive channel properties were taken from the voltage-current relation measured by Reyes et al. (1996) (note that time constants change substantially for developing chick embryos (Kuba et al. 2002), so the age of the modeled chick can be important).
The active potassium channels fall into two main types, KLVA and KHVA, activated, respectively, at low voltage (hyperpolarized) and high voltage (depolarized) thresholds. KLVA is believed to be of the Kv1 subfamily (Carr et al. 2001). It shows strong outward rectification near rest and is dendrotoxin-I (DTX) and 4-amino-pyridine (4-AP) sensitive (Kuba et al. 2002). It presents rapid activation near rest (low threshold), and its deactivation is slow or partial. Conductances and kinetic parameters for KLVA are taken from studies in NM (Rathouz and Trussell 1998), modified to conform to measured behavior in NL (Reyes et al. 1996), e.g., activation threshold and voltage change to small current injections. KHVA is believed to be Kv3.1 (Parameshwaran et al. 2001) and is characterized by positively shifted voltage dependencies and very fast deactivation rates (Parameshwaran et al. 2001; Rudy and McBain 2001). These properties allow this channel to enable fast repolarization of action potentials without compromising spike initiation or height; this leads to quick recovery after an action potential. It is tetraethylammonium (TEA) sensitive. Conductances and kinetic parameters for KLVA are also taken from studies in NM (Rathouz and Trussell 1998), fit to a linear gating model (see Tables 10 and and1111 for parameter values).
Axonal sodium and potassium channels use simple Hodgkin-Huxley models with conductances adjusted to give NL firing rates consistent with Reyes et al. (1996) (see Tables 12 and and1313 for parameter values). Kinetic dependencies on temperature, sorely needed since physiological experiments are conducted at a variety of temperatures, are essentially unknown and are set to the plausible values for all channels.
The time-step interval of 12.5 μs was selected by examining spike rates and timing for a limited number of cases and was chosen as the largest interval below which firing properties did not change. In particular, even at high stimulus frequencies, smaller intervals do not change firing rates, presumably because the low vector strength jitter dominates any artifacts due to time-step artifacts.
Model runs are typically done on the entire two-dimensional array in parallel, with frequencies varying from 350 Hz to 2 kHz in half-octave steps and interaural phase difference varying from 0° to 180° in 22.5° steps (Fig. 1). NL spike times are noted and average firing rates and vector strength calculated. To analyze more properties more quickly, the phase is often limited to the two values of 0° and 180°. Then the “ITD discrimination index” is defined as: one minus the ratio of the inputs-out-of-phase NL firing rate to the inputs-in-phase NL firing rate. This index is near one when the in-phase firing rate is high and the out-of-phase firing rate is low, but drops to zero when the in-phase and out-of-phase firing rates are indistinguishable. Figures 2 and and33 show plots of ITD discrimination index as a function of frequency.
As shown in Fig. 1, the NL neurons fire well when their ipsilateral and contralateral inputs are in phase. At low frequencies the dendritically extensive NL neurons fire poorly as the interaural phase difference (IPD) increases – the hallmark of a good ITD-sensitive coincidence detector. At high frequencies the dendritically compact NL neurons still fire at large IPDs, meaning that their ability to discriminate between ITDs worsens. There is also a significant enhancement of the vector strength of the output over the vector strength of the input, even at nonzero IPD, a potential flaw in the model.
As shown in Fig. 2, ITD discrimination is good until ~2 kHz. The sharp phase tuning at low frequencies, and specifically the firing rate dropping to zero, appears as ceiling effect on the ITD discrimination index at those low frequencies. The lower panel shows the same neurons, but receiving inputs from NM with vector strength enhanced to match that of barn owl. In this case, ITD discrimination is significantly stronger, up to ~4 kHz. By fine-tuning parameters of the model (particularly by adding more dendrites), this can be pushed up to ~6 kHz (not shown), but this is still not as effective as the barn owl, which can discriminate up to ~ 10 kHz (Carr and Konishi 1990).
Synaptic sublinearity occurs due to finite synaptic reversal potential, as predicted by cable theory (Rall 1959; Agmon-Snir et al. 1998), and it also appears functionally as a drop in ITD discrimination when the excitatory reversal potential is depolarized (not shown). The most interesting effect of synaptic sublinearity can be seen when the dendritic length gradient is turned off or, equivalently, neurons are allowed to receive stimuli at frequencies (unphysiologically) far from their BF. As shown in Fig. 3, for every stimulus frequency there is a dendritic length longer than that at which performance no longer increases. The effect is most pronounced at lower frequencies. The optimal length of the dendrite increases with decreasing frequency. Thus, the dendritic length gradient is predicted under the assumption that dendritic length is minimized while optimizing for ITD discrimination. This is in agreement with Agmon-Snir et al. (1998).
The model also demonstrates a subtraction nonlinearity, due to the presence of KLVA channels acting as a current sink, and it also appears functionally as a drop in ITD discrimination when the KLVA conductance is reduced (not shown). This has an interesting additional manifestation: the out-of-phase rate is suppressed relative to monaural rate, even in the absence of spontaneous activity from the contralateral side (Fig. 4). This is because the high input activity from the out-of-phase case activates KLVA and suppresses depolarization sufficiently to reach spike threshold. This result is in agreement with experimental findings (Goldberg and Brown 1968; Yin and Chan 1990) but without the spontaneous activity needed by more linear models (e.g., Colburn et al. 1990). Note that for the monaural case the average current injected via excitatory synapses is roughly half that of the binaural out-of-phase case, even though the spike rate of the former is greater than that of the latter. Furthermore, by comparing the firing rates for monaural inputs among different neurons (with enhanced KLVA or missing the opposite dendrite), it can be seen that KLVA plays an important role in reducing the firing rate, and, in particular, the opposite dendrite can act as a current sink.
A finely detailed, yet still highly functional, model of chick NL has been constructed. Typical parameters allow ITD discrimination up to 2 kHz, and enhancements for barn owl allow ITD discrimination up to 6 kHz. There are two nonlinearities that aid ITD discrimination: intradendritic inputs sum sublinearly, and interaction with KLVA subtractively suppresses outof-phase inputs. The response to monaural input does not require any spontaneous activity from the contralateral side. The dendritic length gradient of NL is predicted from optimization of ITD discrimination.
Thanks to Raul Rodriguez Esteban, Kate MacLeod, Daphne Soares, and Suchitra Parameshwaran. Support was given to JZS by National Institutes of Health Grant No. R03DC04382 and NSF Grant No. 972033 and to CEC by National Institutes of Health Grant No. DCD00426.
The values tabulated here are the values used unless explicitly countermanded by the text.