Abstract
Full text Figures and data Side by side Abstract Editor’s evaluation Introduction Results Discussion Materials and methods Appendix 1 Data availability References Decision letter Author response Article and author information Metrics Abstract Local field potential (LFP) recordings reflect the dynamics of the current source density (CSD) in brain tissue. The synaptic, cellular, and circuit contributions to current sinks and sources are ill-understood. We investigated these in mouse primary visual cortex using public Neuropixels recordings and a detailed circuit model based on simulating the Hodgkin–Huxley dynamics of >50,000 neurons belonging to 17 cell types. The model simultaneously captured spiking and CSD responses and demonstrated a two-way dissociation: firing rates are altered with minor effects on the CSD pattern by adjusting synaptic weights, and CSD is altered with minor effects on firing rates by adjusting synaptic placement on the dendrites. We describe how thalamocortical inputs and recurrent connections sculpt specific sinks and sources early in the visual response, whereas cortical feedback crucially alters them in later stages. These results establish quantitative links between macroscopic brain measurements (LFP/CSD) and microscopic biophysics-based understanding of neuron dynamics and show that CSD analysis provides powerful constraints for modeling beyond those from considering spikes. Editor’s evaluation The study demonstrates that utilizing the LFP and/or the CSD in modeling can facilitate model configuration and implementation by revealing discrepancies between models and experiments. The analysis of the biophysical origin of the canonical CSD using the model is an interesting and worthy line of investigation. The dissection of CSD components is detailed and exhaustive. A key novelty of this article is the addition of CSD patterns as another constraint to more accurately infer the model parameters beyond its prior state. https://doi.org/10.7554/eLife.87169.sa0 Decision letter eLife’s review process Introduction The local field potential (LFP) is the low-frequency component (below a few hundred Hertz) of the extracellular potential recorded in brain tissue that originates from transmembrane currents in the vicinity of the recording electrode (Lindén et al., 2011; Buzsáki et al., 2012; Einevoll et al., 2013; Pesaran et al., 2018; Sinha and Narayanan, 2022). While the high-frequency component of the extracellular potential, the single- or multi-unit activity (MUA), primarily reflects action potentials of one or more nearby neurons, the LFP predominantly stems from currents caused by synaptic inputs (Mitzdorf, 1985; Einevoll et al., 2007) and their associated return currents through the membranes. Thus, cortical LFPs represent aspects of neural activity that are complementary to those reflected in spikes, and as such, they can provide additional information about the underlying circuit dynamics from extracellular recordings. Applications of LFP are diverse and include investigations of sensory processing (Rall and Shepherd, 1968; Di et al., 1990; Victor et al., 1994; Kandel and Buzsáki, 1997; Mehta et al., 2000a; Mehta et al., 2000b; Henrie and Shapley, 2005; Einevoll et al., 2007; Belitski et al., 2008; Montemurro et al., 2008; Niell and Stryker, 2008; Nauhaus et al., 2009; Bastos et al., 2015; Senzai et al., 2019), motor planning (Scherberger et al., 2005; Roux et al., 2006), navigation (Tort et al., 2008; Makarova et al., 2011; Fernández-Ruiz et al., 2012; Watrous et al., 2013; Fernández-Ruiz et al., 2017), and higher cognitive processes (Pesaran et al., 2002; Womelsdorf et al., 2006; Liu and Newsome, 2006; Kreiman et al., 2006; Liebe et al., 2012). The LFP is also a promising signal for steering neuroprosthetic devices (Mehring et al., 2003; Andersen et al., 2004; Rickert et al., 2005; Markowitz et al., 2011; Stavisky et al., 2015) and for monitoring neural activity in human recordings (Mukamel and Fried, 2012) because the LFP is more easily and stably recorded in chronic settings than spikes. Due to the vast number of neurons and multiple neural processes contributing to the LFP, however, it can be challenging to interpret (Buzsáki et al., 2012; Einevoll et al., 2013; Hagen et al., 2016). While we have extensive phenomenological understanding of the LFP, less is known about how different cell and synapse types and connection patterns contribute to the LFP or how these contributions are sculpted by different information processing streams (e.g., feedforward vs. feedback) and brain states. One way to improve its interpretability is to calculate the current source density (CSD) from the LFP, which is a more localized measure of activity, and easier to read in terms of the underlying neural processes. The current sinks and sources indicate where positive ions flow into and out of cells, respectively, and are constrained by Kirchoff’s current law (i.e., currents sum to zero over the total membrane area of a neuron). However, the interpretation of current sinks and sources is inherently ambiguous as several processes can be the origin of a current sink or source (Buzsáki, 2006; Pettersen et al., 2006; Einevoll et al., 2007). For example, a current source may reflect an inhibitory synaptic current or an outflowing return current resulting from excitatory synaptic input elsewhere on the neuron. There is no simple way of knowing which it is from an extracellular recording alone (Buzsáki, 2006). Another approach to uncovering the biophysical origins of current sinks and sources, and by extension the LFP, is to simulate them computationally (Pettersen et al., 2008; Einevoll et al., 2013). Following the classic work by Rall in the 1960s (Rall, 1962), a forward-modeling scheme in which extracellular potentials are calculated from neuron models with detailed morphologies using volume conduction theory under the line source approximation has been established (Holt and Koch, 1999). With this framework, we have achieved a good understanding of the biophysical origins of extracellular action potentials (Koch, 1998; Holt and Koch, 1999; Pettersen and Einevoll, 2008; Hay et al., 2011; Lindén et al., 2010). Expanding on this understanding, models composed of populations of unconnected neurons (e.g., Pettersen et al., 2008; Lindén et al., 2011; Schomburg et al., 2012; Łęski et al., 2013; Sinha and Narayanan, 2015; Hagen et al., 2017; Ness et al., 2018) and recurrent network models (e.g., Traub et al., 2005; Vierling-Claassen et al., 2010; Reimann et al., 2013; Głąbska et al., 2014; Tomsett et al., 2015; Hagen et al., 2016; Hagen et al., 2018; Chatzikalymniou and Skinner, 2018) have been used to study the neural processes underlying LFP. While interesting insights about CSD and LFP were obtained from these computational approaches, establishing a direct relationship between the biological details of the circuit structure and the electrical signal like LFP remains a major unresolved challenge. One reason is that the amount and quality of data available for modeling the circuit architecture in detail have been limited. This situation improved substantially in recent years, and a broad range of data on the composition, connectivity, and physiology of cortical circuits have been integrated systematically (Billeh et al., 2020) in a biophysically detailed model of mouse primary visual cortex (area V1). In addition, significant improvements were achieved in experimental recordings of the LFP and the simultaneous spiking responses. In particular, the Neuropixels probes (Jun et al., 2017) record LFP and hundreds of units across the cortical depth in multiple areas, with 20 μm spacing between recording channels allowing for an unprecedented level of spatial detail. These developments provide unique opportunities to improve our understanding of circuit mechanisms that determine LFP patterns. Here, we analyze spikes and LFP from the publicly available Allen Institute’s Visual Coding survey recorded using Neuropixels probes (https://www.brain-map.org; Siegle et al., 2021) and reproduce these using the mouse V1 model developed by Billeh et al., 2020. The model is comprised of >50,000 biophysically detailed neuron models surrounded by an annulus of almost 180,000 generalized leaky-integrate-and-fire units. The neuron models belong to 17 different cell type classes: one inhibitory class (Htr3a) in layer 1, and four classes in each of the other layers (2/3, 4, 5, and 6) where one is excitatory and three are inhibitory (Pvalb, Sst, Htr3a) in each layer. The visual coding dataset consists of simultaneous recordings from six Neuropixels 1.0 probes across a range of cortical and subcortical structures in 58 mice while they are exposed to a range of visual stimuli (about 100,000 units and 2 billion spikes over 2 hr of recording). In our analysis of this dataset, we identified a canonical CSD pattern that captures the evoked response in mouse V1 to a full-field flash. We then modified the biophysically detailed model of mouse V1 to reproduce both the canonical CSD pattern and laminar population firing rates in V1 simultaneously. We reproduce, in a quantitative manner, the shape and timing of the pattern of current sources and sinks that have been described in considerable detail by experimentalists (e.g., Mitzdorf, 1987; Swadlow et al., 2002; Senzai et al., 2019). This shows that adjustments to synaptic parameters such as weights and placement in addition to a circuit architecture that included feedback are sufficient to reproduce experimental findings on both single-cell measures such as spikes and population-level measures such as CSD. We use this model to explain, in a highly mechanistic manner, the biophysical origins of the various ionic current sinks and sources and their location across the various layers of visual cortex. In the process of obtaining a model that could reproduce both spikes and CSD, we discovered that the model can be modified by adjusting the synaptic weights to reproduce the experimental firing rates with only minor effects on the simulated CSD, and, conversely, that the simulated CSD can be altered with only minor effects on the firing rates by adjusting synaptic placement. Furthermore, we found that comparing the simulated CSD to the experimental CSD revealed discrepancies between model and data that were not apparent from only comparing the firing rates. Additionally, it was not until feedback from higher cortical visual areas (HVAs) was added to the model that simulations reproduced both the experimentally recorded CSD and firing rates, as opposed to only the firing rates. This bio-realistic modeling approach sheds light on specific components of the V1 circuit that contribute to the generation of the major sinks and sources of the CSD in response to abrupt visual stimulation. Our findings demonstrate that utilizing the LFP and/or the CSD in modeling can aid model configuration and implementation by revealing discrepancies between models and experiments and provide additional constraints on model parameters beyond those offered by the spiking activity. The new model obtained here is freely accessible (https://doi.org/10.5061/dryad.k3j9kd5b8) to the community to facilitate further applications of biologically detailed modeling. Results Spikes and LFP were recorded across multiple brain areas, with a focus on six cortical (V1, LM, AL, RL, AM, PM) and two thalamic (LGN, LP) visual areas, using Neuropixels probes in 58 mice (Siegle et al., 2021). A schematic of the six probes used to perform the recordings in individual mice is shown in Figure 1A, and the spikes and LFP recorded in V1 of an exemplar mouse during presentation of a full-field bright flash stimulus are displayed in Figure 1B, C. The CSD can be estimated from the LFP (averaged over 75 trials) using the delta iCSD method to obtain a more localized measure of inflowing (sinks) and outflowing currents (sources) (Pettersen et al., 2006; Einevoll et al., 2013). The biophysically detailed model of mouse V1 used to simulate the neural activity and the recorded potential in response to the full-field flash stimulus is illustrated in Figure 1E. The model contains 230,924 neurons, of which 51,978 are biophysically detailed multicompartment neurons with somatic Hodgkin–Huxley conductances and passive dendrites, and 178,946 are leaky-integrate-and-fire (LIF) neurons. These neuron models are arranged in a cylinder with a radius 845 μm and a height 860 μm. The multicompartment neurons are placed in the ‘core’ with a radius of 400 μm, while the LIF neurons form an annulus surrounding this core. Cellular models belong to 17 different classes: one excitatory class and three inhibitory (Pvalb, Sst, Htr3a) in each of layers 2/3, 4, 5 and 6, and a single Htr3a inhibitory class in layer 1. The extracellular electric field in the model was recorded on an array of simulated point electrodes (Dai et al., 2020) arranged in a straight line (Figure 1D) and separated by 20 μm, consistent with Neuropixels probes, shown in Figure 1E to scale with the model. Figure 1 Download asset Open asset Illustration of experimental data and the biophysical model for mouse primary visual cortex (V1). (A) Schematic of the experimental setup, with six Neuropixels probes inserted into six cortical (V1, latero-medial [LM], rostro-lateral [RL], antero-lateral [AL], postero-medial [PM], AM) and two thalamic areas (LGN, LP). (B) Top: spikes from many simultaneously recorded neurons in V1 during a single trial. Bottom: spikes from a single neuron recorded across multiple trials. In both cases, the stimulus was a full-field bright flash (onset at time 0, offset at 250 ms). (C) Top: local field potential (LFP) across all layers of V1 in response to the full-field bright flash, averaged over 75 trials in a single animal. Bottom: current source density (CSD) computed from the LFP with the delta iCSD method. (D) Histology displaying trace of the Neuropixels probe across layers in V1, subiculum (SUB) and dentate gyrus (DG). (E) Visualization of the V1 model with the Neuropixels probe in situ. (Image made using VND.) Uncovering a canonical visually evoked CSD response We first established a ‘typical’ experimentally recorded CSD pattern to be reproduced with the model. Though there is substantial inter-trial and inter-animal variability in the evoked CSD response, we find that most trials and animals have several salient features in common. In Figure 2A, the trial-averaged evoked CSDs from five individual mice are displayed. In the first four animals (#1–4), we observe an early transient sink arising in layer 4 (L4) ~40 ms after flash onset, followed by a sustained source starting ~60 ms, which covers L4 and parts of layers 2/3 (L2/3) and layer 5 (L5). We also observe a sustained sink covering layers 5 and 6 (L6) emerging around 50 ms, as well as a sustained sink covering layers 1 and 2/3 around 60 ms. An animal that does not fully exhibit what we term the ‘canonical’ pattern is shown in the rightmost plot (#5 in Figure 2A); it has an early L4 sink arising at 40 ms, but this sink is not followed by the sustained sinks and sources from 50 to 60 ms and onward observed in the other animals. The timing and location of sinks and sources are, overall, similar to those described earlier by Givre et al., 1994; Schroeder et al., 1998, Niell and Stryker, 2008, and Senzai et al., 2019. Figure 2 with 3 supplements see all Download asset Open asset Variability in experimentally recorded current source density (CSD). (A) Evoked CSD response to a full-field flash averaged over 75 trials, from five animals in the dataset. (B) The first principal component (PC) computed from the CSD of all n = 44 animals, explaining 50.4% of the variance. (C) Illustration of movement of sinks and sources in the calculation of the Wasserstein distance (WD) between the CSD of two animals in the dataset. The gray lines in the rightmost panels display how the sinks or sources of one animal are moved to match the distribution of sinks or sources of the other animal. (D) Left: WDs from each animal to the PC 1 CSD. Right: pairwise WDs between all 44 animals sorted by their distance to the first PC. (E) CSD from five individual trials in example animal 1. (F) Distribution of pairwise distances between single-trial CSD (red) and pairwise distances between trial-averaged CSD of individual animals (blue). Both are normalized to the maximum pairwise distance between the trial-averaged CSD of individual animals. (G) Pairwise WDs between trials in each of 44 animals (white boxplots), normalized to maximal pairwise WDs between trial-averaged CSD of animals. Gray-colored boxplot shows the distribution of pairwise WDs between trial-averaged CSD of individual animals, and the red stars indicate the n = 5 animals for which the inter-trial variability was greater than the inter-animal variability (assessed with Kolmogorov–Smirnov [KS] tests, p < 0.001 in all cases, see Figure 2—figure supplement 3). To identify the robust features across animals in this dataset, we performed principal component analysis (PCA) on the trial-averaged evoked CSD from all animals. Out of the 58 animals in the dataset, 5 did not have readable recordings of LFP in V1 during the presentation of the full-field flash stimuli, and the exact probe locations in V1 could not be recovered for 9 other animals due to fading of fluorescent dye or artifacts in the optical projection tomography (OPT) volume (see ‘Materials and methods’). The remaining 44 (out of the 58) animals in the dataset were retained for the CSD analysis. The trial-averaged CSD plots of all these 44 animals are displayed in Figure 2—figure supplement 1. The first principal component (PC 1) (Figure 2B) constitutes a sum of weighted contributions of the CSD patterns from all 44 animals and explains half (50.4%) of the variance. The salient features typically observed in individual animals are also prominent in the PC 1 CSD pattern (Figure 2B), that is, the canonical pattern. In Figure 2—figure supplement 2, the first 10 principal components cumulatively explaining 90% of the variance are plotted. Quantifying CSD pattern similarity We use the Wasserstein, or Earth Mover’s, distance (WD) to the in CSD patterns (see 'Materials and which can then be used to how well the simulated CSD the CSD typically observed in experiments. The reflects the of one distribution into another by its around et al., 1998; et al., An to the two as two of where the the amount of work that be to the of one around until its distribution the other et al., In the of CSD the reflects the of the distribution of sinks and sources in one CSD pattern into the distribution of sinks and sources in another with greater between CSD patterns. The WDs are computed between the sinks of two CSD patterns and between the sources of two CSD patterns and then to form a total between the CSD patterns (Figure The sum of all sinks and the sum of all sources in each CSD pattern are normalized to and respectively, the only reflects in and not in the The with in and the WDs between the evoked CSD patterns of individual animals and the canonical we find that the animals with CSD patterns by visual the canonical pattern (Figure 2A, animals are animals with while the animal with the more CSD pattern (Figure 2A, animal is an (Figure The of the evoked response is less in the single-trial CSD due to sinks and sources, but there is a in from 40 to 50 ms onward (Figure with the of spiking responses to full-field in V1 (Siegle et al., 2021). An of sinks and sources with a of ms, that is, in the range is apparent in the from to the of which to be or out by more sustained sinks and sources emerging at about 60 ms. of this activity from the visual flash that covers the visual field and that neurons and in the in an (see the in the firing in Figure Figure 3 with 2 supplements see all Download asset Open asset Variability in experimentally recorded spikes. (A) laminar population firing rates of cells, by and across all layers in response to full-field flash. across all animals. (B) Kolmogorov–Smirnov (see 'Materials and between the trial-averaged firing rates of each individual animal and the firing over from all animals line in at of 250 ms flash evoked response to 60 ms after flash and during the sustained 60 to ms). (C) between trial-averaged firing rates of individual mice and all mice ms after flash (D) evoked firing rates for excitatory in visual areas over trials, neurons, and the in the firing of neurons (blue). (E) of population firing rates during evoked response, and the sustained across trials, neurons and time The inter-trial variability is to the inter-animal variability of the trial-averaged responses. the pairwise Wasserstein distances between single CSDs each animal and comparing it to the pairwise between the trial-averaged CSD of each we find that inter-trial variability in CSD is than the inter-animal variability in trial-averaged CSD distance = (Figure The of animals out of have a to the first principal PC 1, of the CSD that is less than half of the between the CSD of individual animals and the PC 1 CSD (Figure the pairwise WDs between animals are also less than half of the maximum pairwise for most animals out of the total pairwise Figure This the that most animals exhibit the canonical CSD pattern captured by the PC 1 CSD (Figure The total inter-trial variability is than the inter-animal both estimated by pairwise WDs (Figure and there are n = 5 animals for which the inter-trial WDs are than the inter-animal WDs (Figure by red with on the distribution of pairwise WDs between animals and pairwise WDs between trials in each see Figure 2—figure supplement 3). Quantifying firing variability For the we between neurons and excitatory and neurons (see 'Materials and and Figure supplement are into one population across all while the are into populations for each layer (Figure The are across layers because we a of at 10 recorded neurons in one layer comparing the population firing in individual animals to the population firing in all animals, and only one animal 10 or more in layer (Figure supplement This was to have a more of the population firing rates in individual animals. We use the similarity as one the see 'Materials and and to the variability in firing rates. We use the experimental variability as a to the model firing rates typically observed in experiments. The similarity the similarity between the of firing rates across neurons in two populations in time with similarity = 1 such, similarity provides a to the of firing rates in time We the as the over 250 ms the flash onset, the as ms after flash onset, and the as ms after flash The similarity during is during the and The on the other is computed between two population firing rates the ms The a measure of the similarity in the of firing rates in this of We establish the experimental variability in and by these between the population firing rates of each individual animal and the population firing rates of all other animals (averaged over trials for both the individual animals and the over all other (Figure and The population firing rates for neurons are more than as as during and the the firing in is the in all followed by L4 and while has the firing rates (Figure between the model and experimental We simulated the response to a full-field flash stimulus with the biophysical network model of mouse primary visual cortex as in Billeh et al., 2020. input to the we used experimentally recorded (Figure see ‘Materials and methods’). A firing at a of 1 provides additional synaptic input to all cells, the from the of the brain The thalamocortical input consists of from units et al., 2018; Billeh et al., The public Neuropixels data recordings from neurons across mice during 75 trials of full-field bright flash resulting in To the input for each of our 10 trials, we 10 unique of from this until all units been a in each trial. Figure 4 with 1 supplement see all Download asset Open asset Local field potential current source density and spikes from simulations with the model. (A) Top: plot of all in the 400 μm radius ‘core’ all in a of a single with the flash Bottom: plot and of spikes from 10 trials for an example (B) Top: simulated LFP averaged over 10 trials of flash Bottom: CSD calculated from the LFP the delta iCSD method. (C) of experimentally recorded used as input to the model. (D) Wasserstein distance between CSD from the model and PC 1 CSD from experiments with the Wasserstein distances from experimental CSD in animal to PC 1 CSD normalized to maximal distance for animals. (E) recorded firing rates and simulated firing rates (blue). (F) Kolmogorov–Smirnov similarity between firing rates in model or individual animals and firing rates in experiments at evoked response, and during the sustained in Figure 3). (G) between firing rates of model or individual animals in experiments and population firing rates in experiments ms). of model firing rates during evoked response, and the sustained across trials, neurons and time Figure the resulting spiking pattern across all layers with its associated LFP. The CSD a sink in the and by a source both starting at ms after flash (Figure However, the early L4 the later sustained L4 and the sustained sink typically observed in the experimental CSD (Figure 2A, are or to the sink and source in and The from the simulated CSD to the experimental PC 1 CSD is greater than the between the CSD of the animal and the PC 1 CSD = normalized to the between CSD of individual animals and PC 1 Thus, using experimental variability as a the CSD from this is an (Figure The population firing rates of the the and