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: 95/5/3113    most recent 00979.2005v1 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 PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via ISI Web of Science (16) Citing Articles via Google Scholar Google Scholar Articles by Gold, C. Articles by Buzsáki, G. Search for Related Content PubMed PubMed Citation Articles by Gold, C. Articles by Buzsáki, G.

On the Origin of the Extracellular Action Potential Waveform: A Modeling Study

¹Computation and Neural Systems, California Institute of Technology, Pasadena, California; ²Merck Research Laboratories, West Point, Pennsylvania; and ³Center for Molecular and Behavioral Neuroscience, Rutgers University, Newark, New Jersey

Submitted 19 September 2005; accepted in final form 24 January 2006

	ABSTRACT

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
Although extracellular unit recording is typically used for the detection of spike occurrences, it also has the theoretical ability to report about what are typically considered intracellular features of the action potential. We address this theoretical ability by developing a model system that captures features of experimentally recorded simultaneous intracellular and extracellular recordings of CA1 pyramidal neurons. We use the line source approximation method of Holt and Koch to model the extracellular action potential (EAP) voltage resulting from the spiking activity of individual neurons. We compare the simultaneous intracellular and extracellular recordings of CA1 pyramidal neurons recorded in vivo with model predictions for the same cells reconstructed and simulated with compartmental models. The model accurately reproduces both the waveform and the amplitude of the EAPs, although it was difficult to achieve simultaneous good matches on both the intracellular and extracellular waveforms. This suggests that accounting for the EAP waveform provides a considerable constraint on the overall model. The developed model explains how and why the waveform varies with electrode position relative to the recorded cell. Interestingly, each cell's dendritic morphology had very little impact on the EAP waveform. The model also demonstrates that the varied composition of ionic currents in different cells is reflected in the features of the EAP.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
Extracellular action potential (EAP) recordings form one of the primary means for studying the activity of the intact brain. Multielectrode arrays and spike-sorting algorithms have advanced to the point where hundreds of neurons can be reliably recorded in a single experiment (e.g., see Csicsvari et al. 2003; Wilson and McNaughton 1993). Typically, EAP recordings are used only to determine whether and when neurons have spiked, under the assumption that the actual waveform of individual EAPs does not convey any information. At the same time, average EAP waveforms are known to exhibit a range of characteristic features when observed on a millisecond timescale, and these variations can be used to distinguish between different neuronal classes (Csicsvari et al. 1999; Mountcastle et al. 1969) as well as individual neurons within classes (e.g., Quian Quiroga et al. 2004). However, there have been few attempts to systematically study the causes of the variability in EAP waveforms either through experimental work or through computer modeling (Buzsáki et al. 1996; Holt and Koch 1999; Moffit and McIntyre 2005; Quirk et al. 2001; Rall 1962). The paucity of studies probably results from both insufficient data and inadequate techniques with which to model EAPs in a meaningful way. It is generally assumed that the main source of EAP variability across different recordings arises from random positioning of the electrode and the morphology of the neuron—a relationship that does not provide useful information about the intracellular state of the cell. Here we show that an accurate computational model of the EAP can shed light on the source(s) of variability in recorded EAP waveforms and can contribute to the analysis of some outstanding questions in the interpretation of EAP recordings. We investigate the effects of cellular morphology, the cell's spatial dimensions, and differential expression of various ionic channels on the waveform of the EAP. We take advantage of recent developments in computational modeling to predict EAPs resulting from simulated neurons at a high level of detail (Holt 1998) and the massive increase in available computing power since the theory was developed in the 1960s (Plonsey 1969; Rall 1962). We use a unique data set (Henze et al. 2000) to reproduce the precise conditions for the generation of a set of intra- and extracellular action potentials recorded in vivo.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
Computational methods

The extracellular potential induced by a spike in a neuron was calculated in two distinct stages. First, we computed the transmembrane currents for a pyramidal neuron model on the basis of standard one-dimensional (1-D) cable theory (e.g., Koch and Segev 1999). Second, we used those currents to compute the extracellular potentials as described below.

CALCULATION OF EXTRACELLULAR POTENTIALS AND THE LINE SOURCE APPROXIMATION. It was previously demonstrated that the neuropil is well modeled by an isotropic volume conductor in which the capacitive effects of the media are negligible in the frequency range of interest to us (1–3,000 Hz). That is, we can satisfactorily describe the extracellular milieu by a purely ohmic conductivity (units of cm) (e.g., see Holt 1998; Plonsey 1969). Under these circumstances, the electric potential in the extracellular space is governed by Laplace's equation

(1)

where is the extracellular potential. At the boundaries, (1/) = J_m, where J_m is the transmembrane current density and is the extracellular resistivity. For a single point source of amplitude I in an unbounded isotropic volume conductor, the solution is dual to the classical physics problem of point charges in free space (Coulomb's law)

(2)

where I is a point source of current and r is the distance from the source to the measurement. Multiple current sources combine linearly by the superposition principle. In real neurons, membrane currents are distributed over elongated cylindrical processes, whose length considerably exceeds their width. The line source approximation (LSA) (Holt and Koch 1999) makes the simplification of locating the transmembrane net current for each neurite on a line down the center of the neurite. By assuming a line distribution of current, the resulting potential from Eq. 2 has a straightforward analytic two-dimensional (2-D) solution in cylindrical coordinates. For a single linear current source having length s, the resulting potential (r, h) is given by

(3)

where r is the radial distance from the line, h is the longitudinal distance from the end of the line, and l = s + h is the distance from the start of the line. Holt and Koch (1999) analyzed the accuracy of the LSA and found it to be highly accurate except at very close distances (i.e., about 1 µm) to the cable (see also Rosenfalck 1969; Tryanova and Enriquez 1991). Because extracellular recording electrodes are typically many micrometers away from neurons, we can use the LSA to calculate extracellular potentials.

The steps in the model are as follows. First, we computed transmembrane currents for a particular neuron with its complement of ionic currents (see following text) using the NEURON Simulation Environment (Hines and Carnevale 1997), assuming that the extracellular potential was constant and equal to zero. In a second step, we used the LSA to compute the extracellular potential at a select number of locations from the transmembrane currents using a custom-written Matlab program. We assumed that the previously calculated transmembrane currents would not be affected by the small changes in extracellular potential (<<1 mV). One could refine this estimate on the basis of an iterative procedure, although this does not significantly affect the numerical results (Holt 1998).

CALCULATION FOR INHOMOGENEOUS RESISTIVITY. The LSA assumes an extracellular medium that is homogeneous. However, recent measurements of CA1 have found that the pyramidal cell body layer has approximately double the resistivity of the surrounding stratum radiatum and stratum oriens ( = 640, 260, 290 cm, respectively) (López-Aguado et al. 2001). Furthermore, these baseline resistivities may be increased by as much as 50% during periods of high activity. Because the inhomogeneity constitutes an approximately planar layer, we can use the "method of images" (Maxwell 1881; Weber 1950) to calculate its impact. Three layers of differential conductivities (₁, ₂, and ₃), separated by two parallel planes, are exactly solved by an infinite series of images with decreasing magnitudes of the form

(4)

where n and m increase for the more distant images. The magnitude of the images decline rapidly at the same time as the distance to the images increases; in practice, the infinite series is well approximated by the first few terms (results shown here use n, m <5). Including the exact calculation for the impact of the high resistivity cell body layer has only a modest impact on extracellular voltages compared with the assumption of homogeneous conductivity equal to that of the stratum oriens or the stratum radiatum. The amplitude of the extracellular action potential increases by a maximum of around 40% at the center of the region of high resistivity, but the impact falls off rapidly with distance. Because of the modest impact on the results and the complexity of the calculation method, we defer details of the solution and our method of verification to a future publication.

Experimental methods

Simultaneous intracellular and extracellular recordings of CA1 neurons in vivo were reported previously in Henze et al. (2000) and we briefly review the methods here. The extracellular electrodes were of three types: 1) single, 60-µm-diameter wires, 2) "tetrodes" as described in Gray et al. (1995), or 3) silicon electrode arrays with six recording sites spaced 25 µm apart, as described in Henze et al. (2000). During numerous attempts to obtain stable intracellular recordings from cells also recorded by the extracellular electrode, Henze et al. (2000) obtained recordings from 38 neurons: three recorded with single-wire electrodes, 22 recorded with tetrodes, and 13 recorded with silicon probe arrays.

Recordings were wideband filtered at either 1 Hz to 3 kHz or 1 Hz to 5 kHz. Averages of the EAPs were made by sampling from the extracellular recording at times triggered by the intracellular spike. In preparing averages for comparison to the model, we used only recordings from the beginning of the session until the cell started to depolarize significantly (>5–10 mV) arising from the shunt current from electrode impalement. The number of spikes available for the average range from a few hundred to a few thousand. After intracellular recordings were complete, cells were injected with biocytin and the rats were killed and their brains sliced, stained at 60 µm, and preserved in slides. Of the 38 recorded cells, 17 cells were stained well enough for reconstruction. In these cases the complete three-dimensional (3-D) structure of recorded cells was measured using the Neurolucida System and then used as the basis for compartmental simulations. In cases where the extracellular electrode track left some visible mark of its location (i.e., blood or debris) this was also measured and used to estimate the electrode location for comparison with the computer simulations. Visible electrode tracks were found in the CA1 area for seven cells and tracks were found in the overlying cortex only for another three cells.

We also used a larger sample of EAP recordings (n = 307) with no coincident intracellular electrode recordings as a reference set for comparison (method as described in Csicsvari et al. 2003), to more accurately estimate the frequency of EAP features observed in the small sample of simultaneous recordings and to confirm that observed features in the simultaneous recordings were not artifacts of intracellular impalement.

Simulation methods

Single trials of standard 1-D compartmental simulations were performed for each reconstructed cell within NEURON (Hines and Carnevale 1997) and compared with the average simultaneous intra- and extracellular recordings. The average number of compartments was around 250, based on a 3-D reconstruction that contained nearly 2,500 measurements of dendrite diameters and locations. The time steps of the simulation were varied by the CVODE method (Hines and Carnevale 2001). During the simulation, membrane currents for all compartments of the cell were saved at intervals of about 0.05–0.1 ms to calculate extracellular potentials.

PASSIVE PARAMETERS AND SPINES. The intracellular resistivity was set to R_i = 70 cm (Stuart and Spruston 1998). The value of this parameter had an important impact on the resulting extracellular potential amplitude. Simulations with higher values for R_i resulted in potentials that were too small to match the recording and histology data. The membrane resistance was set to R_m = 15 kcm² (Spruston and Johnston 1992), to account for the net effect of in vivo synaptic conductances in reducing the membrane resistance without actually modeling the detailed time dependency of synaptic activity (Destexhe and Paré 1999). The specific capacitance was set to C_m = 1 µF/cm² (Koch 1999). The reversal potential for the passive leak mechanism was set to V_rest = –65 mV.

Dendritic spines are accounted for by adjusting the passive membrane parameters R_m and C_m, decreasing the former and increasing the latter by a factor f given by the normalized spine area (Major et al. 1994). Specific spine density estimates are taken from Megias et al. (2001) (for details, see APPENDIX A).

These choices resulted in average somatic input resistances of 31.8 ± 6.5 M, which is in agreement with previous measurements of CA1 input resistances in vivo (48.4 ± 11 M; Henze and Buzsáki 2001) and compatible with the notion that, because of constant synaptic bombardment in vivo, the input resistance is as much as 80% lower (Destexhe and Paré 1999) than that in vitro (Spruston and Johnston 1992).

ACTIVE IONIC CURRENTS. The model includes Hodgkin–Huxley style kinetic models for 12 different voltage-dependent ionic currents carried by Na⁺, K⁺, and Ca²⁺ ions. The following ionic currents included in the model.

• Fast, inactivating Na⁺ (axonal and soma/dendritic varieties: I_Na⁺Ax, I_Na⁺SD) (Colbert and Pan 2002; Magee and Johnston 1995; Martina and Jonas 1997)
  
• A type K⁺ (proximal and distal varieties, I_K⁺AProx, I_K⁺ADist) (Hoffman et al. 1997; Klee et al. 1995)
  
• AHP (afterhyperpolarization) type Ca²⁺-dependent K⁺ (I_K⁺AHP) (Williamson and Alger 1990)
  
• C type voltage and Ca²⁺-dependent K⁺ (I_K⁺C) (Lancaster and Nicoll 1987; Yoshida et al. 1991)
  
• D type K⁺ (I_K⁺D) (Storm 1988)
  
• K type K⁺ (also known as "DR" type, I_K⁺K) (Klee et al. 1995)
  
• M type K⁺ (I_K⁺M) (Halliwell and Adams 1982)
  
• H type mixed cation (somatic and dendritic varieties, I_HSoma, I_HDend) (Magee 1998)
  
• L type Ca²⁺ (I_Ca²⁺L) (Christie et al. 1995; Fisher et al. 1990; Magee and Johnston 1995)
  
• N type Ca²⁺ (I_Ca²⁺N) (Christie et al. 1995; Fisher et al. 1990; Magee and Johnston 1995)
  
• R type Ca²⁺ (I_Ca²⁺R) (Magee and Johnston 1995)
  
• T type Ca²⁺ (I_Ca²⁺T) (Fisher et al. 1990; Magee and Johnston 1995)
  

To model these currents, we follow the style of Borg-Graham (1999). For complete details on the active ionic currents model see APPENDIXES B and C.

The sodium reversal potential E_NA⁺, was set to 70 mV. Ca²⁺ currents were modeled using a conductance-based formalism (not permeability) with E_Ca = 140 mV. Ca²⁺ accumulation, diffusion, and buffering were modeled using concentric shells (APPENDIX B). We used E_K = –140 mV to accurately account for perfusion of K⁺ from the intracellular electrode. This is based on using the Nernst equation where the extracellular concentration of K⁺ is 5 mM. Normally, the intracellular K⁺ concentration is 140 mM. However, the sharp pipette contains 1 M K⁺ and, given the fact that larger molecules like biocytin effectively perfuse from a pipette during current injection (which was applied throughout recording), the soma and apical shaft will most likely approach this level. Although the concentration will be close to 140 mM in oblique and distal dendrites as a result of active membrane processes restoring the K⁺ balance, it is the soma and apical shaft that make the largest contribution to the generation of observable EAPs. Consequently, we simplify the model by assuming a constant 1 M intracellular K⁺ concentration. Our experiences tuning the model, as described in RESULTS, led us to the conclusion that any inaccuracies resulting from this assumption are small compared with those of other sources of uncertainty in the biophysical parameters; the most likely consequence is slightly different tuned active current conductance densities compared with a model incorporating details of K⁺ perfusion and active K⁺ membrane pumps.

The exact kinetics and densities of the currents were tuned to match the simultaneous intra- and extracellular recordings taken in vivo while remaining faithful to the qualitative properties that were established by in vitro studies. This includes the nonuniform distribution of active ionic current conductances (as detailed in APPENDIX C). In general, the time constants were found to require values significantly faster than those measured in vitro and, when necessary, activation/inactivation curves were modified as well. The peak conductance densities associated with each current were treated as variables that needed to be tuned to match the simultaneous intracellular and extracellular recording, leaving the passive parameters and active ionic current kinetics fixed for all cells. The parameters tuned are those listed in Tables C1, C2, and C3 (see APPENDIX C).

Although Mainen et al. (1995) used a very high density of Na⁺ currents to cause axonal action potential initiation (1,000 x the somatic density), we used a model similar to that suggested by Colbert and Pan (2002): the initial segment and nodes of Ranvier have a Na⁺ current conductance density that is only moderately higher than the soma (1–3x) but the Na⁺ current in the initial segment and nodes of Ranvier (I_Na⁺Ax) has kinetics that activates at less depolarized potentials than the Na⁺ current in the soma and dendrites (I_Na⁺SD). This usually resulted in an action potential that initiates at the first node of Ranvier 0.5–1 ms before spreading to the soma by the initial segment and then to the dendrites. However, as described in Result of dendritic initiation (Fig. 7), our experiments suggest that there can be variability in the pattern of initiation with interesting consequences. For details on the kinetics and density for Na⁺ currents see APPENDIXES B and C.

View larger version (56K): [in this window] [in a new window]   FIG. 7. A: extracellular recording and simulation of neuron D180. Recording was made with a 6-site silicon probe with 25-µm spacing between individual taps. Sites closest to the cell body have a typical waveform, including the presence of the capacitive phase on the deepest recording site (site 1) but less so on the site with the strongest recording, presumably closest to the soma (site 2). In the shallower sites in the stratum oriens (4–6) the waveform is of smaller amplitude and the Na+ phase of the waveform has longer duration (see Electrode position and width of Na+ phase). Two recording sites, most distant from the cell (5–6), show unusual waveforms similar to a letter "W" in appearance. B: membrane potential and current for the axon and indicated compartments of dendrite 1. AP initiates in dendrite 1 and the axon before invading the soma. First negative peak in the W waveform results from the influx of current into dendrite 1 at initiation. C: membrane potential and current for the soma and indicated compartments of dendrite 2. AP invades the soma from both the axon and dendrite 1 and then spreads to dendrite 2. Current profile in the basal dendrite (first 2 rows, 2nd column) exhibits a large capacitive dominant phase that creates the "notch" in the W shape when superimposed on the negative potential generated by Na+ current in dendrite 1, the axon and the soma. Second negative peak of the W is the Na+-dominant phase in dendrite 2 (2nd row, 2nd column). D: extracellular potential waveforms in a transverse section 10 µm caudal to the soma and apical trunk (i.e., soma is at z = –10 µm, with respect to the plane of the EAPs). Assumed position of the recording electrode is indicated (dotted line). "W" waveform results when a recording site is close to the basal dendrite (<20 µm) where the AP initiates and also to a second to dendrite to which the AP backpropagates.

Synaptic input was mimicked by varying the leak resistance and reversal potential for the more distal compartments assumed to be receiving synaptic input. Typically, this meant reducing R_m by a factor of 3 to 5, and applying a reversal potential to the leak current of between –55 and –30 mV to mimic a mixed excitatory and inhibitory synaptic input.

Performance

Simulations were 25 ms long: the first 10 ms were used to establish the stability of the rest potential (between –70 and –60 mV) given by the combination of active currents without synaptic input. Simulated synaptic input was switched on to depolarize the cell, from 1 to 5 ms, until an action potential was triggered and the cell repolarized and returned to a stable rest potential as judged by remaining stable for the final 10 ms of the simulation. In addition a 4-ms extracellular spike trace at 10 locations was computed. All the above took 1–2 min, long enough that a brute-force automatic search for optimal parameters to match the simulations to the recordings was not feasible. Consequently, the primary method for tuning the channel density parameters of the model consisted of manual trials, supplemented by an automated local search of the parameter space (see Estimation of simulated EAP accuracy below).

	RESULTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
Membrane currents and the extracellular potential waveform

The simultaneous intra- and extracellular recording along with the model simulation matching the recording for cell D151 is shown in Fig. 1. As can be seen by comparing the total membrane current (I_total) to the EAP, the shape of the waveform is proportional to the time profile of I_total across the membrane of the perisomatic compartments. This is a straightforward consequence of the basic equation of the extracellular potential (Eq. 2) and the superposition principle. Note that, although the apical trunk compartment in the figure has a net current peak amplitude approximately 1/5 to 1/10 that of the soma, there are several other proximal dendritic compartments with virtually identical current versus time profiles. This gives the proximal dendrites about equal weight to the soma in determining the shape of the waveform. An additional consequence of Eq. 2, that will be further analyzed below, is that distal dendritic compartments make virtually no contribution to the EAP in the perisomatic region (i.e., the waveform detectable by an extracellular recording electrode). This is because the smaller diameter of the individual distal dendrites results in smaller net currents and the greater distance further reduces the impact.

View larger version (57K): [in this window] [in a new window]   FIG. 1. Recording and Simulation D151. A: extracellular action potentials (EAPs) in the transverse section containing the soma and the tip of the electrode track (dotted line). B: enlargement of the EAP at the estimated electrode position, and comparison to the recording (strongest channel of the tetrode). EAP is made up of 3 distinct phases: 1) a brief, positive peak; 2) a much larger negative peak; and 3) a positive period of longer duration and slowly decaying amplitude. C: comparison of the average intracellular recording with the simulated spike in the proximal apical trunk. Lack of pronounced afterhyperpolarization (AHP) suggests the intracellular electrode was not at the soma. D: details of the simulation in the indicated compartments. Shape of the EAP waveform is given by the shape of the net membrane current across the membrane at the soma and proximal dendrites (2nd column). Third column: makeup of the membrane current in terms of Na+, K+, and mixed-ion capacitive current. All 3 currents are simultaneously active throughout the action potential (AP); the 3 phases of the EAP correspond to the current that is dominant at that time: Brief positive peak at the start of the waveform is attributed to the positive capacitive current; the main negative peak is attributed to the influx of Na+ current driving the action potential; the final positive phase results from repolarizing K+ current flowing out of the cell. Last column: corresponding conductivity densities of the A, C, D, K, and M type K+ currents.

Another example of an EAP is illustrated in cell D068 (Fig. 2). As in most recordings, the EAP lacks a prominent capacitive current phase. The capacitive phase peak for cell D151 (Fig. 1) is 9% of the amplitude of the Na⁺ phase peak; in all 38 recordings in Henze et al. (2000) only four have a capacitive phase peak that is comparable or larger. We compared these results to the sample of 307 EAPs recorded without a simultaneous intracellular recording and found that 73% of the waveforms had a capacitive peak that was <5% of the amplitude of the Na⁺ peak, whereas 95% of the recordings had a capacitive peak amplitude that was <10% of the amplitude of the Na⁺ peak. The overall average ratio was 6%, and the distribution had a long tail, including a handful of recordings for which the ratio of the capacitive phase peak to the Na⁺ peak was close to 100% (recordings from Henze et al. 2000 contained two such waveforms).

View larger version (57K): [in this window] [in a new window]   FIG. 2. Recording and Simulation D068. A: EAPs in the transverse section containing the tip of electrode track (dotted line), about 5 µm caudal to the soma and apical trunk (the z-axis is the axis perpendicular to the plane of the section). Close to the soma, there is no initial, positive peak in the EAP. At locations along the apical trunk, the initial peak becomes pronounced. B: comparison of extracellular recording (strongest channel of the tetrode) and simulation at the estimated electrode position. There is only a slight hint of a positive deflection before the negative Na+-dominant phase of the waveform. C: comparison of intracellular recording and the simulation in the apical trunk approximately 120 µm from the soma. D: further simulation details. In the soma (middle row) the positive capacitive current, proportional to the change in membrane potential, coincides with the larger Na+ current (3rd column) because the membrane potential change is driven by local Na+ current. In dendritic compartments (4th and 5th rows) the membrane depolarization is initially driven by Na+ current from the soma, until local Na+ currents are activated and the action potential regenerates. In the brief time before the local Na+ currents activate, the positive capacitive current is the dominant membrane current and a capacitive-dominant phase is visible in the net current (2nd column). Consequently, the EAP has a more pronounced capacitive phase close to the apical trunk. Also shown, the Na+ current starts to enter the axon initial segment (2nd row) in a gradual way, driven by the action potential in the first node of Ranvier (not shown) and before the start of the main action potential at the soma. This results in a slightly negative slope in the extracellular potential before the start of the main EAP.

Cell D068 (Fig. 2) also illustrates the contribution of the axon to the EAP: because of the lower threshold of activation, the axon initial segment Na⁺ current activates before the soma Na⁺ current. This creates a slight depolarization at about 0.5 ms before the start of the more prominent three-phase EAP described above. However, most of the Na⁺ current in the initial segment is coincident with the somatic Na⁺ current and adds to the amplitude of the EAP peak. Although the AP typically initiates at the first node of Ranvier nearly 0.5 ms before the AP spreads to the soma, the myelinated axon and nodes of Ranvier make no contribution to the EAP. The myelinated axon lacks active channels and has only passive and capacitive currents that are very small relative to the active perisomatic currents. The nodes of Ranvier are too small and isolated to significantly contribute to the amplitude of the EAP: the peak Na⁺ current in the node is typically about 0.05 nA, two orders of magnitude less than the peak Na⁺ current in the soma.

Although the model does predict significant variation in the shape of the EAP waveform at distal locations, practically all of the significant variability occurs at amplitudes below the threshold of detectability. In Figs. 1–3 only the red and pink traces are above the typical detection and sorting threshold, about 60 µV, in CA1 (Henze et al. 2000). Aside from the aforementioned increase in the capacitive phase of the waveform along the apical trunk, in the region around the soma where an extracellular electrode would detect the spike above the background noise there is relatively little variation in the shape of the waveform, only an increase in amplitude as the electrode is moved closer to the soma.

View larger version (56K): [in this window] [in a new window]   FIG. 3. Recording and Simulation D112.1. A: EAPs in the transverse section containing the electrode track (dotted line), about 25 µm caudal to the soma and apical trunk (the z-axis is the axis perpendicular to the plane of the section). B: enlargement of the EAP at the estimated electrode position, and comparison to the recording (strongest channel of the tetrode). An initial positive peak is completely absent, and the Na+-dominant negative phase is made up of a narrow peak followed by a period of slower repolarization. K+-dominant phase is of minimal amplitude. C: comparison of intracellular recording and simulation at the soma. D: selected details of the simulation. Duration of significant Na+ channel activation is longer than the Na+-dominant phase seen in the EAP. Width of the Na+-dominant phase in the EAP is therefore determined by the timing of the K+ currents: strong A and C type currents (rightmost column). Very fast A current shapes the sharp peak of the Na+ phase by quickly counteracting Na+ current, although the A current rapidly inactivates. C current continues the repolarization in the somatic compartment, but has a much lower density in the dendritic compartments. This leads to the widening of the Na+ phase and the negligible amplitude of the K+ phase. Shape of the EAP (B) is similar to the current profile of the proximal dendrites (1st, 4th, and 5th rows, 2nd column).

The recording of cell D112.1 (Fig. 3) illustrates two other sources of variability seen in average EAPs of CA1 pyramidal cells: the width of the Na⁺-dominant phase may vary and show inflection points in the transition to the K⁺-dominant phase, and also the amplitude and shape of the K⁺-dominant phase may vary. The model suggests that this variability reflects differences in the conductance densities of active ionic currents, rather than electrode position.

Figure 3 illustrates how the width of the Na⁺-dominant phase of the EAP is influenced by the timing of the K⁺ currents. Variability in the timing of the net K⁺ current can occur because there are at least five different K⁺ currents (types A, C, D, K, and M) that make a contribution to the repolarization of the action potential and each type has unique properties with respect to the onset time and duration of the K⁺ current it generates. The A current is fast and inactivating; the C current is slower to activate and also inactivates slowly as a result of Ca²⁺ dependency; the D current is slow but has a low threshold of activation, being active at rest; the K current is fast but has a relatively high activation threshold; and the M current is the slowest but has a low threshold of activation, although not low enough to be active at rest. Although the C and K currents tend to dominate repolarization at the soma, the M and A currents make a larger contribution in the distal dendrites where the C current is not present (as described in APPENDIX C) and where the AP amplitude is too low to trigger significant K current.

Reproducing the simultaneous intra- and extracellular recordings with the model led us to the conclusion that it was necessary to assume significant variation in the conductance density levels of the K⁺ currents to achieve reasonable matches between recording and simulation. This was particularly the case when trying to match the width of the Na⁺ phase and the shape of the K⁺ phase in the extracellular waveform. Evidence for varying levels of expression of ionic channels in pyramidal cells has been found in the analysis of gene expression profiles (Toledo-Rodriguez et al. 2004). The peak conductance densities used to match the waveforms are given in Tables C1, C2, and C3 in APPENDIX C. Although the variability in the peak Na⁺ conductance density (as measured by a single SD) is on the order of only 25% of the mean, for the primary K⁺ currents the variability is of the same order of magnitude as the average.

For the model to match the narrow Na⁺-dominant phase (0.35 ms, measured as the width at 25% of peak amplitude) and rapidly decaying K⁺ phase in the EAP of cell D068 (Fig. 2) requires a dominant K type K⁺ current to rapidly counteract the Na⁺ current and achieve repolarization. Cell D151 (Fig. 1) has a somewhat longer-duration Na⁺ phase and a K⁺ phase that decays more slowly, and the best match with the model is achieved by assuming a balanced combination of C and K type K⁺ currents. Cell D112.1 (Fig. 3) is relatively unusual in that it has a sharply peaked Na⁺ phase in which the transition to the K⁺ phase slows significantly midway, leading to a relatively long duration Na⁺ phase that lasts 0.95 ms. Only three of 38 of the waveforms (Henze et al. 2000) had a similar shape and in a larger sample of EAP recordings (with no associated intracellular recording); the proportion having such a waveform was similar (27 of 307). This type of waveform was best matched by a repolarization consisting of primarily A and C type currents.

Although the amplitude of the waveform varies (generally getting smaller as a recording session progresses), there is little change in the waveform of the EAP during recording sessions that last up to 2 h. This suggests that the conductance density differences between cells are stable on this timescale.

Estimation of simulated EAP accuracy

The quality of fit values between the simulation and recording achieved by these methods are listed in Table 1. The measure of error used is the square root of the mean square error (SMSE) between the recording and the simulation. The SMSE is then normalized by the peak amplitude to enable comparison of the error on intra- versus extracellular recordings and on cells with EAPs of different amplitude. For comparison of the model with the extracellular recording, we also weighted the mean square error to emphasize an accurate match for the samples at the peak amplitude: the sample for the negative (Na⁺) peak was weighted 10 x that of a nonpeak sample, the points immediately before and after the peak were weighted 5x, and the points two samples distant from the peak were weighted 2.5 x (weights on all sample points were then normalized). Because of the longer duration of the intracellular action potential (IAP) peak, weighting was not necessary to ensure that the SMSE reflected a close match to the peak amplitude. For single-wire and tetrode recordings, we compare the simulated EAP with the largest-amplitude channel. For silicon probe recordings, we calculate the SMSE for all channels with peak amplitudes >20 µV and report the average.

The method used to tune the conductance densities of the active ionic current was a combination of manual tuning and automated search. Manual tuning of the parameters was directed by three principles: 1) the net Na⁺ current for the perisomatic region must produce the observed peak amplitude in the intracellular and extracellular APs; 2) K⁺ conductance densities in the perisomatic region are tuned to match the shape of the K⁺-dominant phase of the EAP waveform; 3) distal K⁺ currents are tuned so that the IAP in the simulation repolarizes at the rate given by the intracellular recording. The flow of significant axial currents in the cell (an order of magnitude greater than the membrane current in a typical dendritic compartment) means that M and A type K⁺ currents in the distal dendrites can result in significant repolarization to the soma and proximal dendrites while having virtually no contribution to the routinely recorded EAP.

After a reasonably good fit was achieved by manual exploration of the parameter space (i.e., normalized SMSE close to 5% for a recording with amplitude >20 µV), an additional local search in the parameter space was performed using the Nelder–Mead simplex algorithm. Improvements to the match were relatively modest, typically <1% SMSE for both the intracellular and extracellular waveforms. The parameter values tested by the search were typically within a few percent of the values at the initial starting point. When we attempted to tune the conductance density parameters solely with the automated search (starting from average values for the parameters) the resulting error was invariably inferior to that achieved by manual tuning. This suggests that the performance landscape in the parameter space is nonconvex, with many local optima.

Despite this obstacle and the high degree of variability in our data, we succeeded to model both intra- and extracellular data to within 10%. This could be achieved because the degree of nonconvexity arising from parameter variations at a fine scale is less significant than the relatively constrained choices available for a few dominant parameters. For a given recording, the choice of Na⁺ conductance density is well constrained by the requirement that the net Na⁺ current reproduces the peak amplitude of both the intra- and extracellular APs. Only the precise balance between somatic and dendritic densities leads to local optima in performance. At the same time, the choice of the dominant K⁺ conductance is constrained by the shape of the EAP: in all cases we examined cells whose K⁺ phase begins and ends rapidly were best matched by a dominant K type conductance, whereas cells with a K⁺ phase that begins and ends more slowly were best matched by a dominant C type conductance. Given that the timing of onset and offset, activation threshold, and density distributions of the K⁺ conductances are fixed to be consistent with in vitro studies (as described in APPENDIXES B and C) the model results suggest that the choice for which K⁺ conductance(s) dominate in a given cell are well constrained: it is fine-tuning the exact trade-off between the K⁺ currents that are less significant that creates many local optima. These observations lead us to conclude that the tuned parameters for conductance densities do, in fact, represent relatively unique combinations of dominant current types. Better results from the automated search may result from methods designed to escape from very small local optima, such as simulated annealing.

Although it is possible to match simultaneous intra- and extracellular recordings with a compartmental simulation, it is actually exceedingly difficult to simultaneously achieve low errors on both. Interestingly, we find that matching the EAP seems to be a tighter constraint on the model parameters than matching the IAP. Table 1 shows that lower errors were achieved for the match to the intracellular recording than the extracellular recording in nearly every case, despite equal weight being assigned to both errors in the manual and automated parameter tuning: the average error for all intracellular recordings is 5.5%, whereas the average error for all extracellular recordings is 9.8%. Figure 4 illustrates this point further by demonstrating that a model cell can have virtually identical membrane potentials despite significant changes in the distribution of active current densities. The extracellular waveform varies to a much greater degree when the distribution of active currents is changed. Consequently, it would be very unlikely to have significant changes in the IAP (peak amplitude, duration, or magnitude of AHP) attributed to conductance density changes without having a significant impact on the EAP. This results from the fact that the membrane potential in any given compartment combines membrane currents with axial currents, whereas the EAP results directly from the active membrane currents in the perisomatic region alone. This gives the membrane potential additional degrees of freedom compared with the extracellular potential.

View larger version (17K): [in this window] [in a new window]   FIG. 4. A: an EAP at (x = –20, y = 20, z = 0) for cell D151 simulated with 3 different distributions of active currents. EAP changes significantly as a result of the change in parameters: the EAPs each have normalized error [square root of the mean square error (SMSE)] of about 12.5% with each other. B: intracellular membrane potentials at the soma for the same simulations. Parameters have been tuned so that the membrane potentials are virtually identical, despite the very different EAPs: average SMSE is 3.5%. C: intracellular membrane potentials in the apical trunk at a distance of nearly 100 microns from the soma: the average SMSE is 5.5%. Although at this location the membrane potential reflects the difference in the active current conductance densities more than at the soma, the introduction of a second membrane potential measurement as a constraint does not provide as much of a measure of error as a single EAP measurement. Because of the high dependency of the EAP on the precise active current conductance density the reverse result—similar EAPs for different intracellular membrane potentials—would be exceedingly unlikely.

Earlier we presented evidence that the main impact of electrode position (other than on the amplitude of the EAP) was on the presence or absence of a positive capacitive-dominant peak in the waveform, and that the main determinant of the width of the Na⁺-dominant phase is the balance of ionic currents. However, the position of the extracellular electrode makes a secondary contribution to the duration of the Na⁺-dominant phase of the waveform, as shown in Fig. 5: the more distant the electrode from the cell, the longer the duration of the Na⁺-dominant phase. It follows that the duration of the Na⁺ phase of the EAP should be negatively correlated to the amplitude of the EAP. This relationship was analyzed in the sample of 307 EAP recordings made by silicon probes. The resulting correlation coefficient between the duration of the Na⁺ phase (as measured as the width at 25% peak amplitude) and the amplitude was –0.246 (P = 1.4 x 10^–5).

View larger version (20K): [in this window] [in a new window]   FIG. 5. A: duration of EAP Na+ phase vs. the simulated electrode position when moving away from the soma in a direction perpendicular to the apical trunk. Duration is measured as the width at 25% peak amplitude of the Na+-dominant (negative) phase. There is a nearly linear relationship between the duration of the Na+-dominant phase and the distance from the electrode to the soma. B: same data as in A, plotted vs. the amplitude of the waveform. C: membrane current vs. time at different locations illustrates the time required for the action potential to propagate back into the dendrite. Combined inward current phase of the soma and proximal dendrites makes a reasonably close prediction for the apparent duration of the Na+ phase observed at 45 µm, around 0.4 ms.

Cell morphology

A salient aspect of the reconstruction and simulations that had very little impact on the EAP waveform was the details of each cell's dendritic morphology. To demonstrate this relative indifference, we performed simulations using the reconstructed morphology of a CA1 basket cell (Gulyas et al. 1999) with the membrane parameters tuned to match one of our recordings from a pyramidal cell (Fig. 6). Although a pyramidal cell has a single apical dendrite and a skirt of basal dendrites extending from opposite ends of a conical soma, this basket cell has an approximately spherical soma with two thick dendrites and an assortment of narrow dendrites extending at more irregular angles. However, the resulting EAP very closely matches the original recording. We also observe that given the freedom to choose the electrode distance from the soma (to match the amplitude of the waveform), we can match virtually any EAP waveform using an arbitrarily chosen CA1 pyramidal cell for the simulation. This results from the fact that the observable waveform is determined by the net current profiles of the perisomatic compartments: if the distribution of active currents in the perisomatic region is similar for two cells, then the EAP waveform will be similar regardless of the precise location and branching pattern of the dendrites. These observations suggest that the major determinant of the EAP waveform depends primarily on the unique distribution of active current densities and considerably less on the morphological details of the cell.

View larger version (53K): [in this window] [in a new window]   FIG. 6. CA1 Basket cell PV08D (Gulyas et al. 1999) simulated with the conductance densities tuned to match the recording for the pyramidal cell D112.1 (Fig. 3). A: extracellular waveform is virtually the same as that for neuron D112.1, which suggests that the dendritic morphology makes little direct contribution to the waveform. B: EAPs around the simulated basket cell in the plane of the soma. EAPs below threshold are somewhat different from those for D112.1, but the above threshold EAPs are very similar.

Result of dendritic initiation

As mentioned earlier, the model predicts a variety of EAP waveforms in the subthreshold region of the basal dendrites, as in Fig. 7. Note that the pyramidal cell recorded in the session was not stained sufficiently for reconstruction, and so simulations were performed with a cell reconstructed from another session. As the simulation demonstrates, a "W"-shaped waveform can result from an AP in a single cell and does not require coincident firing in two nearby cells. The EAP is produced by superposition of the waveform resulting from an AP that initiates in a single dendrite and then propagates to the soma and then to a second dendrite close to the recording site. (The AP in the axon initiates before the soma but slightly after the first dendrite.) The two negative peaks in the waveform occur from two distinct Na⁺ influxes into two dendrites near the recording site. The dendritic initiation arises from a combination of strong synaptic input in the basal dendrite, below average Na⁺ density in the axon (but still greater than that in the soma), and a single dendrite having higher input resistance, as a result of smaller size, relative to the other dendrites. The amplitude of the Na⁺ phase peak in the "W" waveform is <10 µV and, although they appear clearly in an average of several hundred spikes, such waveforms would not be detected if it were not for the fact that the six-site probe also had recording sites close to the soma. Five of 13 silicon probe recordings from Henze et al. (2000) revealed similar waveforms on distal recordings sites. Recordings and modeling experiments prove that a large variety of waveforms can occur depending on the exact proximity to the local basal dendrite(s) and the exact pattern of initiation. In rare cases "W" waveforms occur at high amplitude. One recording in the data from Henze et al. (2000) was of a large-amplitude (about 40 µV negative peak amplitude) "W" waveform at a proximal site, which we have not been able to reproduce in the model.¹

	DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
Our main finding is that a detailed compartmental model combined with the LSA can accurately reproduce both the waveform and amplitude of recorded EAPs. Further, the model provides precise explanations for many of the features observed experimentally in EAPs. Despite the high degree of variability typical of biological measurements, the complexity of the system under investigation (a single, highly nonlinear nerve cell embedded into nonuniform neuronal tissue), the remoteness of distal dendritic events, and uncertainty about the exact electrode position, our model yields <10% average normalized error (Table 1). This reasonably close match validates the simplified assumptions we made to describe the extracellular milieu in our domain of interest (i.e., at frequencies <3 kHz) using ohmic but no capacitive components. It is known that action potentials in CA1 tend to fire in synchrony with local field potential (LFP) oscillations (Csicsvari et al. 1999). As a result, the EAP waveforms do in fact seem to contain an average LFP component of a few (typically <10) microvolts. These factors may significantly increase the error in matching low-amplitude EAPs.

Variability of conductance density

An important observation of the present experiments is that exact channel distributions, rather than the details of the dendritic arbors, are the main determinants of the EAP waveform. The implication of this finding is that CA1 pyramidal cells that are usually regarded as a homogeneous set of neurons do have substantial variability in their biophysical properties and that such variability may be identified with proper analysis of the extracellularly recorded spikes.

We demonstrated that for neurons to produce significantly different EAP waveforms in the perisomatic region, it is necessary for the net current profiles in the soma and proximal dendrites to be different, and that the source of this difference is most likely a varied balance of repolarizing K⁺ currents. We assumed that variability of the K⁺ current contributions made by different ionic current types on different cells is attributable solely to differences in conductance densities (i.e., ionic channel densities) and not to variability in channel behaviors (i.e., kinetics). However, it is known that the activation curves of some ionic channels may shift by several millivolts under various conditions in slice preparations (e.g., LTP; Frick et al. 2004). It seems possible that if we allowed the kinetics of ion current activation/inactivation to vary from cell to cell, we could match the simultaneous intra- and extracellular recordings with less variability in the conductance density levels.

In this light, we consider the choice to fix active current kinetics to be a pragmatic decision based on the need to manage the complexity of the model. In fact, several other important aspects of ionic current behavior were also effectively "lumped" into the conductance density. For example, Na⁺ channels and K and D type K⁺ channels in CA1 pyramidal cells are all known to undergo slow, activity-dependent inactivation (Colbert et al. 1997; Klee et al. 1995; Storm 1988). Because an average recording is modeled with a single simulated spike, what may in fact be slow inactivation of these channels is represented in the model as a reduced peak conductance density level. Similarly, the activity of the I_K⁺C current is strongly dependent on the intracellular [Ca²⁺], which may vary from cell to cell in real neurons but does not vary in the model. Consequently, what may in fact be variations in [Ca²⁺] are also reflected in the model as differences in the density of the I_K⁺C conductance. Finally, a variety of neuromodulators affect the activation and inactivation of ionic channels in pyramidal neurons (for review see, e.g., Tsubokawa 2000). Any impacts of these various factors are all folded into the parameter controlling the peak conductance density. Consequently, we consider conductance density as a single parameter that actually reflects a large number of factors contributing to the currently available number of any particular ionic channel and not only the absolute number of ionic channels present in the membrane.

Extracellular recording as a model constraint

We performed extensive and laborious manual and automated tuning of the peak conductance densities to minimize the error between the model and the recordings. However, these results do not prove that we have found a globally optimal solution for all of the parameters. In particular, we have not exhaustively searched the large number of parameters controlling the kinetics of ionic current activation and inactivation. Even without relaxing our assumption that active current kinetics are constant across different cells, it may be that an alternative formulation of the channel kinetics may provide a better match to the recording data.

Although the example of Fig. 4 is admittedly contrived, it suggests that EAPs have been an underutilized resource for constraining compartmental models. Combining extracellular recordings with intracellular recordings is a useful technique, even for in vitro physiology experiments, because it will add precision to the analysis of events involving active ionic currents. If the general properties of the IAP are known for a given cell class (i.e., typical amplitude and duration), then EAP measurements alone may provide sufficient information for constraining a compartmental model to match a specific experimental protocol. In turn, this feature may provide a useful technique for studying the properties of active ionic currents in vivo.

Impact of cell morphology on EAP

Our results predict that, although sizes of the soma and proximal dendrites constitute an important factor determining the amplitude of the EAP, the details of cell morphology make relatively little impact on the EAP waveform. However, a few caveats apply to this result: In the example of Fig. 6, not only the conductance densities but also the axonal and spine properties appropriate to a CA1 pyramidal cell were mapped onto a basket neuron. This is inexact because basket cells are spineless and have unmyelinated axons. A spineless cell with an unmyelinated axon would have significantly different input resistance and AP initiation properties compared with those of a pyramidal cell. To simulate such a cell the conductance density parameters that were tuned to match a pyramidal cell would need to be significantly modified to produce a stable action potential—let alone a particular EAP waveform—which is the reason we have used the pyramidal cell properties for this illustrative example. Therefore the example shows only that the shape of the soma and the location and size of the dendrites do not have an important impact on the EAP, but does not suggest that other properties of specific cell types will not be significant. Also, we have not ruled out the possibility that the cell morphology may be indirectly related to the shape of the EAP waveform by a feedback mechanism through which morphological details control the conductance density levels.

Pattern of AP initiation

Our results simulating W-shaped waveforms suggest that in specific circumstances the pattern of action potential initiation can be reflected in the EAP. In the simulation of Fig. 7 the dendritic AP initiation was caused by higher input resistance in one dendrite under conditions of uniform dendritic stimulation. However, it is plausible that, in vivo, different patterns of dendritic initiation may result from the precise pattern of synaptic input to different dendrites. Alternatively, if the active current conductance densities were significantly different in two dendrites, a similar result may follow (our model always assumed identical properties for all apical and all basal dendrites in a given cell). Consequently, monitoring EAPs may provide insights into the pattern of AP initiation in vivo. Initial experiments simulating the waveform of interneurons suggest that similar W-shaped EAPs may be generated when the axon extends from a proximal dendrite (Martina et al. 2000); in that situation the two negative peaks of the waveform are generated by the AP initiation in the axon that then spreads to the soma, but separated by propagation through a connecting dendrite. However, our failure to reproduce one W-shaped waveform in our recording set suggests that there may be other causes for the phenomena.

Expected developments

Because the model accurately predicts the amplitude of the EAP, we can use it to predict the distance at which an electrode would detect cells not only in the hippocampus but also in other brain regions, where sufficient data are available to constrain the model. Our LSA method is currently being applied to investigate the relative importance of different factors in determining at what distance EAPs may be detected by a recording electrode and to analyze how sampling bias arising from cell size may influence EAP recordings.

To date, extracellular unit recordings have been used mainly to monitor the spike output of neurons. However, multiple-site high-resolution extracellular recordings can also provide useful information about intracellular features of the recorded neurons. For example, the extent of somadendritic backpropagation of action potentials can be monitored in the intact brain as a function of behavior (e.g., Buzsáki et al. 1996; Quirk et al. 2001). Our findings suggest that accurate and high-resolution monitoring of EAP can provide information about alteration of conductance densities in single neurons as a function of state changes and plasticity. Because subcortical neuromodulators and use-dependent activation of neurons are known to modulate the various conductances, our results imply that by monitoring the waveforms, EAPs can provide access to this valuable information in behaving animals. Finally, models consisting of large numbers of geometrically arranged neurons can provide the needed tool for the perfection of unit-classifying algorithms and to the development of ideal recording electrode configurations without the need of simultaneous extracellular and intracellular recordings from various cell types and brain regions.

	APPENDIX

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION APPENDIX GRANTS ACKNOWLEDGMENTS REFERENCES

 
A. Passive parameters and spines

Spines are accounted for by adjusting the passive membrane parameters only for spiny compartments according to a folding factor f, similar to that described in Major et al. (1994). But here we modify the passive properties of the compartment directly without modifying the compartment length or diameter as that would also impact the properties of active ionic currents. That is

where A_s is the estimated spine surface area for the compartment and A_c is the actual surface are