|
|
||||||||
1Department of Mathematical Sciences and Technology and Center for Integrative Genetics, Norwegian University of Life Sciences, Ås, Norway; 2Athinoula A. Martinos Center, Massachusetts General Hospital, Charlestown, Massachusetts; 3Departments of Neurosciences and Radiology, University of California, San Diego, La Jolla, California; and 4Institute for Psychology of the Hungarian Academy of Sciences, Budapest, Hungary
Submitted ; accepted in final form
| ABSTRACT |
|---|
|
|
|---|
| INTRODUCTION |
|---|
|
|
|---|
Technology for large-scale electrical recording of neural activity is rapidly improving, enabling routine measurement of electric potentials from large numbers of closely spaced locations (Buzsaki 2004
; Nadasdy et al. 1998
). Local field potentials (LFP), i.e., the low-frequency band of extracellularly recorded potentials, are thought to predominantly stem from dendritic processing of synaptic inputs, but a direct interpretation in terms of the underlying neural activity is difficult (Nicholson and Freeman 1975
; Pettersen et al. 2006
). A standard measurement procedure has been to record the field potential at equidistant, linearly positioned electrode contacts using a laminar electrode vertically penetrating the cortical layers (Mitzdorf 1985
; Nadasdy et al. 1998
; Rappelsberger et al. 1981
; Ulbert et al. 2001
). The current source density (CSD) can be estimated from the second spatial derivative of the recorded LFP (Nicholson and Freeman 1975
; Rappelsberger et al. 1981
) or by using the newly developed inverse CSD (iCSD) method (Pettersen et al. 2006
). However, the CSD measures the net transmembrane current provided by all the underlying neural populations and the contributions from individual populations are difficult to assess. In a series of papers, Barth and collaborators pursued such an assessment by the application of principal component analysis (PCA) of the CSD (Barth and Di 1991
; Barth et al. 1989
, 1990
; Di et al. 1990
), and in their study of stimulus-evoked data from the rat barrel cortex, they identified putative cortical populations in supra- and infragranular layers of the barrel column (Di et al. 1990
).
PCA is one of several statistical methods where functions of two variables (here electrode contact positions and time) are expanded into sums over spatiotemporally separable functions, i.e., functions that can be written as a product of a spatial function and a temporal function. Such an expansion can be done in an infinite number of ways, and additional constraints are needed to make the expansion unique. In PCA, the functions (i.e., components) are constrained to be orthogonal both in space and time. In this paper, we instead use physiological constraints to specify the expansion. As a consequence, the expansion will not only be compatible with the physiology, but each component will also have a clear physiological interpretation.
The basic underlying constraint inherent in the laminar population analysis (LPA) presented here is that the observed LFP is evoked by the firing of action potentials in the modeled neuronal populations. An experimental measure of this population firing is obtained from the high-frequency component of the laminarly recorded extracellular potentials, i.e., the multiunit activity (MUA) (Ulbert et al. 2001
). The outcomes of the laminar population analysis are 1) identification of the relevant laminar cortical populations and their vertical spatial position and extent, 2) estimates of the firing-rates of these populations, and 3) estimates of the spatiotemporal LFP signature (and CSD signature) after action potential firing in the individual populations.
Here we apply LPA to stimulus-evoked laminar-electrode data from rat somatosensory (barrel) cortex, a well-studied example of topographic mapping where each of the large facial vibrissa (whiskers) is mapped onto a specific barrel column (Woolsey and Van der Loos 1970
). Several sample data sets, distinguished by stimulation paradigm, anesthesia conditions, and electrical boundary conditions at the cortical surface, are considered. Despite the different experimental situations, LPA seems to identify the same set of cortical populations (typically 1 supragranular, 1 granular, and 2 infragranular). Furthermore, the CSD signatures after population firing are generally found to be similar between experiments, in particular for the three topmost populations accounting for most of the observed LFP.
While these population CSDs give more information about the underlying neural circuit activity than traditional CSD analysis, they still do not provide a complete understanding of the synaptic connection pattern between the various populations. To estimate this pattern in this experimental situation, we use a new template-fitting analysis where LFP population templates are fitted to the observed LFP. These LFP population templates are calculated using the electrostatic forward solution with transmembrane currents obtained from compartmental modeling of morphologically reconstructed neurons. The results are found to agree with previous findings about the projection from the granular population onto the layer 2/3 supragranular population (Lübke et al. 2000
; Petersen and Sakmann 2001
) and provide predictions about other synaptic connections.
Preliminary results from this project were presented earlier in poster format (Einevoll et al. 2004
).
| METHODS |
|---|
|
|
|---|
All experimental procedures were approved by the Massachusetts General Hospital Subcommittee on Research Animal Care.
Male Sprague-Dawley rats (250350 g, Taconic) were used in the experiments. Glycopyrrolate (0.5 mg/kg, im) was administered 10 min before the initiation of anesthesia. Two anesthesia regimens were used. 1) Rats were anesthetized with 1.5% halothane in oxygen for surgery. After surgery, halothane was discontinued, and anesthesia was maintained with 50 mg/kg intravenous bolus of
-chloralose followed by continuous intravenous infusion at 40 mg/kg/h. 2) Rats were anesthetized with a single dose of ketamine hydrochloride (50 mg/kg, ip), supplemented with sodium pentobarbital (12 mg/kg) during the surgery. During the experiment, animals were supplemented by constant intraperitoneal infusion of 50 mg/kg/h of ketamine hydrochloride. During surgery, a tracheotomy was performed, and cannulas were inserted in the femoral artery and vein. All incisions were infiltrated with 2% lidocaine. After tracheotomy, rats were mechanically ventilated with 30% O2 in air. Ventilation parameters were adjusted to maintain PaCO2 between 35 and 45 mmHg, PaO2 between 140 and 180 mmHg, and pH between 7.35 and 7.45. Heart rate, ECG, and blood pressure were monitored continuously. Body temperature was maintained at 37.0 ± 0.5°C with a homeothermic blanket (Harvard Apparatus, Holliston, MA). The animal was fixed in a stereotaxic frame. An area of skull overlying the primary somatosensory cortex was exposed and thinned with a dental burr. The thinned skull was removed, and the dura matter was dissected to expose the cortical surface. A barrier of dental acrylic was built around the border of the exposure and filled with either saline or mineral oil. Mapping with a single metal microelectrode (FHC, 25 M
) was done to determine the positioning of the laminar electrode array. The optimal position was identified by listening to an audio monitor while stimulating different whiskers. After the mapping procedure, the electrode was withdrawn, and the array was slowly introduced at the same location perpendicular to the cortical surface. Contact 1 was positioned at the cortical surface using visual control. A linear array multielectrode with 23 contacts with diameter 0.04 mm spaced at 0.1 mm was used (Ulbert et al. 2001
), and the data from contacts 223 were used in the modeling.
Single whiskers were deflected upward by a wire loop coupled to a computer-controlled piezoelectric stimulator. We used a fast, randomized event-related stimulus presentation paradigm. The stimulus sequence was optimized for event-related response-estimation efficiency using the approach described by Dale (1999)
. The stimulation paradigm consisted of single deflections of varying amplitude with interstimulus interval (ISI) of 1 s.
We used three stimulus conditions. 1) In the nine-stimulus condition case, the vertical displacement of the whisker varied from 0.12 (amplitude 1) to 1.08 mm (amplitude 9). Intervening amplitudes had whisker displacements with equal increments on a linear scale. The stimulus rise time (time from onset to maximum displacement) was 30 ms, and the (mean) stimulus angular velocity increased from 76 (amplitude 1) to 660°/s (amplitude 9). 2) In the second stimulus condition, with altogether 27 different stimuli, three stimulus amplitudes [vertical displacements 0.40 (a1), 0.80 (a2), and 1.2 mm (a3)] each with nine stimulus rise times [20 (t1), 30 (t2), . ., 100 ms (t9)] were applied. The stimulus angular velocity varied between 76 (a1,t9) and 1,090°/s (a3,t1). 3) In the third stimulus condition, a stimulus rise time of 30 ms was used, with 27 different, linearly incrementing, vertical displacements
1.2 mm.
The recorded potential was amplified and filtered into two signals: a low-frequency part (0.1500 Hz, sampled at 2 kHz with 16 bits) and a high-frequency part (1505,000 Hz, sampled at 20 kHz with 12 bits; see Ulbert et al. 2001
for details). The low-frequency part is referred to as the LFP. The high-frequency part was further filtered digitally between 750 and 5,000 Hz using a zero phase-shift second-order Butterworth filter and rectified along the time axis to provide the MUA. Finally, the MUA data were decimated by a factor 10 along the time axis to provide the same time resolution (0.5 ms) as the LFP data. Note that the temporal form of the MUA was found to be essentially the same without decimation or when 500 or 1,000 Hz was used as the lower cut-off frequency in the band-pass filtering. Example traces of LFP and MUA (before rectification) for a single trial are shown in Fig. 1.
|
In this paper we consider four experiments. In experiment 1,
-chloralose was used as anesthesia, the stimulus-condition paradigm with three amplitudes and nine rise times was applied, and saline covered the exposed cortex around the laminar electrode. In experiment 2, ketamine was used as anesthesia, the nine stimulus-condition paradigm was applied, and oil covered the exposed cortex. In experiment 3,
-chloralose, saline, and the 27-amplitude stimulus condition was used, and in experiment 4,
-chloralose, saline and the three amplitude, nine rise-time stimulus condition was used. Experiment 1b corresponds to a separate 27-amplitude experiment performed on the rat in experiment 1. Eight of 5,400 trials were removed because of artifacts in the MUA. In experiments 24, the MUA signal from contact 9 was erratic, and for this channel, the arithmetic mean of the MUA signals from the neighboring contacts (contacts 8 and 10) was used instead.
The stimulus-averaged, baseline-subtracted LFP and MUA data used in the modeling are denoted as
L(zi,tj) and
M(zi,tj), respectively. Here zi refers to electrode contact position and tj to time. The time resolution is 0.5 ms.
LPA
Two-dimensional data like our LFP [
L(zi,tj)] and MUA [
M(zi,tj)] data can be expanded in a series of spatiotemporally separable functions
![]() | (1) |
In PCA, the expansion functions fn(zi) (called spatial loading) and gn(tj) (called temporal scores) are chosen so that the first component picks up most of the data variance, the second component most of the remaining variance, and so on. Furthermore, the PCA algorithm forces both the spatial loadings [fn(zi)] and the temporal scores [gn(tj)] to be orthogonal, i.e.,
, for n
m (Gershenfeld 1999
).
In LPA, the MUA and LFP data are modeled jointly. First, we assume the MUA data to be modeled as a sum over spatiotemporally separable contributions from several neuronal populations
![]() | (2) |
Because the MUA data are rectified, we impose the constraint that the spatial functions Mn(zi) should be nonnegative. The interpretation of rn(tj) as the population firing rate should imply that these functions should also be nonnegative. However, in our analysis of stimulus-averaged data, we subtract the baseline, i.e., mean MUA before stimulus onset. This means that rn(tj) should be interpreted as the population firing rate measured relative to baseline and should thus be allowed to have small negative values.
Next, we assume that the LFP data are driven by the population firing activities rn(tj) seen in the MUA; firing of an action potential of a neuron in population n will lead to postsynaptic transmembrane currents (including both the ligand-gated synaptic currents and consequent return currents). These transmembrane currents will in turn contribute to the LFP (Nicholson and Freeman 1975
). Consequently, we assume that the LFP data can be decomposed into contributions from each of the neuronal populations, and we suggest the following form
![]() | (3) |
rn)(tj) is the temporal convolution given by
![]() | (4) |
Mm(zi,tj) and
Lm(zi,tj) and the experimental data
M(zi,tj) and
L(zi,tj) as small as possible.
FORM OF SPATIAL POPULATION PROFILES Ln(zi) and Mn(zi).
The spatial profiles Ln(zi) are allowed to vary freely. However, for the MUA spatial profiles Mn(zi), we assume particular parameterized functional forms. The exact profiles will depend both on the morphology of the neurons in the population and the spatial distribution of the neurons constituting the population and are difficult to know in detail a priori. However, the results of Somogyvari et al. (2005)
for the vertical extent of the extracellular potential after an action potential (measured with the same type of laminar electrode) indicate that the MUA signal is mainly caused by sums of action potential contributions from neurons with their somas positioned within a distance 0.1 mm above or below the recording position. With the further assumption that each population is spatially continuous and covers a certain depth range in cortex, this implies that the MUA spatial profiles of the different populations will be largely nonoverlapping and that the individual spatial profiles will be spatially confined (bump-like) with rather sharp interface slopes and limited overlap to the neighboring populations. We thus model Mn(zi) as a trapezoidal symmetric profile centered at the position z0n. The profile has a constant height between z0n an/2 and z0n + an/2 and decays linearly to zero over a length bn on each side. We here constrain the flat maxima of the spatial profiles Mn(zi) to be nonoverlapping, (e.g., z0n + an/2 < z0N aN/2 if population N is positioned immediately below population n in cortex). We further constrain bn to be smaller than 0.1 mm to enforce a rather sharp interface between neighboring populations.
FORM OF TEMPORAL COUPLING KERNEL hn(tk).
This kernel can be thought of as the time-course of the effect of an action potential on the LFP, and one could in principle derive the functional form from knowledge of the electrophysiological effects of synaptic inputs onto neuronal dendrites. However, as described in the Results, our sample experiment 2 offers a way to estimate a suitable form directly from the experimental data. The resulting form is found to be
![]() | (5) |
(t) is the Heaviside unit step function [
(t
0) = 1,
(t < 0) = 0], and
n and
n are the time constant and a time-delay parameter, respectively, of population n. In this application, we make the additional simplifying assumption that all temporal kernels are identical, i.e., same
n and
n for all populations so that hn(tk) = h(tk).
NUMERICAL OPTIMIZATION PROCEDURE.
In the first step of the LPA procedure, we fit the MUA data to the expression in Eq. 2. Each spatial MUA profile Mn(z) is described by three parameters (z0n, an, bn). Thus in the numerical optimization scheme 3 x Npop parameters are varied to minimize the deviation eM between the experimental MUA data and the model expression. Here eM is the relative mean square error
![]() | (6) |
In the second step of the procedure, we fit the LFP data given the fitted firing rates rn(tj) from the MUA fit. The spatial LFP profiles Ln(zi) are fitted nonparametrically, so the parameter fit involves only the two model parameters of h(tj), the time constant
, and the delay
. The optimized solution minimizes the relative mean square error, defined in analogy to Eq. 6. The details of the numerical optimization procedure are given in APPENDIX 1.
CSD estimation
The recently developed iCSD method was used to estimate the population CSD Cn(zi) from the population LFP profiles Ln(zi) (Pettersen et al. 2006
). This method is based on the explicit inversion of the electrostatic solutions and can incorporate finite lateral extents of the underlying CSD and discontinuities in the extracellular conductivity at the cortical surface. Here we used the
-source iCSD method, where the CSD sources are assumed positioned at the electrode contacts and distributed in discs of infinitesimal thickness with lateral extent corresponding to the barrel diameter, combined with a three-point Hamming spatial filter (Pettersen et al. 2006
). For the case where the discs have infinite lateral extent and the extracellular conductivity is homogeneous, the
-source iCSD method was shown in Pettersen et al. (2006)
to correspond to the standard CSD method involving the double spatial derivative (Nicholson and Freeman 1975
).
Template-fitting analysis of LFP
More specific information than the CSD can be obtained by taking advantage of the fact that the forward electrostatic solutions, i.e., the LFP generated by synaptic inputs onto a particular dendritic structure (of a morphologically reconstructed neuron), can be evaluated: The transmembrane currents, i.e., CSD, after onset of synaptic input follow immediately from compartmental modeling of neurons, and given the extracellular conductivity, the LFP at any spatial position can be calculated (Gold et al. 2006
; Holt and Koch 1999
). This approach can be generalized to obtain corresponding LFP signatures for synaptic activation of a population of neurons. The observed LFP can be decomposed as a linear sum over such LFP population templates, and information about the synaptic connection pattern can be extracted from the fitted weights of the linear sum. Here we use this approach, denoted LFP template-fitting analysis, to obtain information about the synaptic connections between the laminar populations identified by LPA. For a detailed description, see APPENDIX 2.
| RESULTS |
|---|
|
|
|---|
In Fig. 2, we show the stimulus-averaged LFP and MUA data for experiment 1 (
-chloralose/saline/27 stimulus conditions) from the time of stimulus onset until 100 ms after onset. Results for 9 of the 27 stimulus conditions are shown: three different stimulus rise times (t9,t5,t1) for all three amplitudes (a1a3). The shortest rise time (t1) gives the largest responses, whereas the longest (t9) gives the smallest responses. From Fig. 2, we observe that the stimulus rise time has a stronger influence on the response than does the stimulus amplitude. The corresponding stimulus-averaged LFP and MUA data for experiment 2 (ketamine/oil/9 stimulus conditions) are shown in Fig. 3. Here, results for all nine amplitude stimulus conditions are shown (a1a9).
|
|
|
0.60.8 mm below the cortical surface. This MUA signal is expected to be dominated by neurons with a vertical displacement of less than
0.1 mm from the electrode contact position (Somogyvári et al. 2005Comparison of Figs. 2 and 3 also reveals that the maximum LFP and MUA signals are larger in experiment 1 than in experiment 2. This may reflect differences in electrode placement and type and depth of anesthesia. We further see that the MUA in experiment 1 has a longer duration and extends somewhat deeper down in cortex than in experiment 2.
Principal component analysis
Di et al. (1990)
used PCA in the interpretation of laminar-electrode LFP data from the barrel cortex of anesthetized rats using a similar stimulus paradigm. From stimulus-averaged LFP data, they derived an estimate of the CSD by performing a smoothed discrete spatial second derivative (Nicholson and Freeman 1975
). They found that two principal components (PCs) accounted for essentially all (98%) of the observed variance in CSD data obtained by averaging results from five animals. Although the direct interpretation of the PCs is difficult, PCA is a useful first analysis step for elucidating the effective dimensionality of the data, providing a hint at the model complexity needed and warranted in the further analysis.
In Fig. 4, we show the results from applying PCA on the stimulus-averaged LFP and MUA data in experiment 1 including all 27 stimulus conditions. For the LFP data, the first PC accounts for 90.6% of the variance and the second PC for 7.5%, and together they thus account for
98% of the data variance. For the MUA, the contributions from the first two PCs are 92.0% and 2.6%, respectively. To facilitate a comparison with the CSD analysis of Di et al. (1990)
, we also did PCA of CSD data obtained by using a smoothed double spatial derivative (cf. Ulbert et al. 2001
) on the LFP data. Here the two first PCs were found to account for 70.8% and 18.3% of the data variance, respectively. The spatial loadings and temporal scores of the first two PCs resulting from the application of PCA on the LFP, MUA, and CSD data are shown in Fig. 4.
In Fig. 5, we show the corresponding results from applying PCA on the data in experiment 2 including data for all nine stimulus conditions. Here the first two PCs account for 88.0% and 8.7% of the variance for the LFP data, 87.5% and 6.0% for the MUA data, and 89.7% and 7.2% for the CSD data. Comparison of the spatial loadings and temporal scores of the LFP and MUA data with the corresponding results for experiment 1 in Fig. 4 show many similar features. The main differences in the spatial loadings occur high up in the cortex, consistent with the different electrical boundary conditions at the cortical surface. The differences in temporal scores are somewhat larger. For example, the first PC for the LFP in experiment 1 has a shorter time-course than in experiment 2, whereas the situation is somewhat reversed for the MUA. For the CSD, the spatial loadings of the first PC have a similar shape for experiments 1 and 2, but the dominant sink extends deeper into cortex for experiment 1. For the second CSD PC, the difference is even larger, a difference also reflected in the difference in variance accounted for by this PC (18.3% for experiment 1 vs. 7.2% for experiment 2).
|
LPA
ESTIMATION OF TEMPORAL COUPLING KERNEL.
In LPA, the functional form of the temporal coupling kernel hn(tj) used in Eq. 3 must be specified. Data from experiment 2 offer a way to estimate it; close inspection of the MUA contour plots in Fig. 3B shows that the MUA response for the smallest stimulus (a1) is not only weaker and delayed compared with the higher-stimulus responses, it also seems to be spatially more confined. A detailed comparison with the corresponding LFP response in Fig. 3A further shows that the onset of the LFP response is delayed compared with the MUA onset. The spatial position of the MUA response (maximum response at
0.7 mm cortical depth) supports an interpretation that this response stems from layer 4 cells, and this suggests that the subsequent LFP is caused by synaptic inputs from these layer 4 cells. The absence of any significant MUA response after this LFP indicates that these synaptic inputs are not sufficiently strong to evoke significant firing of the postsynaptic cells for this particularly weak stimulus (Brecht et al. 2003
). The lack of any observed LFP before the onset of MUA in layer 4 indicates that the magnitude and spatial distribution of the thalamic synaptic input (as well as the morphological properties of the postsynaptic cells) are such that their contribution to the present LFP is small.
Our hypothesis is thus that for the stimulus condition a1 firing in the layer 4 population dominates the MUA and that this firing drives the subsequently observed LFP. The time-courses of these can be estimated by performing a PCA decomposition keeping only the first and thus dominant PCs. The temporal scores of these PCs, denoted gM1(tj) for the MUA and gL1(tj) for the LFP, can be used to obtain a direct estimate of the temporal coupling kernel h(tj). Figure 6 shows that the temporal coupling kernel for this situation is excellently modeled by the delayed exponential kernel in Eq. 5 with the parameters
= 16.4 ms and
= 0.6 ms. We also note that the peak in the temporal score of the LFP PC gL1(tj) occurs about 5 ms after the peak in corresponding MUA temporal score gM1(tj). These results seem to support the fundamental assumption of the model that the presently observed LFP predominantly results from MUA and that their coupling can be described using the simple delayed exponential temporal coupling kernel hn(tj) given in Eq. 5. We will therefore use this in the further analysis of the data from all stimuli.
|
|
|
= 13.4 ms and
= 2.0 ms.
The fitted LFP spatial profiles Ln(zi), n = 1, 2, 3, 4, are shown in Fig. 8B, and in Fig. 8D we show the time-courses of the postsynaptic LFP activity, (h
rn)(tj), for the four populations. These time-courses are obtained by convolving the population firing rates rn(tj) in Fig. 8C with the delayed, exponentially decaying kernel in Eq. 5 using the fitted parameter values for
and
. This convolution effectively corresponds to a low-pass filtering of the population firing rates, and comparison of (h
rn)(tj) in Fig. 8D with rn(tj) in Fig. 8C thus shows much smoother time-courses. The contribution to the LFP from each population n is determined by (h
rn)(tj) multiplied by the corresponding Ln(zi). From Fig. 8, B and D, we thus see that populations 2 (blue) and 3 (red) give large contributions to the LFP, whereas population 4 (green) generally provides a meager contribution.
In Fig. 9, C and A, we show results for the fitted MUA and LFP models, respectively, for experiment 2, assuming four cortical populations. The same procedure as for the model fitting to experiment 1 was used, and the fitting errors were found to be eM = 0.066 and eL = 0.087, respectively. The corresponding difference between the experimental data and the fitted model results are shown in Fig. 9, D and B. Overall, a very good agreement is observed. The corresponding fitted spatial profiles for the MUA and LFP [Mn(zi) and Ln(zi), n = 1, 2, 3, 4] are shown in Fig. 10, A and B, whereas the fitted population firing rates rn(tj) and the postsynaptic time-courses (h
rn)(tj) are shown in Fig. 10, C and D, respectively.
|
|
With the exception of population 4 (green), the predicted LFP spatial profiles Ln(zi) for the two experiments are seen to share most of the qualitative features (cf. Figs. 8B and 10B). The main difference occurs at the top where the influence on the shape of the LFP profiles from the different electrical boundary conditions at the cortical surface is strongest. The fitted model parameters for the temporal coupling kernel are
= 21.2 ms and
= 1.2 ms.
Population CSD
The interpretation of the predicted LFP spatial profiles Ln(zi) (shown in Figs. 8B and 10B) is less obvious than for the MUA spatial profiles Mn(zi). These LFP profiles represent the spatial form of the LFP after firing of action potentials of neurons in population n and will consist of a sum over contributions to the LFP from the CSD generated by the synaptic projections from population n onto all other populations. By analogy with previous CSD analyses, we can thus separately estimate the CSD resulting from activity in each population, i.e., the population CSD denoted Cn(zi). This measure differs from the traditional estimation of CSD that gives the total CSD across all populations.
The estimation of the population CSD Cn(zi) from the LFP profiles Ln(zi) raises the same methodological issues as the estimation of the total CSD from laminar LFP recordings. In general, the LFP stemming from a particular CSD distribution will depend on both the spatial extent of the CSD and the extracellular conductivity (Pettersen et al. 2006
). Here we use the iCSD method (Pettersen et al. 2006
) to incorporate these effects in the CSD estimation. The effective lateral sizes of the CSDs [Cn(zi)] underlying the LFP profiles [Ln(zi)] are not known, however. We therefore consider two different situations: 1) infinite lateral extent of the CSD (as assumed in the standard CSD method, see Pettersen et al. 2006
) and 2) cylindrical CSD distribution of diameter 0.5 mm centered at the axis of the laminar electrode, with constant CSD in the horizontal direction. In addition to the common choice of assuming a constant electrical conductivity
, we also consider other electrical boundary conditions at the cortical surface. For experiment 1 with highly conductive saline added at the cortical surface, we also evaluate CSD estimates for the situation where
/
0 = 0, i.e., the conductivity
0 above the cortical surface (z < 0) is infinitely large. For experiment 2 with poorly conducting oil added at the cortical surface, we also evaluate CSD estimates for the situation where
0 = 0, i.e., the conductivity above the cortical surface is zero. Although none of these situations can be completely realistic, qualitative features of the population CSD [Cn(zi)] present in all sets of estimates are more likely to be genuine than features present only for a particular choice of CSD diameter and/or electrical boundary condition.
In Fig. 11, we show the resulting estimates for the population CSDs for experiments 1 and 2. A first feature to note is that effects of different electrical boundary conditions in the estimation procedure are mainly seen at the electrode contacts positioned 0.1 and 0.2 mm below the cortical surface. The effects of lateral spatial confinement of the CSD sources are seen at all electrode positions, however.
|
0.50.6 mm depth, with a dominant source above. For population 3 (red), this pattern is reversed. In addition, population 3 seems to generate a sink at
1 mm. For population 1 (black), both experiments predict a dominant sink at
0.3 mm, but the source below is predicted to extend deeper into cortex in experiment 2. Furthermore, in experiment 2, a source is predicted at the top electrode position (0.1 mm), whereas for experiment 1, such a source is only seen if a constant conductivity is assumed (solid black). The CSD profiles for the deepest population, population 4 (green), are seen to be quite different in the two experiments, in particular in the upper 1 mm of cortex, but this population also accounts for a rather small part of the total LFP, especially for experiment 2. Most of the same qualitative features of the population CSDs are observed when a population diameter of d = 0.5 mm is assumed. Even though the sinks and sources generally are predicted to be somewhat wider and smoother compared with the infinite-diameter assumption, the main predictions described above seem to be robust.
In Fig. 12 B, we show the population CSD profiles, assuming infinite population diameters and homogeneous electrical conductivity, for the three topmost populations for all our five data sets from four rats. As observed, the salient features of these population CSDs seem to be quite robust across the experiments, in particular for the two topmost populations. We note, however, that for two of the experiments (3 and 4), the population CSD profiles seem to have been shifted toward lower electrode contacts. Because the same is observed for the MUA profiles in Fig. 12A, it seems that the laminar-electrode position was shifted in the vertical direction (relative to the cortical layers) compared with the other experiments.
|
The population LFPs [Ln(z); and CSDs, Cn(z)] estimated in the LPA procedure will reflect the LFP (CSD) generated by synaptic inputs from a particular presynaptic population, but it does not specify the postsynaptic neuronal populations receiving these synaptic projections. Such information can be obtained by taking advantage of the fact that the forward electrostatic solution, i.e., the LFP generated by synaptic activation onto a population of similar neurons, can be evaluated. We here apply an LFP template-fitting analysis (see METHODS and APPENDIX 2) to estimate the synaptic connection pattern between the neuronal populations for our two sample data sets, experiments 1 and 2.
FORWARD CALCULATION OF CSD PROFILES.
In the left columns of Fig. 13, A and B, we show the steady-state CSD profiles of a layer 5 and a layer 3 cortical neuron (morphology taken from Mainen and Sejnowski 1996
) after injection of a sink-like steady-state current (corresponding to excitatory inputs) into the dendritic trees at different vertical positions measured from the soma. Current has been injected into all dendritic branches crossing a horizontal plane at the particular vertical position, and the total injected current reflects the number and circumference of the intersections with the dendritic branches. From Fig. 13, we see that for the layer 5 cell, excitatory input onto the lower basal dendrites gives a deep sink accompanied by a spatially more distributed source from the return current around the soma and at the top of the apical dendrite. Conversely, excitatory input at the apical dendrite gives a high sink accompanied by a spatially distributed return-current source with a notable contribution from the soma area. A similar, but spatially more compressed pattern is seen in the CSD profile for the layer 3 neuron.
|
ESTIMATION OF LFP POPULATION TEMPLATES. The simple but characteristic dipolar LFP patterns observed in Fig. 13, A and B, for synaptic activation for these two example pyramidal neurons suggests that the population LFPs [Ln(z)] can be expanded as a weighted sum over such LFP patterns characteristic for the neuronal populations in barrel cortex. Furthermore, the standardized dipolar pattern forms suggest that the LFP from such synaptic activation can be well described by the first PC (found in a singular-value decomposition) of the LFP patterns in Fig. 13.
As detailed in the Discussion, the LPA results suggests the following identification of populations: population 1 (black in Fig. 8A) corresponds to a population of layer 2/3 pyramidal cells, population 2 (blue) corresponds to a population of predominantly stellate layer 4 cells, population 3 (red) corresponds to a population of layer 5 pyramidal cells, and population 4 (green) correspond to a population of layer 6 (or deep layer 5) pyramidal cells. Contributions from projections onto the (predominantly) stellate cells (population 2) are assumed negligible because of the closed-field structure of the stellate cells in this layer. As seen in Figs. 2 and 3, the electrical boundary conditions significantly affect the LFP at the topmost electrode and is thus included in the forward LFP modeling.
As described in APPENDIX 2, the characteristic spatial LFP patterns, called LFP templates, for synaptic activation of each of the three pyramidal populations (L2/3, L5, and L6) were obtained in a numerical optimization procedure identifying the population diameters and electrical conductivity-jump parameter Wim that gives the templates providing the best fit to the total LFP. The soma positions zn0 were kept fixed at the values given by the MUA fit, whereas the population diameters and the electrical conductivity jump at the cortical surface were allowed to vary freely. The LFP profiles resulting from minimizing the value of eL for the two experiments are shown in Fig. 13C. For experiment 1, the temporal sum over the three optimized LFP templates was found to account for 98.6% of the LFP data, whereas for experiment 2, 96.8% of the LFP data were accounted for. Thus our LFP templates are able to account for most of the LFP data.
TEMPLATE-FITTING ANALYSIS. As described in APPENDIX 2 the population LFPs [Ln(z)] can now be projected onto the three fitted LFP templates, produced by activation of the different cell populations: the magnitudes of the fitted weight coefficients measure the proportion of the population LFPs caused by projections from particular presynaptic and onto particular postsynaptic populations, whereas the sign of the fitted weight coefficients tells about the spatial positioning of the synaptic projections (i.e., basal or apical). For example, for both sample experiments, the layer 4 population is found to have a significant positive weight of the template corresponding to the layer 2/3 population (seen in Fig. 13C). This means that the layer 4 population is predicted to have a significant excitatory synaptic projection onto the lower basal dendrites of the layer 2/3 population. Likewise, the analysis predicts the layer 4 population to have a weaker excitatory synaptic projection onto the basal dendrites of the layer 5 population (alternatively, an inhibitory projection onto the apical dendrites of the same population). If we assume the dominant synaptic projections to be excitatory, some further predictions from the analysis of our two experiments are as follows: the layer 2/3 population projects to its own basal dendrites, to the apical dendrites of the layer 5 population, and (weakly) to the apical dendrites of the layer 6 population. The layer 5 population projects to its own basal dendrites and to the apical dendrites of the layer 2/3 population. The predicted projections from the layer 6 population are, however, found be conflicting in the two experiments.
| DISCUSSION |
|---|
|
|
|---|
Inspection of the predicted spatial and temporal MUA profiles (Figs. 8 and 10) and population CSDs in Fig. 11 suggests the following interpretation: population 1 (black) positioned above
0.5 mm cortical depth corresponds to a population of layer 2/3 pyramidal cells, population 2 (blue) centered at
0.7 mm corresponds to a population of predominantly stellate layer 4 cells, population 3 (red) centered at
1.11.2 mm corresponds to a population of layer 5 pyramidal cells, and population 4 (green) positioned below
1.5 mm corresponds to a deep infragranular population of layer 6 (or deep layer 5) pyramidal cells. In the following, we will for clarity use this suggested identification when we denote the populations. This population identification is supported by several observations:
CORTICAL DEPTH POSITIONS OF ESTIMATED POPULATION BORDERS.
The depth positions extracted by LPA are compatible with previously suggested borders between cortical layers in barrel cortex of similarly sized rats (Armstrong-James et al. 1992
; Di et al. 1990
). In Armstrong-James et al. (1992)
, the border between layers 3 and 4 were estimated to be around 0.5 mm, the border between layers 4 and 5 around 0.750.8 mm, and the border between layers 5 and 6 at 1.5 mm.
TEMPORAL ORDERING OF POPULATION FIRING.
Armstrong-James et al. (1992)
measured the response latency after single whisker flicks. In the principal barrel, they found the following approximate sequence of first firing of action potentials: 1) layer 4, 2) layer 5b, 3) layer 3 (slightly later), 4) layer 5a and 2, and 5) layer 6. The temporal resolution of our predicted firing rates rn(t) is lower, but for experiment 1 (Fig. 8B), it is clearly seen that the layer 4 population initiates firing before the layer 2