Uncovering circuit mechanisms of current sinks and sources with biophysical simulations of primary visual cortex

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 numerous cortical 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. Our findings show that CSD analysis provides powerful constraints for modeling beyond those from considering spikes.


36
The local field potential (LFP) is the low-frequency component (below a few hundred Hertz) 37 of the extracellular potential recorded in brain tissue that originates from the transmembrane 38 currents in the vicinity of the recording electrode (Lindén et al., 2011; sum to zero over the total membrane area of a neuron). However, the interpretation of 71 current sinks and sources is inherently ambiguous, as several processes can be the origin of 72 a current sink or source (Buzsáki, 2006;Pettersen et al., 2006;Einevoll et al., 2007). For 73 example, a current source may reflect an inhibitory synaptic current or an outflowing return 74 current resulting from excitatory synaptic input elsewhere on the neuron, and there is no 75 simple way of knowing which it is from an extracellular recording alone (Buzsáki, 2006).

77
Another approach to uncovering the biophysical origins of current sinks and sources, and by 78 extension the LFP, is to simulate them computationally (Pettersen,Hagen,and Einevoll,79 2008; Einevoll et al., 2013). Following the classic work by Rall in the 1960s (Rall, 1962), a 80 forward-modeling scheme in which extracellular potentials are calculated from neuron 81 models with detailed morphologies using volume conduction theory under the line source 82 approximation has been established (Holt and Koch, 1999). With this framework, we have

93
While interesting insights about CSD and LFP were obtained from these computational 94 approaches, establishing a direct relationship between the biological details of the circuit 95 structure and the electrical signal like LFP remains a major unresolved challenge. One 96 reason is that the amount and quality of data available for modeling the circuit architecture in 97 detail has been limited. This situation improved substantially in recent years, and a broad 98 range of data on the composition, connectivity, and physiology of cortical circuits have been 99 integrated systematically (Billeh et al., 2020) in a biophysically detailed model of mouse

144
Spikes and LFP were recorded across multiple brain areas, with a focus on six cortical (V1, 145 LM, AL, RL, AM, PM) and two thalamic (LGN, LP) visual areas, using Neuropixels probes in 146 58 mice (Siegle et al., 2021).  to scale with the model.

Uncovering a canonical visually evoked CSD response 195
We first established a "typical" experimentally recorded CSD pattern, to be reproduced with 196 the model. Though there is substantial inter-trial and inter-animal variability in the evoked 197 CSD response, we find that most trials and animals have several salient features in 198 common. In Fig. 2A, the trial-averaged evoked CSDs from five individual mice are displayed.

199
In the first four animals (# 1-4), we observe an early transient sink arising in layer 4 (L4) 200 around 40 ms after flash onset, followed by a sustained source starting at about 60 ms which 201 covers L4 and parts of layers 2/3 (L2/3) and layer 5 (L5). We also observe a sustained sink 202 covering layers 5 and 6 (L6) emerging at around 50 ms, as well as a sustained sink covering 203 layers 1 and 2/3 from about 60 ms. An animal that does not fully exhibit what we term 204 the "canonical" pattern is shown in the rightmost plot (# 5 in Fig. 2A); it has an early L4 sink 205 arising at 40 ms, but this sink is not followed by the sustained sinks and sources from 50-60 206 ms and onwards observed in the other animals. The timing and location of sinks and sources 207 are, overall, similar to those described earlier by Niell & Stryker (2008), and Senzai,

253
To identify the robust features across animals in this dataset, we performed Principal

Quantifying CSD pattern similarity 265
We use the Wasserstein or Earth Mover's distance (WD), to quantify the differences in CSD 266 patterns (see Methods). The WD reflects the cost of transforming one distribution into 267 another by moving its "distribution mass" around (Rubner et al., 1998;Arjovsky et al., 2017).

268
An often-used analogy refers to the two distributions as two piles of dirt, where the WD tells 269 us the minimal amount of work that must be done to move the mass of one pile around until 270 its distribution matches the other pile (Rubner et al., 1998). In the context of CSD patterns,

286
The onset of the evoked response is less conspicuous in the single-trial CSD, due to 287 pronounced, ongoing sinks and sources, but there is still a visible increase in magnitude 288 from 40-50 ms and onwards (Fig. 2E), compatible with the latency of spiking responses to 289 full-field flashes in V1 (Siegle et al., 2021). An oscillation of sinks and sources with a 290 periodicity of ~20 ms, i.e., in the gamma range is apparent in the region stretching from L2/3 291 to the top of L5, which appears to be either partially interrupted or drowned out by more

297
The inter-trial variability is roughly comparable to the inter-animal variability of the trial-298 averaged responses. By computing the pairwise Wasserstein distances between single trial 299 CSDs within each animal, and comparing it to the pairwise WD between the trial-averaged 300 CSD of each animal, we find that inter-trial variability in CSD is significantly lower than the 301 inter-animal variability in trial-averaged CSD (Kolmogorov-Smirnov distance = 0.33; p < 302 0.001) (Fig. 2F).

304
The majority of animals (39 out of 44) have a WD to the 1 st principal component, PC 1, of the 305 CSD that is less than half of the greatest WD between the CSD of individual animals and the 306 PC 1 CSD ( Fig 2D); the pairwise WDs between animals are also less than half of the 307 maximum pairwise WD for most animals (921 out of the total 946 pairwise WDs; Fig 2E).

308
This supports the view that most animals exhibit the canonical CSD pattern captured by the 309 PC 1 CSD (Fig. 2B). The total inter-trial variability is smaller than the inter-animal variability,

371
The population firing rates for FS neurons are more than twice as high than RS cells during 372 baseline, peak and sustained. Among the RS populations, the firing rate in L5 is the highest 373 at the peak and baseline, followed by L4 and L6, while L2/3 has the lowest firing rates (

444
Thus, using experimental variability as a reference, the CSD from this simulation is an outlier 445 (Fig. 4C).

447
The population firing rates of the model, the KS similarities and correlation between the

491
The resulting spatial pattern (but not the magnitude) of the CSD, however, was largely 492 unchanged (Fig. 5D) compared to the original CSD (Fig. 4B). The overall magnitude was 493 reduced, and there were some traces of a sink arising at 40 ms after flash onset, and a L2/3 494 (and L1) sink after 60 ms, but they were substantially weaker relative to the L5/L6 dipole

Two-way dissociation between spikes and CSD 548
Simulations demonstrate that the LFP, and the associated CSD, can be significantly altered  leftmost plot), while the magnitude of its firing rates were reduced (Fig. 5G). In this case, the 571 pattern displayed a sink in the middle with a source above and below it.

573
These results indicate a two-way dissociation that can occur between CSD and firing rates of

Identifying the biophysical origins of the canonical CSD 670
With the canonical CSD (Fig. 2B) reproduced, we can use the model to probe the 671 biophysical origins of its sinks and sources. We began by removing all recurrent connections 672 and only feeding the LGN input to the model to find the contribution from the thalamocortical 673 synapses onto excitatory and inhibitory neurons (Fig. 7A, top). The main thalamic 674 contribution to the CSD is from synapses onto excitatory neurons, in line with the 675 expectation that neurons with a spatial separation between synaptic input currents and the 20 return currents dominate the cortical LFP generation . (Neurons without 677 apical dendrites will have largely overlapping synaptic input currents and return currents, 678 resulting in a cancellation of current sinks and sources.)

680
We further observed that the early L4 and the sustained L5/L6 sinks are present in the CSD 681 contributions of excitatory neurons, though the magnitude of the L5/L6 sink is substantially 682 reduced compared to its magnitude when the model is configured with recurrent synapses 683 intact (Fig. 4B, 5D, and 6D). The sustained L2/3 sink and L4 source, on the other hand, were 684 not visible. This suggests that the early L4 sink and the L5/L6 sink are at least partly 685 generated by thalamocortical synapses. However, the substantially diminished magnitude of 686 the L5/L6 sink indicates that recurrent synapses also contribute significantly to the 687 generation of this sink.

689
We then removed the LGN input and added the feedback (while keeping the recurrent 690 connections cut), which resulted in a prominent upper layer dipole, with the sink residing in 691 L1 and L2/3, and the source residing in L4 (Fig. 7A, bottom). Together with their absence 692 when input came from LGN only (Fig. 7A, top), this suggests that the sustained L2/3 sink 693 and the L4 source in the canonical pattern originate at least in part from the feedback 694 synapses onto the apical dendrites of L2/3 and L5 pyramidal cells and the activity this input 695 generates.

697
To assess the extent to which active channels at the somata contributed to the CSD pattern, in the model. The most prominent discrepancy between the CSD with and without active 703 channels is the magnitude of the L4 source and the L2/3 sink (Fig. 7B). In this all-passive 704 setting, the L4 source is significantly attenuated, and the L2/3 sink is either absent or 705 dominated by a source in the same region.

707
We explored whether the contributions from currents in recurrent connections come primarily 708 from excitatory or inhibitory synapses by removing all connections from inhibitory (Pvalb, 709 Sst, Htr3a) neurons to all other neurons, so that all postsynaptic currents stem from 710 excitatory thalamocortical synapses, excitatory synapses from higher visual areas, or 711 recurrent excitatory synapses in V1 (Fig. 7D and Fig. S5). Note that inhibitory synaptic 712 currents give rise to sources, while excitatory synaptic currents give rise to a sink. Of course, 713 21 without inhibition, the network is unbalanced, which limits the conclusions that can be drawn 714 from this simulation. However, the fact that the major sinks and sources are still present is 715 an indication that the currents from excitatory input account for the majority of the sinks and 716 sources observed in the experimental CSD.

759
The contributions from each population to the total CSD in the final model ( Fig. 6D) with both

760
LGN and feedback input and intact recurrent connections are displayed in Fig. 7D. From

798
(2019). The L4 source and L2/3 sink were also observed in that study, but emerge 799 somewhat later than in our data -just after 100 ms as opposed to ~60 ms in our canonical 800 pattern. This discrepancy in onset might simply be due to differences in stimuli. In Senzai,

801
Fernandez-Ruiz, and Buzsáki (2019) the animals were exposed to 100 ms light pulses, while 802 the animals in our data were presented with 250 ms whole-field flashes of white screens.

803
Nonetheless, the canonical CSD pattern exhibits good overall agreement with the pattern 804 seen in Senzai, Fernandez-Ruiz, and Buzsáki (2019).

806
We introduced the Wasserstein distance as a method to evaluate the difference between 807 two CSD patterns and used it to quantify the variability in trial-averaged CSD between 808 animals (Fig. 2D), the trial-to-trial variability in CSD within animals (   feedback can play a significant role in shaping the sustained sinks and sources (Fig. 6B-D).

865
The resulting model CSD reproduced the major sinks and sources identified in the canonical 866 CSD pattern and was no longer among the outliers when compared to the experimental 867 variability (Fig. 6E). Interestingly, absence of the feedback was not apparent from analysis of

886
In line with observations made by Reimann et al. (2013), we found that the somatic voltage-887 dependent membrane currents significantly shape the CSD signature ( Fig. 5H and 7B). Even 888 so, our findings still emphasize the importance of synaptic inputs in sculpting the CSD, as 889 the addition of synaptic input (Fig. 6D) and changes to synaptic placement (Fig. 5F) 890 substantially altered the CSD pattern.

892
This investigation into the biophysical origins of sinks and sources is limited by the fact that 893 the contributions from recurrent connections are difficult to estimate precisely due to the non-894 26 linear effects of these connections in terms of how they contribute to spiking. That is, their 895 contribution cannot simply be found by subtracting the CSD from thalamocortical and 896 feedback synapses with all recurrent connections removed (Fig. 7A) from the total CSD with 897 the same input and recurrent connections intact (Fig. 6D). Still, this analysis provides an 898 initial estimate into the biophysical origins of the sinks and sources observed experimentally 899 and demonstrates the insights that can be obtained from modeling of extracellular signals.

901
There is ample evidence that firing rates and LFP are modulated by the behavioral state of 902 the animal, including measures like the pupil size (considered to be a proxy for arousal level) stimulus (out-class generalization). This is a long-term goal, and can be called "the holy 930 grail" of visual system modeling.

931
27 932 In this study, we developed a systematic framework to quantify experimental variability in 933 both LFP/CSD and spikes and to evaluate model performance. We identified a canonical 934 CSD pattern observed during presentations of full-field flash stimuli and obtained a bio-935 realistic model that reproduced both the canonical CSD pattern and spikes simultaneously.

936
We utilized this validated model to identify the biophysical origins of the canonical sinks and