## Abstract

With rapidly advancing multi-electrode recording technology, the local field potential (LFP) has again become a popular measure of neuronal activity in both research and clinical applications. Proper understanding of the LFP requires detailed mathematical modeling incorporating the anatomical and electrophysiological features of neurons near the recording electrode, as well as synaptic inputs from the entire network. Here we propose a hybrid modeling scheme combining efficient point-neuron network models with biophysical principles underlying LFP generation by real neurons. The LFP predictions rely on populations of network-equivalent multicompartment neuron models with layer-specific synaptic connectivity, can be used with an arbitrary number of point-neuron network populations, and allows for a full separation of simulated network dynamics and LFPs. We apply the scheme to a full-scale cortical network model for a ∼1 mm^{2} patch of primary visual cortex, predict laminar LFPs for different network states, assess the relative LFP contribution from different laminar populations, and investigate effects of input correlations and neuron density on the LFP. The generic nature of the hybrid scheme and its public implementation in

## Introduction

The local field potential (LFP), the low-frequency component ($\u2272500Hz$) of the extracellular potential recorded in the brain, is commonly used as a measure of neuronal activity (Buzsáki et al. 2012; Einevoll et al. 2013). The LFP originates from transmembrane currents (Nicholson and Freeman 1975), and at the single-cell level the biophysical origin of such extracellular potentials is well understood (see, e.g., Rall and Shepherd 1968; Holt and Koch 1999; Buzsáki et al. 2012; Einevoll et al. 2013). However, the interpretation of the LFP remains difficult due to the large number of neurons contributing to the recorded signal. In neocortex, for example, the measured LFP is typically generated by thousands or even millions of neurons near the recording electrode (Kajikawa and Schroeder 2011; Lindén et al. 2011; Łe¸ski et al. 2013). Moreover, the LFP reflects also synaptic input generated by remote populations, for example, inputs from other cortical or subcortical areas in addition to local network interactions (Herreras et al. 2015). A thorough theoretical description of the LFP, therefore, needs to account not only for the anatomical and electrophysiological features of neurons in the vicinity of the recording electrode, but also for the entire large-scale neuronal circuitry generating synaptic input to these cells.

Modeling large-scale neural-network dynamics with individual spiking neurons is challenging due to the memory required to represent the large number of synapses. With current technology and using the largest supercomputers available today, simulations of neural networks comprising up to 10^{9} neurons and 10^{13} synapses (roughly corresponding to the size of a cat brain) are feasible for simplified model neurons (Diesmann 2013; Kunkel et al. 2014). Typically, these simplified models neglect the spatial aspects of neuronal morphologies and describe neurons as points in space (point-neuron models). Despite their simplicity, point-neuron-network models explain a variety of salient features of neural activity observed in vivo, such as spike-train irregularity (Softky and Koch 1993; van Vreeswijk and Sompolinsky 1996; Amit and Brunel 1997; Shadlen and Newsome 1998), membrane-potential fluctuations (Destexhe and Paré 1999), asynchronous firing (Ecker et al. 2010; Renart et al. 2010; Ostojic 2014), correlations in neural activity (Gentet et al. 2010; Okun and Lampl 2008; Helias et al. 2013), self-sustained activity (Ohbayashi et al. 2003; Kriener et al. 2014), and realistic firing rates across laminar cortical populations (Potjans and Diesmann 2014). Point-neuron networks are amenable to mathematical analysis (see, e.g., Brunel 2000; Deco et al. 2008; Tetzlaff et al. 2012; Helias et al. 2013; de Kamps 2013; Schuecker et al. 2015; Bos et al. 2016) and can be efficiently evaluated numerically (Brette et al. 2007; Plesser et al. 2007; Helias et al. 2012; Kunkel et al. 2014). The mechanisms governing networks of biophysically detailed multicompartment model neurons, in contrast, are less accessible to analysis and these models are more prone to overfitting. Existing multicompartment neuron network models accounting for realistic cell morphologies are restricted to sizes of ∼10^{4}–10^{5} neurons (Hines et al. 2008; Reimann et al. 2013; Migliore et al. 2014; Markram et al. 2015). Large-scale models are, however, necessary to include contributions to the LFP from distant populations in situations where the spatial reach of the LFP is known to be large (Lindén et al. 2011; Łe¸ski et al. 2013).

Although point-neuron networks capture many features of in vivo spiking activity, they fail to predict extracellular potentials that result from transmembrane currents distributed across the cell surface. According to Kirchhoff's law of current conservation, the sum of all transmembrane currents, including all ionic and capacitive currents, must be zero for each neuron. In a point-neuron model, all transmembrane currents are collapsed in a single point in space. The net transmembrane current, and hence the extracellular potential, therefore vanishes. Only the spatial separation between current sinks and sources leads to a nonzero extracellular potential (Pettersen et al. 2012; Einevoll et al. 2013). A priori, the prediction of extracellular potentials, therefore, requires spatially extended neuron models accounting for the spatial distribution of transmembrane currents, commonly handled using multicompartment neuron models (De Schutter and Van Geit 2009). Note that the principle of current conservation implies a current sum rule for multicompartment neuron models as well: the sum of all single-cell transmembrane currents remains zero, also across neuron populations and the whole column. In several previous studies (Bazhenov et al. 2001; Hill and Tononi 2005; Ursino and La Cara 2006; Mazzoni et al. 2008; 2010; 2011), the activity of point-neuron networks (e.g., population firing rates, synaptic currents and membrane potentials) has nevertheless been used as a proxy for the LFP when comparing with experiments. In a recent study comparing different candidate proxies, it was found that a suitably chosen sum of synaptic currents could provide a good LFP proxy, but only for the case when the LFP is generated from transmembrane currents of a single population of pyramidal neurons (Mazzoni et al. 2015). In cortex, however, several populations in general contribute to the LFP, and there are spatial cancellation effects when positive LFP contributions from one population overlap in space with negative LFP contributions from other populations. This effect cannot be accounted for by a simple LFP proxy.

In this article, we present a hybrid modeling scheme that combines the simplicity and efficiency of point-neuron network models and the biophysical principles underlying LFP generation captured by multicompartment neuron models with anatomically reconstructed morphologies. The scheme allows for arbitrary numbers of LFP-contributing populations, and directly incorporates spatial cancellation effects. Furthermore, the spatially extended LFP-generating neurons assure that effects from intrinsic dendritic filtering of synaptic inputs are included in the predicted LFP (Lindén et al. 2010). The scheme assumes that the spiking activity of the neural network (Fig. 1*B*) generating the synaptic input reflected in the LFP is well described by a point-neuron network model (Fig. 1*A*). The network spiking activity serves as synaptic input to a population of mutually unconnected multicompartment model neurons with realistic morphologies positioned in 3-dimensional (3D) space (Fig. 1*C*) and is thereby translated into a distribution of transmembrane currents and, hence, an LFP (Fig. 1*D*). Thus each multicompartment model neuron has its equivalent in the point-neuron network and receives input spikes from the same presynaptic neurons as this point-neuron equivalent.

In the proposed hybrid modeling scheme, the LFP stems from the presynaptic spiking activity, but does not affect the spike-generation dynamics. Thus, the modeling of the spike trains and the LFP generation are separated so that the effects of the spatial and electrophysiological properties of the postsynaptic (multicompartment) neurons on the LFP can be investigated independently of the spike-generation dynamics. Due to the linearity of Maxwell's equations and volume conduction theory linking transmembrane currents to extracellular potentials (Pettersen et al. 2012; Einevoll et al. 2013), the compound LFP results from the linear superposition of all single-cell LFPs generated by the collection of neurons in the multicompartment model neuron population (Einevoll et al. 2013). Note that this linear superposition principle applies even for nonlinear cell dynamics (e.g., nonlinear synaptic integration, action-potential generation and active conductances) as in Reimann et al. 2013. As ephaptic interactions (Anastassiou et al. 2011) are neglected, the LFP contribution from each multicompartment model neuron can be treated independently from the others. The computational hybrid LFP scheme proposed here exploits the methodological and conceptual advantages due to the independence of the contributions to the LFP from each multicompartment model neuron: the evaluation of the LFPs becomes “embarrassingly parallel” (see Foster 1995) and simulations of the multicompartment model neuron dynamics can be easily distributed in parallel across many compute units (i.e., CPUs). Although tailored toward use on high-performance computing facilities, the hybrid simulation can in principle be run on a single laptop.

The hybrid scheme predicts spatially and temporally resolved neural activity at various scales: spikes, synaptic currents, membrane potentials, current-source densities (CSD, see e.g., Nicholson and Freeman 1975; Pettersen et al. 2006; 2008), and LFPs. It therefore allows for investigation of relationships between different measures of neural activity. Thus, although point-neuron networks until now only have connected to in vivo experiments via measurement of spikes, single-neuron membrane potentials and currents, the present hybrid scheme allows for comparison of model predictions also with measured LFPs (and associated CSDs).

As an illustration, we apply the hybrid scheme to a multilayered point-neuron network model of an early sensory cortical microcircuit (Potjans and Diesmann 2014). We thereby demonstrate how to obtain LFP predictions from point-neuron network models using additional spatial connectivity information from anatomical data (Binzegger et al. 2004; Izhikevich and Edelman 2008). The example illustrates how the hybrid scheme can be used to examine the relation between single-neuron and population signals, that is, spikes and LFPs, the effect of network dynamics on the LFP, and the interpretation of the LFP in terms of underlying laminar neuron populations. We further use the example to demonstrate that synaptic-input correlations result in a nontrivial dependence of the LFP on the neuron density. Correct LFP predictions can therefore only be obtained by accounting for realistic neuron densities.

The network model of Potjans and Diesmann (2014) is chosen here since it has a minimum level of detail in the sense that individual neurons have simplified leaky integrate-and-fire (LIF) dynamics, but still represents a cortical column with full density of neurons and connections. The connectivity in such a full-scale circuit alone suffices to explain realistic firing rates across populations as well as propagation of activity through layers (Potjans and Diesmann 2014). Applicability of the scheme is, however, not restricted to this model as it in principle can be used for all network models generating spikes.

In Methods and Materials, we detail the components of the hybrid scheme and their application to the cortical microcircuit model: the point-neuron-network model, the populations of multicompartment neurons, the synaptic connectivity of the point-neuron network and the multicompartment model neuron populations, and the biophysical forward-modeling scheme of extracellular potentials. We further describe the analysis of the data generated by the simulations, as well as the

## Methods and Materials: Hybrid LFP Modeling Scheme

### Point-Neuron Network Model

The point-neuron network model is a key component of the hybrid scheme. The hybrid scheme enables LFP predictions from network models with an arbitrary number of populations and thus permits application to a large class of networks with arbitrarily complex single-neuron and synapse dynamics. The example network of “spike-generators” used here is, except for some minor adjustments (see below), the multilayered model of a cortical microcircuit published by Potjans and Diesmann (2014). The model is implemented and included in

The network model describes 1 mm^{2} of primary sensory cortex and consists of 4 layers with one excitatory (E) and one inhibitory (I) neuron population each, as illustrated in Fig. 1*A*. The network receives modulated thalamic input in addition to stationary external input. While the neuron (LIF) and synapse (static, exponential-current-based) model are intentionally left simple, the focus of this network implementation is on the complex connectivity which integrates multiple sources of anatomical and electrophysiological data (Potjans and Diesmann 2014). Apart from the layer identity, the model does not explicitly account for cell positions. For the full network description, see Tables 1, 2, and 5. The microcircuit model reproduced experimentally observed distributions of firing rates across populations and propagation of activity across layers (Potjans and Diesmann 2014). It thus forms a suitable starting point for LFP predictions in a cortical column.

A | Model summary |
---|---|

Structure | Multilayered excitatory-inhibitory (E–I) network |

Populations | 8 cortical in 4 layers, 1 thalamic (TC) |

Connectivity | Random, independent, population-specific, fixed number of connections |

External input | Cortico-cortical: constant current with population-specific strength |

Neuron model | Cortex: LIF; TC: point process |

Synapse model | Exponential postsynaptic currents, static weights, population-specific weight distributions |

Measurements | Spike activity, input currents, membrane potential of each neuron |

B | Network model |

Connectivity | Connection probability C_{YX} ($X,Y\u2208{L2/3,L4,L5,L6}\xd7{E,I}\u222aTC$, $CYX=0$ for $Y=TC$) |

•Fixed number of synapses K_{YX} between populations X and Y | |

•Binomial in-/out-degrees | |

Input | Cortico-cortical direct current $IYext$ |

C | Neuron model |

Cortex | |

Type | Leaky integrate-and-fire neuron (LIF) |

Description | Dynamics of membrane potential V_{i}(t) (neuron $i\u2208[1,N]$): |

•Spike emission at times t^{i}_{l} with $Vi(tli)\u2a7e\theta $ | |

•Subthreshold dynamics: $\tau mVi\u02d9=\u2212Vi+RmIi(t)$ if $\u2200l:t\u2209(tli,tli+\tau ref]$ | |

•Reset + refractoriness: $Vi(t)=Vreset$ if $\u2200l:t\u2208(tli,tli+\tau ref]$ | |

Exact integration with temporal resolution dt (Rotter and Diesmann 1999) | |

Uniform distribution of membrane potentials at t = 0 | |

Thalamus | |

Type | •DC current for constant background input |

•Nonstationary Poisson process for modulation | |

Description | DC current included in external DC input |

Types of thalamic input modulation: | |

•Spontaneous activity: no modulation in activation of thalamic neurons | |

•Thalamic pulses: fixed-interval coherent activation of all $NTC$ thalamic neurons | |

•AC modulation: Poisson spike trains with sinusoidally modulated rate profile (discretized with time resolution dt): (18) $\nu th(t)=\nu \xafTC+\Delta \nu TCsin(2\pi tfTC)$ |

A | Model summary |
---|---|

Structure | Multilayered excitatory-inhibitory (E–I) network |

Populations | 8 cortical in 4 layers, 1 thalamic (TC) |

Connectivity | Random, independent, population-specific, fixed number of connections |

External input | Cortico-cortical: constant current with population-specific strength |

Neuron model | Cortex: LIF; TC: point process |

Synapse model | Exponential postsynaptic currents, static weights, population-specific weight distributions |

Measurements | Spike activity, input currents, membrane potential of each neuron |

B | Network model |

Connectivity | Connection probability C_{YX} ($X,Y\u2208{L2/3,L4,L5,L6}\xd7{E,I}\u222aTC$, $CYX=0$ for $Y=TC$) |

•Fixed number of synapses K_{YX} between populations X and Y | |

•Binomial in-/out-degrees | |

Input | Cortico-cortical direct current $IYext$ |

C | Neuron model |

Cortex | |

Type | Leaky integrate-and-fire neuron (LIF) |

Description | Dynamics of membrane potential V_{i}(t) (neuron $i\u2208[1,N]$): |

•Spike emission at times t^{i}_{l} with $Vi(tli)\u2a7e\theta $ | |

•Subthreshold dynamics: $\tau mVi\u02d9=\u2212Vi+RmIi(t)$ if $\u2200l:t\u2209(tli,tli+\tau ref]$ | |

•Reset + refractoriness: $Vi(t)=Vreset$ if $\u2200l:t\u2208(tli,tli+\tau ref]$ | |

Exact integration with temporal resolution dt (Rotter and Diesmann 1999) | |

Uniform distribution of membrane potentials at t = 0 | |

Thalamus | |

Type | •DC current for constant background input |

•Nonstationary Poisson process for modulation | |

Description | DC current included in external DC input |

Types of thalamic input modulation: | |

•Spontaneous activity: no modulation in activation of thalamic neurons | |

•Thalamic pulses: fixed-interval coherent activation of all $NTC$ thalamic neurons | |

•AC modulation: Poisson spike trains with sinusoidally modulated rate profile (discretized with time resolution dt): (18) $\nu th(t)=\nu \xafTC+\Delta \nu TCsin(2\pi tfTC)$ |

A | Model summary |
---|---|

Topology | Cortical column under $1mm2$ of cortical surface |

Populations | 8 excitatory and 8 inhibitory cell types |

Input | Spiking activity of thalamic and cortical populations as modeled by point-neuron network |

Neuron model | Multicompartment, passive cable formalism |

Synapse model | Exponential postsynaptic current, static weights |

Measurements | Current source density (CSD), local field potential (LFP) |

B | Topology |

Type | Cylindrical volume with layer-specific distribution of cell types and synapses |

Description | Cylinder radius r |

Laminar, defining upper/lower boundaries of layers 1, 2/3, 4, 5, 6 | |

C | Populations |

Type | Each cell type y assigned to population Y, $y\u2208Y$ |

Description | Populations $Y\u2208{L2/3,L4,L5,L6}\xd7{E,I}$ (population size N_{Y}, cell types $y\u2208Y$) |

(e.g., $L4E={p4,ss4(L23),ss4(L4)}$, cf., Fig. 2). | |

Cell types y: | |

•Size $Ny=Fy\u2211YNY$, F_{y} is the occurrence of cell type y in the full model | |

•Morphology M_{y} | |

•Extrapolated according to spatial connectivity data (Table 7) | |

Somatic placement, population Y: | |

•Random soma placement in cylindrical volumes with radius r, thickness h | |

•Volumes centered between boundaries of layers 2/3–6 | |

Morphologies | |

Type | 3D histological reconstructions from slice preparations (see Jacobs et al. 2009; De Schutter and Van Geit 2009) of cat visual and somatosensory cortices |

Description | One morphology M_{y} per cell type: |

•Excitatory and inhibitory cells in layers 2/3–6 | |

•For all cells $j\u2208y:$M_{j} = M_{y} | |

•For some cell types $y,y\u2032:$$My=My\u2032$ (limited availability) | |

Orientations: | |

•Pyramidal cells: apical dendrites oriented along depth axis with random depth-axis rotation | |

•Interneurons, stellate cells: random rotation around all axes | |

Corrections: | |

•Apical dendrites of pyramidal cells elongated to accommodate spatial connectivity | |

•Axons removed if present | |

Reconstructed morphologies (cf., Fig. 2): | |

•Cat visual cortex (Kisvárday and Eysel 1992; Mainen and Sejnowski 1996; Contreras et al. 1997; Stepanyants et al. 2008) | |

•Cat somatosensory cortex from NeuroMorpho.org (Contreras et al. 1997; Ascoli et al. 2007). |

A | Model summary |
---|---|

Topology | Cortical column under $1mm2$ of cortical surface |

Populations | 8 excitatory and 8 inhibitory cell types |

Input | Spiking activity of thalamic and cortical populations as modeled by point-neuron network |

Neuron model | Multicompartment, passive cable formalism |

Synapse model | Exponential postsynaptic current, static weights |

Measurements | Current source density (CSD), local field potential (LFP) |

B | Topology |

Type | Cylindrical volume with layer-specific distribution of cell types and synapses |

Description | Cylinder radius r |

Laminar, defining upper/lower boundaries of layers 1, 2/3, 4, 5, 6 | |

C | Populations |

Type | Each cell type y assigned to population Y, $y\u2208Y$ |

Description | Populations $Y\u2208{L2/3,L4,L5,L6}\xd7{E,I}$ (population size N_{Y}, cell types $y\u2208Y$) |

(e.g., $L4E={p4,ss4(L23),ss4(L4)}$, cf., Fig. 2). | |

Cell types y: | |

•Size $Ny=Fy\u2211YNY$, F_{y} is the occurrence of cell type y in the full model | |

•Morphology M_{y} | |

•Extrapolated according to spatial connectivity data (Table 7) | |

Somatic placement, population Y: | |

•Random soma placement in cylindrical volumes with radius r, thickness h | |

•Volumes centered between boundaries of layers 2/3–6 | |

Morphologies | |

Type | 3D histological reconstructions from slice preparations (see Jacobs et al. 2009; De Schutter and Van Geit 2009) of cat visual and somatosensory cortices |

Description | One morphology M_{y} per cell type: |

•Excitatory and inhibitory cells in layers 2/3–6 | |

•For all cells $j\u2208y:$M_{j} = M_{y} | |

•For some cell types $y,y\u2032:$$My=My\u2032$ (limited availability) | |

Orientations: | |

•Pyramidal cells: apical dendrites oriented along depth axis with random depth-axis rotation | |

•Interneurons, stellate cells: random rotation around all axes | |

Corrections: | |

•Apical dendrites of pyramidal cells elongated to accommodate spatial connectivity | |

•Axons removed if present | |

Reconstructed morphologies (cf., Fig. 2): | |

•Cat visual cortex (Kisvárday and Eysel 1992; Mainen and Sejnowski 1996; Contreras et al. 1997; Stepanyants et al. 2008) | |

•Cat somatosensory cortex from NeuroMorpho.org (Contreras et al. 1997; Ascoli et al. 2007). |

D | Neuron models |
---|---|

Type | Passive, multicompartment, reconstructed morphologies |

Description | Compartment n membrane potential $Vmjn$ of cell j having length l_{jn}, diameter d_{jn} and surface area A_{jn}: (22) $CmjndVmjndt=\u2211k=1mIajkn\u2212GLjn(Vmjn\u2212EL)\u2212\u2211iIjin,$ (23) $Cmjn=cmAjn,$ (24) $Iajkn=Gajkn(Vmjk\u2212Vmjn),$ (25) $Gajkn=\pi (djk2+djn2)/(4ra(ljk+ljn)),$ (26) $GLjn=Ajn/rm,$ (27) $Imjn=CmjndVmjndt+GLjn(Vmjn\u2212EL)+\u2211iIjin.$ n and neighboring compartment k (out of m compartments), $Gajkn$ axial conductance between n and k, I_{jin} synaptic currents, and $Imjn$ transmembrane current of compartment n. For specific parameter values, see Table 6. Membrane potentials and transmembrane currents are computed using NEURON through LFPy (Carnevale and Hines 2006; Lindén et al. 2014), assuming the extracellular potential to be zero everywhere on the outside of the neuron, that is, an infinite extracellular conductivity. |

E | Synapse model |

Type | Exponential postsynaptic current, static weights |

Description | Neuron j input current of synapse formed with presynaptic neuron i: (28) $Iji(t)=Iji,max\u2211lexp(\u2212(t\u2212tli\u2212di)/\tau s)H(t\u2212tli\u2212di),$ (29) $Iji,max=\u2212Cm\mu YXofpoint\u2212neuronnetwork,$ (30) $H(t)=1fort\u2a7e0,and0elsewhere.$ |

•Delays d_{i} from Gaussian distribution with mean d_{X} ($i\u2208X$), relative standard deviation $\sigma d,rel$ | |

•Synapse activation times: network spike trains plus delay | |

•No cortico-cortical connections: $Iext=0$ (cf., Table 5) | |

F | Input |

Type | Spike times t^{i}_{l} of spiking neuron network (including thalamic input spikes), no cortico-cortical input |

Description | Synapse placement, postsynaptic cell $j\u2208y,y\u2208Y$ (see Methods): |

•Number of synapses from presynaptic population X in layer L: k_{yXL} (Eq. 9) | |

•Compartment specificity of connections: $Ajn/\u2211n\u2208LAjn$, compartment $n\u2208L$ | |

•Synapse locations within layers are chosen randomly among dendritic compartments only | |

G | Measurements |

Type | Local field potential (LFP) and current source density (CSD) |

Description | Laminar multi-electrode, see parameter values in Table 6: |

•Axis perpendicular to pial surface | |

•$ncontacts:$ number of contacts | |

•$hcontacts:$ intercontact distance | |

•$rcontact:$ contact surface radius |

D | Neuron models |
---|---|

Type | Passive, multicompartment, reconstructed morphologies |

Description | Compartment n membrane potential $Vmjn$ of cell j having length l_{jn}, diameter d_{jn} and surface area A_{jn}: (22) $CmjndVmjndt=\u2211k=1mIajkn\u2212GLjn(Vmjn\u2212EL)\u2212\u2211iIjin,$ (23) $Cmjn=cmAjn,$ (24) $Iajkn=Gajkn(Vmjk\u2212Vmjn),$ (25) $Gajkn=\pi (djk2+djn2)/(4ra(ljk+ljn)),$ (26) $GLjn=Ajn/rm,$ (27) $Imjn=CmjndVmjndt+GLjn(Vmjn\u2212EL)+\u2211iIjin.$ n and neighboring compartment k (out of m compartments), $Gajkn$ axial conductance between n and k, I_{jin} synaptic currents, and $Imjn$ transmembrane current of compartment n. For specific parameter values, see Table 6. Membrane potentials and transmembrane currents are computed using NEURON through LFPy (Carnevale and Hines 2006; Lindén et al. 2014), assuming the extracellular potential to be zero everywhere on the outside of the neuron, that is, an infinite extracellular conductivity. |

E | Synapse model |

Type | Exponential postsynaptic current, static weights |

Description | Neuron j input current of synapse formed with presynaptic neuron i: (28) $Iji(t)=Iji,max\u2211lexp(\u2212(t\u2212tli\u2212di)/\tau s)H(t\u2212tli\u2212di),$ (29) $Iji,max=\u2212Cm\mu YXofpoint\u2212neuronnetwork,$ (30) $H(t)=1fort\u2a7e0,and0elsewhere.$ |

•Delays d_{i} from Gaussian distribution with mean d_{X} ($i\u2208X$), relative standard deviation $\sigma d,rel$ | |

•Synapse activation times: network spike trains plus delay | |

•No cortico-cortical connections: $Iext=0$ (cf., Table 5) | |

F | Input |

Type | Spike times t^{i}_{l} of spiking neuron network (including thalamic input spikes), no cortico-cortical input |

Description | Synapse placement, postsynaptic cell $j\u2208y,y\u2208Y$ (see Methods): |

•Number of synapses from presynaptic population X in layer L: k_{yXL} (Eq. 9) | |

•Compartment specificity of connections: $Ajn/\u2211n\u2208LAjn$, compartment $n\u2208L$ | |

•Synapse locations within layers are chosen randomly among dendritic compartments only | |

G | Measurements |

Type | Local field potential (LFP) and current source density (CSD) |

Description | Laminar multi-electrode, see parameter values in Table 6: |

•Axis perpendicular to pial surface | |

•$ncontacts:$ number of contacts | |

•$hcontacts:$ intercontact distance | |

•$rcontact:$ contact surface radius |

A | Global simulation parameters | |||||||||
---|---|---|---|---|---|---|---|---|---|---|

Symbol | Value | Description | ||||||||

T | 5200 ms | Simulation duration | ||||||||

dt | 0.1 ms | Temporal resolution | ||||||||

B | Point-neuron network | |||||||||

Populations and external input | ||||||||||

Symbol | Value | Description | ||||||||

X | L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | Name |

N_{X} | 20,683 | 5834 | 21,915 | 5479 | 4850 | 1065 | 14,395 | 2948 | 902 | Size |

$kXext$ | 1600 | 1500 | 2100 | 1900 | 2000 | 1900 | 2900 | 2100 | Ext. in-degree per neuron | |

$Iext$ | $\tau syn\nu bgJ,\nu bg=8Hz$ | DC ampl. per ext. input | ||||||||

Connectivity | ||||||||||

C_{YX} | from X | |||||||||

L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | ||

toY | L23E | 0.101 | 0.169 | 0.044 | 0.082 | 0.032 | 0.0 | 0.008 | 0.0 | 0.0 |

L23I | 0.135 | 0.137 | 0.032 | 0.052 | 0.075 | 0.0 | 0.004 | 0.0 | 0.0 | |

L4E | 0.008 | 0.006 | 0.050 | 0.135 | 0.007 | 0.0003 | 0.045 | 0.0 | 0.0983 | |

L4I | 0.069 | 0.003 | 0.079 | 0.160 | 0.003 | 0.0 | 0.106 | 0.0 | 0.0619 | |

L5E | 0.100 | 0.062 | 0.051 | 0.006 | 0.083 | 0.373 | 0.020 | 0.0 | 0.0 | |

L5I | 0.055 | 0.027 | 0.026 | 0.002 | 0.060 | 0.316 | 0.009 | 0.0 | 0.0 | |

L6E | 0.016 | 0.007 | 0.021 | 0.017 | 0.057 | 0.020 | 0.040 | 0.225 | 0.0512 | |

L6I | 0.036 | 0.001 | 0.003 | 0.001 | 0.028 | 0.008 | 0.066 | 0.144 | 0.0196 | |

Connection parameters | ||||||||||

Symbol | Value | Description | ||||||||

J | 87.81 pA | Reference synaptic strength. All synapse weights are measured in units of J. | ||||||||

$\sigma J,rel$ | Relative width of synaptic strength distribution | |||||||||

3 | •for lognormal distribution | |||||||||

0.1 | •for Gaussian distribution | |||||||||

g_{YX} | Relative synaptic strength: | |||||||||

1 | $X\u2208{TC,L23E,L4E,L5E,L6E},$ | |||||||||

$\u22124$ | $X\u2208{L23I,L4I,L5I,L6I}$, except for: | |||||||||

2 | $(X,Y)=(L4E,L23E)$ | |||||||||

$\u22124.5$ | $(X,Y)=(L4I,L4E)$ | |||||||||

$dE$ | 1.5 ms | Mean excitatory spike transmission delay | ||||||||

$dI$ | 0.75 ms | Mean inhibitory spike transmission delay | ||||||||

$\sigma d,rel$ | 0.5 | Relative width (stdev/mean) of transmission delay distributions | ||||||||

Neuron model | ||||||||||

Symbol | Value | Description | ||||||||

$Rm$ | 40 MΩ | Membrane resistance | ||||||||

$Cm$ | 250 pF | Membrane capacitance | ||||||||

$\tau m$ | $RmCm$ (10 ms) | Membrane time constant | ||||||||

$EL$ | −65 mV | Resting potential | ||||||||

θ | −50 mV | Fixed firing threshold | ||||||||

$Vm(t=0)$ | $[\u221265,\u221250]mV$ | Uniformly distributed initial membrane potential | ||||||||

$Vreset$ | $EL$ | Reset potential | ||||||||

$\tau ref$ | 2 ms | Absolute refractory period | ||||||||

$\tau syn$ | 0.5 ms | Postsynaptic current time constant | ||||||||

Thalamocortical (TC) input | ||||||||||

Symbol | Value | Description | ||||||||

$\nu \xafTC$ | 30 s^{−1} | Mean firing rate per TC neuron | ||||||||

$\Delta \nu TC$ | 30 s^{−1} | Firing-rate modulation amplitude per TC neuron | ||||||||

$fTC$ | 15 Hz | Frequency of sinusoidal firing-rate modulation |

A | Global simulation parameters | |||||||||
---|---|---|---|---|---|---|---|---|---|---|

Symbol | Value | Description | ||||||||

T | 5200 ms | Simulation duration | ||||||||

dt | 0.1 ms | Temporal resolution | ||||||||

B | Point-neuron network | |||||||||

Populations and external input | ||||||||||

Symbol | Value | Description | ||||||||

X | L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | Name |

N_{X} | 20,683 | 5834 | 21,915 | 5479 | 4850 | 1065 | 14,395 | 2948 | 902 | Size |

$kXext$ | 1600 | 1500 | 2100 | 1900 | 2000 | 1900 | 2900 | 2100 | Ext. in-degree per neuron | |

$Iext$ | $\tau syn\nu bgJ,\nu bg=8Hz$ | DC ampl. per ext. input | ||||||||

Connectivity | ||||||||||

C_{YX} | from X | |||||||||

L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | ||

toY | L23E | 0.101 | 0.169 | 0.044 | 0.082 | 0.032 | 0.0 | 0.008 | 0.0 | 0.0 |

L23I | 0.135 | 0.137 | 0.032 | 0.052 | 0.075 | 0.0 | 0.004 | 0.0 | 0.0 | |

L4E | 0.008 | 0.006 | 0.050 | 0.135 | 0.007 | 0.0003 | 0.045 | 0.0 | 0.0983 | |

L4I | 0.069 | 0.003 | 0.079 | 0.160 | 0.003 | 0.0 | 0.106 | 0.0 | 0.0619 | |

L5E | 0.100 | 0.062 | 0.051 | 0.006 | 0.083 | 0.373 | 0.020 | 0.0 | 0.0 | |

L5I | 0.055 | 0.027 | 0.026 | 0.002 | 0.060 | 0.316 | 0.009 | 0.0 | 0.0 | |

L6E | 0.016 | 0.007 | 0.021 | 0.017 | 0.057 | 0.020 | 0.040 | 0.225 | 0.0512 | |

L6I | 0.036 | 0.001 | 0.003 | 0.001 | 0.028 | 0.008 | 0.066 | 0.144 | 0.0196 | |

Connection parameters | ||||||||||

Symbol | Value | Description | ||||||||

J | 87.81 pA | Reference synaptic strength. All synapse weights are measured in units of J. | ||||||||

$\sigma J,rel$ | Relative width of synaptic strength distribution | |||||||||

3 | •for lognormal distribution | |||||||||

0.1 | •for Gaussian distribution | |||||||||

g_{YX} | Relative synaptic strength: | |||||||||

1 | $X\u2208{TC,L23E,L4E,L5E,L6E},$ | |||||||||

$\u22124$ | $X\u2208{L23I,L4I,L5I,L6I}$, except for: | |||||||||

2 | $(X,Y)=(L4E,L23E)$ | |||||||||

$\u22124.5$ | $(X,Y)=(L4I,L4E)$ | |||||||||

$dE$ | 1.5 ms | Mean excitatory spike transmission delay | ||||||||

$dI$ | 0.75 ms | Mean inhibitory spike transmission delay | ||||||||

$\sigma d,rel$ | 0.5 | Relative width (stdev/mean) of transmission delay distributions | ||||||||

Neuron model | ||||||||||

Symbol | Value | Description | ||||||||

$Rm$ | 40 MΩ | Membrane resistance | ||||||||

$Cm$ | 250 pF | Membrane capacitance | ||||||||

$\tau m$ | $RmCm$ (10 ms) | Membrane time constant | ||||||||

$EL$ | −65 mV | Resting potential | ||||||||

θ | −50 mV | Fixed firing threshold | ||||||||

$Vm(t=0)$ | $[\u221265,\u221250]mV$ | Uniformly distributed initial membrane potential | ||||||||

$Vreset$ | $EL$ | Reset potential | ||||||||

$\tau ref$ | 2 ms | Absolute refractory period | ||||||||

$\tau syn$ | 0.5 ms | Postsynaptic current time constant | ||||||||

Thalamocortical (TC) input | ||||||||||

Symbol | Value | Description | ||||||||

$\nu \xafTC$ | 30 s^{−1} | Mean firing rate per TC neuron | ||||||||

$\Delta \nu TC$ | 30 s^{−1} | Firing-rate modulation amplitude per TC neuron | ||||||||

$fTC$ | 15 Hz | Frequency of sinusoidal firing-rate modulation |

The stationary thalamic Poisson input and cortico-cortical input to the microcircuit present in the original model of Potjans and Diesmann (2014) are here replaced by DC currents. DC input slightly increases the degree of synchrony (see e.g., Brunel (2000)), but retains network dynamics and firing rate distributions across populations as in Potjans and Diesmann (2014).

The network of Potjans and Diesmann (2014) shows slightly synchronous behavior due to the E–I network of layer 4 being close to the synchronous irregular (SI) regime (Brunel 2000; Bos et al. 2016). In order to reduce synchrony, we here increased the average synaptic weight from neurons in population L4I (inhibitory) to L4E (excitatory) neurons by 12.5%, resulting in attenuated oscillations in layer 4. Taking advantage of the fact that point-neuron networks are amenable for theoretical analysis, we derived modified weights based on predictions from dynamical mean-field theory applied to the microcircuit model (Bos et al. 2016). Moreover, we found that high-frequency network oscillations seen for Gaussian synaptic weight distributions are reduced when using lognormally distributed synaptic weights (Sarid 2007; Iyer et al. 2013; Teramae and Fukai 2014). This made the dynamics more similar to experimental observations (Song et al. 2005; Buzsáki and Mizuseki 2014), and we thus also chose this for our network. Henceforth, we refer to our modified network as the “reference network”. Modulated activity of each thalamo-cortical (TC) neuron in the external thalamic population was modeled as synchronous spikes or as independent non-stationary Poisson processes with sinusoidally oscillating rate profiles (cf. Eq. (18)).

### Populations of Multicompartment Model Neurons

Cancellation effects from positive and negative contributions to extracellular potentials and effects of intrinsic dendritic filtering can only be captured with spatially extended multicompartment neuron models (Einevoll et al. 2013). In the hybrid scheme, extracellular potentials are estimated from the spiking activity in the point-neuron network through synaptic activation of populations of multicompartment model neurons (“LFP generators”). In principle, these mutually unconnected model neurons mirror their network counterparts and receive inputs from exactly the same point neurons.

In addition to the description of the point-neuron network model, different types of spatial information are thus needed to predict LFPs. For one, detailed dendritic morphologies are required for each individual network population (Fig. 2). Furthermore, the positions of neurons and synaptic connections must also be specified, as well as the separation of network populations into morphologically distinct cell types (Figures 2 and 4).

Availability of detailed cell-type specific connectivity of neural circuits, especially including information about synapse positions, is limited due to the substantial experimental effort involved. However, several ongoing large-scale neuroscience projects (Kandel et al. 2013) address this issue, and detailed connectomes are beginning to become publicly available (Jiang et al. 2015; Reimann et al. 2015; Markram et al. 2015). In the present example application, we used the connection probabilities as given by Izhikevich and Edelman (2008) derived from Binzegger et al. (2004). Note that the point-neuron network connectivity was partially derived from the same data (Potjans and Diesmann 2014). Quantitative data were provided for the number of connections in 5 cortical layers (layer 1 (L1), layers 2 and 3 grouped into a joint layer 2/3 (L2/3), and layers 4 (L4), 5 (L5), and 6 (L6)) between 17 cortical cell types, cortico-cortical connections from other areas, and 2 TC relay cell types. We follow the nomenclature of Izhikevich and Edelman (2008), where $y=p23$ denotes pyramidal cell types in layer 2/3, $y=b23$ and $y=nb23$ basket interneurons and non-basket interneurons within the same layer, $y=ss4(L23)$ spiny stellate cells in layer 4 with targets mainly within layer 2/3, $y=p4$ layer 4 pyramidal cells and so forth. Out of the 17 covered intracortical cell types only the $y=nb1$ cell type is not associated with any point-neuron network population in our scheme. To account for the lack of layer 1 neurons in our model, we renormalized the connection probabilities for the remaining 16 cortical cell types including the 2 TC relay-cell types, such that the occurrences *F*_{y} of all cell types *y* summed to 100% as given in Table 8.

Furthermore, we assumed that the excitatory point-neuron network populations within one layer are composed of pyramidal cells and spiny stellate cells if both are present in the layer, and that inhibitory network populations encompass both types of interneurons. This results in the grouping of cell types *y* into postsynaptic populations *Y* illustrated in Figure 2. The neuron count *N*_{y} of each cell type is then trivially computed from the frequency of occurrence *F*_{y} as given in Table 8 and Figure 2.

Inclusion of cell-type and layer-specific connections in the present hybrid scheme has some implications for how we proceed with setting up equivalent populations consisting of morphologically detailed model neurons. Different cell types belonging to a particular population may have different spatial distributions of synapses, or the populations may consist of different morphological classes of neurons (Kisvárday and Eysel 1992; Nowak 2003; Stepanyants et al. 2008). An example is layer 4 in which spiny stellate cells lack apical dendrites, while pyramidal cells have apical dendrites extending into layer 1. To incorporate some of this morphological diversity, we considered altogether 16 cell types for the 8 cortical network populations.

For each of the 16 cell types, we acquired representative morphological reconstructions of predominantly cat visual cortex neurons from several sources (Contreras et al. 1997; Kisvárday and Eysel 1992; Mainen and Sejnowski 1996; Stepanyants et al. 2008) (cf. Fig. 2, Table 7). Morphology files were obtained either from NeuroMorpho.org (Ascoli et al. 2007) or through personal communication with the authors. Constrained by layer boundary depths (Stepanyants et al. 2008) and laminar connectivities (see below) we applied an intermediate preprocessing step to our pyramidal cell morphologies. Assuming that the soma compartments of each cell type were centered in their corresponding layer, and noting that the layer-specific connectivity (cf. Table 8) implies connections to layer 1, we stretched the apical dendrites along the axis perpendicular to the cortical surface such that they reached the pial surface. The only exception was the p6(L4) morphology, which we extended to reach the center of layer 2/3 in accordance with Stepanyants et al. (2008) and the observation that Table 8 predicts zero connections within layer 1 and very few connections in layer 2/3 to the p6(L4) cell type. Due to lack of available morphologies of sufficient reconstruction quality, certain cell types were represented by the same neuron morphology. Interneuron types and spiny stellate cells in a given layer shared morphologies, the same interneuron morphology was reused in both layers 5 and 6, and finally the p5(L23) and p6(L56) cell morphologies were similar except for the stretching of the apical dendrites.

Morphology files | ||||
---|---|---|---|---|

Cell type y | Morphology $My$ | File | Source | Online source |

p23 | p23 | oi24rpy1.hoc | (Kisvárday and Eysel 1992) | #NMO_00851 (#NMO_10045) |

b23 | i23 | oi38lbc1.hoc | (Stepanyants et al. 2008) | — |

nb23 | i23 | oi38lbc1.hoc | (Stepanyants et al. 2008) | — |

p4 | p4 | oi53rpy1.hoc | (Kisvárday and Eysel 1992) | #NMO_00855 (#NMO_10040) |

ss4(L23) | ss4 | j7_L4ste.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 (#NMO_00905) |

ss4(L4) | ss4 | j7_L4ste.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 (#NMO_00905) |

b4 | i4 | oi26rbc1.hoc | (Stepanyants et al. 2008) | — |

nb4 | i4 | oi26rbc1.hoc | (Stepanyants et al. 2008) | — |

p5(L23) | p5v1 | oi15rpy4.hoc | (Kisvárday and Eysel 1992) | #NMO_00850 (#NMO_10046) |

p5(L56) | p5v2 | j4a.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 |

b5 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

nb5 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

p6(L4) | p6 | 51–2a.CN.hoc | (Contreras et al. 1997) | #NMO_00879 |

p6(L56) | p5v1 | oi15rpy4.hoc | (Kisvárday and Eysel 1992) | #NMO_00850 (#NMO_10046) |

b6 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

nb6 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

Morphology files | ||||
---|---|---|---|---|

Cell type y | Morphology $My$ | File | Source | Online source |

p23 | p23 | oi24rpy1.hoc | (Kisvárday and Eysel 1992) | #NMO_00851 (#NMO_10045) |

b23 | i23 | oi38lbc1.hoc | (Stepanyants et al. 2008) | — |

nb23 | i23 | oi38lbc1.hoc | (Stepanyants et al. 2008) | — |

p4 | p4 | oi53rpy1.hoc | (Kisvárday and Eysel 1992) | #NMO_00855 (#NMO_10040) |

ss4(L23) | ss4 | j7_L4ste.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 (#NMO_00905) |

ss4(L4) | ss4 | j7_L4ste.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 (#NMO_00905) |

b4 | i4 | oi26rbc1.hoc | (Stepanyants et al. 2008) | — |

nb4 | i4 | oi26rbc1.hoc | (Stepanyants et al. 2008) | — |

p5(L23) | p5v1 | oi15rpy4.hoc | (Kisvárday and Eysel 1992) | #NMO_00850 (#NMO_10046) |

p5(L56) | p5v2 | j4a.hoc | (Mainen and Sejnowski 1996) | #MDB_2488 |

b5 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

nb5 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

p6(L4) | p6 | 51–2a.CN.hoc | (Contreras et al. 1997) | #NMO_00879 |

p6(L56) | p5v1 | oi15rpy4.hoc | (Kisvárday and Eysel 1992) | #NMO_00850 (#NMO_10046) |

b6 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

nb6 | i5 | oi15rbc1.hoc | (Stepanyants et al. 2008) | — |

Online source numbers #NMO_* refer to NeuroMorpho.org identifiers, #MDB_* refer to ModelDB identifiers.

Postsynaptic | Layer | Occurr. | Num. | Presynaptic populations X and cell types x | ||||||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

pop. | cell type | syn. | L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | |||||||||||

p23 | b23 | nb23 | ss4(L4) | ss4(L23) | p4 | b4 | nb4 | p5(L23) | p5(L56) | b5 | nb5 | p6(L4) | p6(L56) | b6 | nb6 | TCs | TCn | |||||

Y | y | L | $Fy$ (%) | $kyL$ | $pyxL$ (%) | |||||||||||||||||

L23E | p23 | 2/3 | 26.7 | 5800 | 59.9 | 9.1 | 4.4 | 0.6 | 6.9 | 7.7 | — | 0.8 | 7.4 | — | — | — | 2.3 | — | — | 0.8 | — | — |

1 | 1306 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L23I | b23 | 2/3 | 3.2 | 3854 | 51.6 | 10.6 | 3.4 | 0.5 | 5.8 | 6.6 | — | 0.8 | 6.3 | — | — | — | 2.1 | — | — | 0.7 | — | 0.5 |

nb23 | 2/3 | 4.3 | 3307 | 48.6 | 11.4 | 3.3 | 0.5 | 5.5 | 6.2 | — | 0.8 | 5.9 | — | — | — | 1.8 | — | — | 0.6 | — | 0.7 | |

L4E | ss4(L4) | 4 | 9.4 | 5792 | 2.7 | 0.2 | 0.6 | 11.9 | 3.7 | 4.1 | 7.1 | 2 | 0.8 | 0.1 | — | — | 32.7 | — | — | 5.8 | 1.7 | 1.3 |

ss4(L23) | 4 | 9.4 | 4989 | 5.6 | 0.4 | 0.8 | 11.3 | 3.8 | 4.3 | 7.2 | 2.1 | 1.1 | 0.1 | — | — | 31.1 | — | — | 5.5 | 1.7 | 1.3 | |

p4 | 4 | 9.4 | 5031 | 4.3 | 0.2 | 0.6 | 11.5 | 3.6 | 4.2 | 7.2 | 2.1 | 1.2 | 0.1 | — | — | 31.4 | 0.1 | — | 5.9 | 1.7 | 1.3 | |

2/3 | 866 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 806 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L4I | b4 | 4 | 5.5 | 3230 | 5.8 | 0.5 | 0.8 | 11 | 3.8 | 4.2 | 8.4 | 2.4 | 1.1 | — | — | — | 30.3 | — | — | 5.4 | 1.6 | 1.2 |

nb4 | 4 | 1.5 | 3688 | 2.7 | 0.2 | 0.6 | 11.7 | 3.6 | 4 | 8.2 | 2.3 | 0.8 | 0.1 | — | — | 32.2 | — | — | 5.7 | 1.7 | 1.3 | |

L5E | p5(L23) | 5 | 4.8 | 4316 | 45.9 | 1.8 | 0.3 | 3.3 | 2 | 7.5 | — | 0.9 | 11.7 | 1 | 0.8 | 1.1 | 2.3 | 2.1 | — | 11.5 | 0.1 | 0.4 |

4 | 283 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 412 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 185 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

p5(L56) | 5 | 1.3 | 5101 | 44.3 | 1.7 | 0.2 | 3.2 | 2 | 7.3 | — | 0.8 | 11.3 | 1.2 | 0.8 | 1.1 | 2.3 | 2.5 | 0.3 | 11.3 | 0.2 | 0.5 | |

4 | 949 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 1367 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 5658 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L5I | b5 | 5 | 0.6 | 2981 | 45.5 | 2.3 | 0.2 | 3.3 | 2 | 7.5 | — | 1.1 | 11.6 | 1 | 0.9 | 1.3 | 2.3 | 2 | — | 11.4 | 0.1 | 0.4 |

nb5 | 5 | 0.8 | 2981 | 45.5 | 2.3 | 0.2 | 3.3 | 2 | 7.5 | — | 1.1 | 11.6 | 1 | 0.9 | 1.3 | 2.3 | 2 | — | 11.4 | 0.1 | 0.4 | |

L6E | p6(L4) | 6 | 14.0 | 3261 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.3 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

5 | 1066 | 46.8 | 0.8 | 0.3 | 3.4 | 2.1 | 7.7 | — | 0.6 | 11.9 | 1 | 0.6 | 0.8 | 2.3 | 2.1 | — | 11.7 | 0.1 | 0.4 | |||

4 | 1915 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 121 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

p6(L56) | 6 | 4.6 | 5573 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.3 | 1.2 | 13.2 | 7.8 | 7.8 | 0.6 | 2.9 | |

5 | 257 | 46.8 | 0.8 | 0.3 | 3.4 | 2.1 | 7.7 | — | 0.6 | 11.9 | 1 | 0.6 | 0.8 | 2.3 | 2.1 | — | 11.7 | 0.1 | 0.4 | |||

4 | 243 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 286 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 62 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L6I | b6 | 6 | 2.0 | 3230 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.4 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

nb6 | 6 | 2.0 | 3230 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.4 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

Postsynaptic | Layer | Occurr. | Num. | Presynaptic populations X and cell types x | ||||||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

pop. | cell type | syn. | L23E | L23I | L4E | L4I | L5E | L5I | L6E | L6I | TC | |||||||||||

p23 | b23 | nb23 | ss4(L4) | ss4(L23) | p4 | b4 | nb4 | p5(L23) | p5(L56) | b5 | nb5 | p6(L4) | p6(L56) | b6 | nb6 | TCs | TCn | |||||

Y | y | L | $Fy$ (%) | $kyL$ | $pyxL$ (%) | |||||||||||||||||

L23E | p23 | 2/3 | 26.7 | 5800 | 59.9 | 9.1 | 4.4 | 0.6 | 6.9 | 7.7 | — | 0.8 | 7.4 | — | — | — | 2.3 | — | — | 0.8 | — | — |

1 | 1306 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L23I | b23 | 2/3 | 3.2 | 3854 | 51.6 | 10.6 | 3.4 | 0.5 | 5.8 | 6.6 | — | 0.8 | 6.3 | — | — | — | 2.1 | — | — | 0.7 | — | 0.5 |

nb23 | 2/3 | 4.3 | 3307 | 48.6 | 11.4 | 3.3 | 0.5 | 5.5 | 6.2 | — | 0.8 | 5.9 | — | — | — | 1.8 | — | — | 0.6 | — | 0.7 | |

L4E | ss4(L4) | 4 | 9.4 | 5792 | 2.7 | 0.2 | 0.6 | 11.9 | 3.7 | 4.1 | 7.1 | 2 | 0.8 | 0.1 | — | — | 32.7 | — | — | 5.8 | 1.7 | 1.3 |

ss4(L23) | 4 | 9.4 | 4989 | 5.6 | 0.4 | 0.8 | 11.3 | 3.8 | 4.3 | 7.2 | 2.1 | 1.1 | 0.1 | — | — | 31.1 | — | — | 5.5 | 1.7 | 1.3 | |

p4 | 4 | 9.4 | 5031 | 4.3 | 0.2 | 0.6 | 11.5 | 3.6 | 4.2 | 7.2 | 2.1 | 1.2 | 0.1 | — | — | 31.4 | 0.1 | — | 5.9 | 1.7 | 1.3 | |

2/3 | 866 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 806 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L4I | b4 | 4 | 5.5 | 3230 | 5.8 | 0.5 | 0.8 | 11 | 3.8 | 4.2 | 8.4 | 2.4 | 1.1 | — | — | — | 30.3 | — | — | 5.4 | 1.6 | 1.2 |

nb4 | 4 | 1.5 | 3688 | 2.7 | 0.2 | 0.6 | 11.7 | 3.6 | 4 | 8.2 | 2.3 | 0.8 | 0.1 | — | — | 32.2 | — | — | 5.7 | 1.7 | 1.3 | |

L5E | p5(L23) | 5 | 4.8 | 4316 | 45.9 | 1.8 | 0.3 | 3.3 | 2 | 7.5 | — | 0.9 | 11.7 | 1 | 0.8 | 1.1 | 2.3 | 2.1 | — | 11.5 | 0.1 | 0.4 |

4 | 283 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 412 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 185 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

p5(L56) | 5 | 1.3 | 5101 | 44.3 | 1.7 | 0.2 | 3.2 | 2 | 7.3 | — | 0.8 | 11.3 | 1.2 | 0.8 | 1.1 | 2.3 | 2.5 | 0.3 | 11.3 | 0.2 | 0.5 | |

4 | 949 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 1367 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 5658 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L5I | b5 | 5 | 0.6 | 2981 | 45.5 | 2.3 | 0.2 | 3.3 | 2 | 7.5 | — | 1.1 | 11.6 | 1 | 0.9 | 1.3 | 2.3 | 2 | — | 11.4 | 0.1 | 0.4 |

nb5 | 5 | 0.8 | 2981 | 45.5 | 2.3 | 0.2 | 3.3 | 2 | 7.5 | — | 1.1 | 11.6 | 1 | 0.9 | 1.3 | 2.3 | 2 | — | 11.4 | 0.1 | 0.4 | |

L6E | p6(L4) | 6 | 14.0 | 3261 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.3 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

5 | 1066 | 46.8 | 0.8 | 0.3 | 3.4 | 2.1 | 7.7 | — | 0.6 | 11.9 | 1 | 0.6 | 0.8 | 2.3 | 2.1 | — | 11.7 | 0.1 | 0.4 | |||

4 | 1915 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 121 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

p6(L56) | 6 | 4.6 | 5573 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.3 | 1.2 | 13.2 | 7.8 | 7.8 | 0.6 | 2.9 | |

5 | 257 | 46.8 | 0.8 | 0.3 | 3.4 | 2.1 | 7.7 | — | 0.6 | 11.9 | 1 | 0.6 | 0.8 | 2.3 | 2.1 | — | 11.7 | 0.1 | 0.4 | |||

4 | 243 | 2.8 | 0.1 | 0.7 | 12.2 | 3.8 | 4.2 | 5.2 | 1.5 | 0.8 | 0.1 | — | — | 33.7 | — | — | 5.9 | 1.8 | 1.4 | |||

2/3 | 286 | 63.1 | 5.1 | 4.1 | 0.6 | 7.2 | 8.1 | — | 0.6 | 7.8 | — | — | — | 2.5 | — | — | 0.8 | — | — | |||

1 | 62 | 6.3 | 0.1 | 1.1 | — | — | 0.1 | — | — | 0.1 | — | — | — | — | — | — | — | — | 4.1 | |||

L6I | b6 | 6 | 2.0 | 3230 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.4 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

nb6 | 6 | 2.0 | 3230 | 2.5 | 0.1 | 0.1 | 0.7 | 0.9 | 1.3 | — | 0.1 | 0.1 | 4.9 | — | 0.4 | 1.2 | 13.2 | 7.7 | 7.7 | 0.6 | 2.9 |

Preserving the laminar cell density under 1 mm^{2} surface area of the point-neuron network model, we created for each postsynaptic cell type *y* model populations where somas were assigned random locations in 3D within cylindrical slabs with radius $r=564\mu m$ and thickness $h=50\mu m$, each centered in their respective layer (illustrated in Fig. 1*D*, see also Fig. 6*E*). Regardless of the vertical offset of the soma of pyramidal cells, postsynaptic target dendrites were therefore present within the $\u223c80\mu m$ thick (Stepanyants et al. 2008) uppermost layer 1 except for cell type p6(L4). For simplicity, each cell type was represented by a single reconstructed morphology in the present application. The full specification of the populations is given in Table 3.

Each neuron is modeled using the multicompartmental, passive cable formalism (Rall 1964; 2009; De Schutter and Van Geit 2009), describing the changes in membrane voltage and the associated transmembrane currents throughout all parts of the neuron geometry (Table 4). We used (non-plastic) exponential current-based synapses as in the point-neuron network model (Table 2). Synapse locations were randomly assigned onto cell compartments assuming a probability proportional to the compartment's surface area divided by the total surface area of the same cell within the target layer. Tables 4–6 summarize parameters relevant for the synapse models, synapse locations, and passive parameters of the multicompartment models.

Multicompartment model neurons | ||
---|---|---|

Symbol | Value | Description |

$cm$ | 1.0 μFcm^{−2} | Membrane capacity |

$rm$ | $\tau m/cm$ | Membrane resistivity |

$ra$ | 150 Ωcm | Axial resistivity |

$EL$ | $EL$ | Passive leak reversal potential |

$Vinit$ | $EL$ | Membrane potential at t = 0 ms |

$\lambda f$ | 100 Hz | Frequency of AC length constant |

$\lambda d$ | 0.1 | Factor for d _lambda rule (Hines and Carnevale 2001) |

$\sigma e$ | $0.3Sm\u22121$ | Extracellular conductivity |

r | $10002/\pi \mu m$ | Population radius |

h | 50 μm | Soma layer thickness |

$ncontact$ | 16 | Number of electrode contacts |

$helec$ | 100 μm | Laminar-electrode intercontact distance |

$rcontact$ | 7.5 μm | Electrode contact-point radius |

Multicompartment model neurons | ||
---|---|---|

Symbol | Value | Description |

$cm$ | 1.0 μFcm^{−2} | Membrane capacity |

$rm$ | $\tau m/cm$ | Membrane resistivity |

$ra$ | 150 Ωcm | Axial resistivity |

$EL$ | $EL$ | Passive leak reversal potential |

$Vinit$ | $EL$ | Membrane potential at t = 0 ms |

$\lambda f$ | 100 Hz | Frequency of AC length constant |

$\lambda d$ | 0.1 | Factor for d _lambda rule (Hines and Carnevale 2001) |

$\sigma e$ | $0.3Sm\u22121$ | Extracellular conductivity |

r | $10002/\pi \mu m$ | Population radius |

h | 50 μm | Soma layer thickness |

$ncontact$ | 16 | Number of electrode contacts |

$helec$ | 100 μm | Laminar-electrode intercontact distance |

$rcontact$ | 7.5 μm | Electrode contact-point radius |

Values for $\tau mandEL$ are inherited from network parameters in Table 5.

### Spatial Synaptic Connectivity

A full description of the connectivity in networks of multicompartment model neurons requires a 3D representation, for example in the form of sparse $NX\xd7NY\xd7ncomp$ matrices of synaptic weights and spike-transmission delays between presynaptic neurons $i\u2208[1,NX]$ and compartments $n\u2208[1,ncomp]$ of postsynaptic cells $j\u2208[1,NY]$. Here, *N*_{X} and *N*_{Y} denote the number of presynaptic and postsynaptic neurons in populations *X* and *Y*, respectively, and $ncomp$ the number of compartments of the postsynaptic cell.

In point-neuron networks, in contrast, connectivity is by definition only 2-dimensional (2D) as the cell morphology is collapsed into a single point and, consequently, the specificity of synapse locations on the postsynaptic morphology is ignored. In the proposed hybrid modeling scheme, the connectivity within the point-neuron network is consistent with the connectivity between point neurons and multicompartment model neurons. Ideally, each multicompartment model neuron has its equivalent in the point-neuron network and receives inputs from exactly the same presynaptic sources as its point-neuron counterpart. Synapses should be positioned on the dendritic tree according to anatomical data, and synaptic weights and time constants should be adapted such that the somatic membrane potential or somatic current match the point-neuron counterparts. Such mapping between point neurons and passive multicompartment neurons is feasible (Koch and Poggio 1985; Wybo et al. 2013; 2015).

In the current application of the hybrid scheme to the cortical microcircuit model, we make the simplest approximation to the mapping problem and fixed the current amplitudes $Iji,max$ and synaptic time constants as in the network model, with compartment specificity of connections dependent on compartment surface area (see Table 4). We further preserve only the statistics of connections (average number of inputs, distribution of spike-transmission delays) for each pair of presynaptic and postsynaptic neuron populations, exploiting that connections between network populations are drawn randomly with fixed probabilities. Finally, we simplify the positioning of synapses to a layer specificity of connections. The activation times of each synapse are then given by the spike train of a randomly drawn point neuron in the network model, with random delays consistent with the delay distribution in the network (Tables 2 and 4).

We first show how to derive a 2D point-neuron connectivity from a given 3D multicompartment-neuron connectivity and describe the case where the complexity of the point-neuron network is further reduced by pooling cell types. Then we describe the opposite procedure, connecting an existing (published) point-neuron network with a predefined 2D connectivity to a population of multicompartment model neurons such that the resulting 3D connectivity is as consistent as possible with anatomical data sets accounting for the compartment specificity of connections (e.g., the layer specificity of connections as in the anatomical data published by Binzegger et al. 2004). The procedures outlined below allow a reduction of complexity within the point-neuron network while accounting for the full diversity in cell types and synapse locations for multicompartment-neuron populations which is essential for predicting extracellular potentials (Fig. 3).

#### Construction of Point-Neuron Network Connectivity

For our example point-neuron network model, the cortical microcircuit model by Potjans and Diesmann (2014), the connectivity is to a large extent based on anatomical data from cat visual cortex (Binzegger et al. 2004; Izhikevich and Edelman 2008) (cf. Table 8).

From Table 8 we obtain 1) the number *N*_{y} of neurons belonging to cell type *y*, 2) the average total number *k*_{yL} of synapses on all compartments in layer *L* (input layer) of a single postsynaptic neuron of type *y*, and 3) the fraction *p*_{yxL} of the *k*_{yL} synapses formed with presynaptic neurons of cell type *x*. The quantity

*x*and a single postsynaptic cell of type

*y*in input layer

*L*(cf. network connectivity in Izhikevich and Edelman 2008). The number of synapses between all neurons in

*x*and all neurons in

*y*, irrespective of the input layer

*L*, is given by

The number *K*_{yx} of connections in combination with a chosen connectivity model (e.g., random graphs with binomially distributed (Erdős and Rényi 1959) or fixed in-/out-degree (Newman 2003) or random graphs with defined higher order statistics (Song et al. 2005; Zhao et al. 2011)) is sufficient for setting up the point-neuron network. Assuming independently drawn synapses (allowing multiple connections between neurons), the probability *C*_{yx} of at least one connection between a neuron of type *x* and a neuron of type *y* can be obtained from *K*_{yx} as (Potjans and Diesmann 2014)

In our case, the point-neuron microcircuit model consists of excitatory and inhibitory populations $X,Y$ (see Tables 1 and 2) pooling different presynaptic and postsynaptic cell types $x\u2208X$ and $y\u2208Y$ (cf. Fig. 4). Given a single multicompartment model neuron of type *y*, we compute the number *k*_{yXL} of incoming connections (in-degree) from cell types *x* in each presynaptic population *X* in a given layer *L* by pooling all connections as illustrated in Figure 4*A* as

The total number of connections onto postsynaptic cells *y* from cells in *X* is then

*C*

_{yXL}(Fig. 5

*B*) can be derived from Equation (5) analogous to Equation (3) for a presynaptic population size

*N*

_{X}(here, $NX=\u2211x\u2208XNx$).

In order to obtain the connectivity within the point-neuron network, that is, between populations *X* and *Y*, we also need to pool over all synapses of input layers *L* and cell types *y* within the postsynaptic population *Y* (dashed/dotted lines in Fig. 4*B*). Thus

*C*

_{YX}(cf. Eq. 3, Fig. 5

*A*).

#### From Pooled to Specific Network Connectivity

In case of an already existing point-neuron network model, the reverse task of creating a spatial connectivity *C*_{yXL} from a given point-neuron network connectivity *C*_{YX} is necessary. This inverse procedure compared with pooling over cell types and input layers entails introducing the cell-type specificity

*X*and

*Y*that are formed with a specific postsynaptic cell type

*y*(Fig. 5

*C*), and the layer specificity of connections

*X*and all cells of cell type

*y*formed in a particular layer

*L*(Fig. 5

*D*). The product $TyXLyXL$ defines the probability of a synapse between populations

*X*and

*Y*formed with a specific postsynaptic cell type

*y*in a particular layer

*L*(Fig. 4

*B*). Thus, if

*K*

_{YX}is given, the total number of connections in layer

*L*onto postsynaptic cells

*y*from cells in

*X*is

If *K*_{YX} is constructed from the same data as $TyX$ and $LyXL$, Equations (7–9) are fully consistent. However, *K*_{YX} can also be computed from any given point-neuron network connectivity *C*_{YX}. This is particularly relevant for the network connectivity *C*_{YX} (Fig. 5*A*) of Potjans and Diesmann (2014) that includes additional data sets for which spatial information on synapse locations is not available. Here the number of synapses *k*_{yXL} from population *X* established in layer *L* on each multicompartment model neuron of type *y* is obtained from Equation (9) as $kyXL=KyXL/Ny$.

### Forward Modeling of Extracellular Potentials

The LFP signal reflects transmembrane currents weighted according to the distance from the source to the measurement location (Einevoll et al. 2013), and here we compute the LFP from the model neurons using a now well-established forward-modeling scheme combining multicompartment neuron modeling and electrostatic (volume-conduction) theory (Holt and Koch 1999; Gold et al. 2006; Pettersen et al. 2008; Lindén et al. 2010; 2011; Reimann et al. 2013; Lindén et al. 2014; Tomsett et al. 2014).

Each morphology was spatially discretized into compartments using the

*f*= 100 Hz.

In this forward modeling scheme, localized synaptic activation of a morphologically detailed neuron results in spatially distributed transmembrane currents across the neuronal membrane as calculated using standard cable theory, see, for example, De Schutter and Van Geit (2009). The extracellular potentials, including the LFP, are in turn given as a weighted sum of transmembrane currents as described by volume conduction theory (Holt and Koch 1999; Einevoll et al. 2013).

The cable-equation description is summarized in box D in Table 4. Equation (27) relates synaptic-input currents *I*_{jin} onto compartment *n* in a neuron *j* from presynaptic neurons *i*, the membrane voltages $Vmjn$ and transmembrane currents $Imjn$, and is derived from the assumption (Kirchoff's law) of current balance in the intracellular node of the equivalent electrical circuit of a cylindrical compartment *n* with *m* neighboring compartments.

We use the standard convention that a positive membrane current is a positive current from the intracellular to the extracellular space across the membrane. $Imjn$ is assumed to be homogeneously distributed across the outer surface of the cylindrical compartment, and the calculation assumes that the electrical potential on the outside boundary of the membrane is zero at all times. Hence, there are no mutual interactions (i.e., ephaptic coupling (Anastassiou et al. 2011)) between the extracellular potential estimated using volume conduction theory, the transmembrane currents, and intracellular potentials.

The associated extracellular potential resulting from the transmembrane currents is calculated based on volume conduction theory (Nunez and Srinivasan 2006; Einevoll et al. 2013). In the present application where the signal frequencies are well below 1000 Hz, this calculation is simplified by applying the quasistatic approximation to Maxwell's equations, that is, terms with time derivatives of the electrical and magnetic fields are omitted, cf. (Hämäläinen et al. 1993). Furthermore, we assume the extracellular medium to be linear, isotropic, homogeneous, and ohmic (Pettersen et al. 2012; Einevoll et al. 2013) and represented by a scalar extracellular conductivity $\sigma e$.

Given a time-varying point current source with magnitude *I*(*t*) at position $r\u2032$, the scalar extracellular potential $\varphi (r,t)$ at position **r** and time *t* is then given by (Nunez and Srinivasan 2006; Lindén et al. 2014)

Contributions to the extracellular potential from multiple current sources, that is, transmembrane currents of all individual compartments *n* from all cells *j* in a population of *N* cells sum linearly. In accordance with the assumed homogeneous current distribution along each cylindrical compartment, the *line-source* approximation is used for dendritic compartments (Holt and Koch 1999). The line-source forward-modeling formula is obtained by integrating Equation (10) along the cylindrical axis of each compartment *n*, and summing the contributions from all $ncomp$ compartments (Holt and Koch 1999; Pettersen et al. 2008; Lindén et al. 2014):

Presently, we approximate the thick soma compartments as spherical current sources, and thus combine the point-source equation (Eq. 10) with the line-source formula (Eq. 11) for dendrite compartments, obtaining (Lindén et al. 2014):

Here, $\Delta sjn$ denotes compartment length, $r\u22a5jn$ the perpendicular distance from the electrode point contact to the axis of the line compartment, *h*_{jn} the longitudinal distance measured from the start of the compartment, and $ljn=\Delta sjn+hjn$ the longitudinal distance from the other end of the compartment. If the distance between electrode contacts and dendritic current sources becomes smaller than the radius of the dendritic segment, an unphysical singularity in our extracellular potential may occur. In these cases, singularities are avoided by setting $\u2223r\u2212rj,soma\u2223$ or $r\u22a5jn$ equal to the compartment radius.

Electrode contacts of real recording devices have finite spatial extent and are not point contacts as assumed above. However, the recorded signal can be well approximated as the mean of the potential averaged across the uninsulated surface (Robinson 1968; Nelson et al. 2008; Nelson and Pouget 2010; Ness et al. 2015), at least for current sources positioned further away than an electrode radius or so (Ness et al. 2015). Here we employed the *disc-electrode* approximation to the potential (Camuñas Mesa and Quiroga 2013; Lindén et al. 2014; Ness et al. 2015):

We further considered circular electrode contacts with a radius of $rcontact=7.5\mu m$, and we averaged the point-contact potential in Equation (12) over *m* = 50 random locations $rh$ across the contact surface *S*, *A*_{S} being the surface area. The chosen locations were distributed with uniform probability on circular discs representing each contact surface, with surface vectors oriented perpendicular to the electrode axis (Lindén et al. 2014). Calculations of extracellular potentials were facilitated by

### Data Analysis and Software

#### Model Measurements

The main simulation output of the hybrid scheme consists of spike trains of each neuron in the point-neuron network, “ground-truth” CSD and LFP of each neuron in the morphologically detailed postsynaptic model populations (see Tables 3 and 4). Here the term “ground-truth” refers to the fact that the CSD is computed from transmembrane currents rather than estimated from the LFP.

As the transmembrane current of each compartment (Eq. 27) is known at each simulation time step, we follow the procedure of Pettersen et al. (2008) to compute the “ground-truth” CSD in addition to the LFP. From *N* model neurons with $ncomp$ compartments having membrane currents $Imjn(t)$ and lengths $\Delta sjn$, we calculate the CSD signals $\rho (r,t)$ inside cylinder elements $V\rho (r)$ centered around each electrode contact as

*r*and heights equal to the electrode separation $helec$ (cf. Table 6).

In the present example application, the extracellular potential is computed at locations corresponding to a laminar multi-electrode array with 16 recording electrodes with an interelectrode distance of $helec=100\mu m$, positioned at the cylindrical axis of the model column with the topmost contact at the pial surface (cf. Fig. 1*D*, see Table 6 for details). Each electrode contact is set to have a radius of $7.5\mu m$ (cf. Eq. (13)). In the network, we also record membrane voltages and input currents from a subset of cells in each of the eight cortical populations (see Table 9).

Measurements and derived signals | ||
---|---|---|

Measurements | ||

Symbol | Description | Number of recorded units |

$Iiex(t)$ | Excitatory synaptic input current of neuron i | 100 per population X |

$Iiin(t)$ | Inhibitory synaptic input current of neuron i | 100 per population X |

V_{i}(t) | Somatic voltage of neuron i | 100 per population X |

s_{i}(t) | Spike train of neuron i | N_{X} neurons |

$\varphi i(r,t)$ | Single-cell LFP generated by neuron i | N_{X} neurons |

$\rho i(r,t)$ | Single-cell CSD generated by neuron i | N_{X} neurons |

Derived signals | ||

Symbol | Definition | Description |

I_{i}(t) | $Iiex(t)+Iiin(t)$ | Total synaptic input current of neuron i |

$I\xafXex(t)$ | $1nav\u2211i\u2208XnavIiex(t)$ | Average excitatory synaptic input current of population X ($nav=100$) |

$I\xafXin(t)$ | $1nav\u2211i\u2208XnavIiin(t)$ | Average inhibitory synaptic input current of population X ($nav=100$) |

$I\xafX(t)$ | $1nav\u2211i\u2208XnavIi(t)$ | Average total synaptic input current of population X ($nav=100$) |

$V\xafX(t)$ | $1nav\u2211i\u2208XnavVi(t)$ | Average membrane voltage of population X ($nav=100$) |

$\nu X(t)$ | $nXs(t)tbin$ | Instantaneous population (X) firing rate, with $nXs(t)$ being the number of spikes in $(t,t+tbin)$ of all cells in population X, $tbin=1ms$ |

$\nu \xafX(t)$ | $\nu X(t)NX$ | Average instantaneous firing rate of population X |

$\varphi X(r,t)$ | $\u2211i\u2208X\varphi i(r,t)$ | Population LFP of population X |

$\rho X(r,t)$ | $\u2211i\u2208X\rho i(r,t)$ | Population CSD of population X |

$\varphi (r,t)$ | $\u2211X\varphi X(r,t)$ | Compound LFP of all cells |

$\rho (r,t)$ | $\u2211X\rho X(r,t)$ | Compound CSD of all cells |

Rescaled signals | ||

$\varphi \gamma (r,t)$ | $\u2211X\u2211i\u2208X\u2032\u2282X\varphi i(r,t)$ | Compound LFP signal from subset of neurons $i\u2208X\u2032\u2282X$ with $NX\u2032=\gamma NX$, $\gamma \u2208[0,1]$ |

$\varphi \gamma \xi (r,t)$ | $\xi \varphi \gamma (r,t)(31)$ | Compound LFP signal from subset of neurons $i\u2208X\u2032\u2282X$ with $NX\u2032=\gamma NX$, $\gamma \u2208[0,1]$, rescaled by factor ξ |

Measurements and derived signals | ||
---|---|---|

Measurements | ||

Symbol | Description | Number of recorded units |

$Iiex(t)$ | Excitatory synaptic input current of neuron i | 100 per population X |

$Iiin(t)$ | Inhibitory synaptic input current of neuron i | 100 per population X |

V_{i}(t) | Somatic voltage of neuron i | 100 per population X |

s_{i}(t) | Spike train of neuron i | N_{X} neurons |

$\varphi i(r,t)$ | Single-cell LFP generated by neuron i | N_{X} neurons |

$\rho i(r,t)$ | Single-cell CSD generated by neuron i | N_{X} neurons |

Derived signals | ||

Symbol | Definition | Description |

I_{i}(t) | $Iiex(t)+Iiin(t)$ | Total synaptic input current of neuron i |

$I\xafXex(t)$ | $1nav\u2211i\u2208XnavIiex(t)$ | Average excitatory synaptic input current of population X ($nav=100$) |

$I\xafXin(t)$ | $1nav\u2211i\u2208XnavIiin(t)$ | Average inhibitory synaptic input current of population X ($nav=100$) |

$I\xafX(t)$ | $1nav\u2211i\u2208XnavIi(t)$ | Average total synaptic input current of population X ($nav=100$) |

$V\xafX(t)$ | $1nav\u2211i\u2208XnavVi(t)$ | Average membrane voltage of population X ($nav=100$) |

$\nu X(t)$ | $nXs(t)tbin$ | Instantaneous population (X) firing rate, with $nXs(t)$ being the number of spikes in $(t,t+tbin)$ of all cells in population X, $tbin=1ms$ |

$\nu \xafX(t)$ | $\nu X(t)NX$ | Average instantaneous firing rate of population X |

$\varphi X(r,t)$ | $\u2211i\u2208X\varphi i(r,t)$ | Population LFP of population X |

$\rho X(r,t)$ | $\u2211i\u2208X\rho i(r,t)$ | Population CSD of population X |

$\varphi (r,t)$ | $\u2211X\varphi X(r,t)$ | Compound LFP of all cells |

$\rho (r,t)$ | $\u2211X\rho X(r,t)$ | Compound CSD of all cells |

Rescaled signals | ||

$\varphi \gamma (r,t)$ | $\u2211X\u2211i\u2208X\u2032\u2282X\varphi i(r,t)$ | Compound LFP signal from subset of neurons $i\u2208X\u2032\u2282X$ with $NX\u2032=\gamma NX$, $\gamma \u2208[0,1]$ |

$\varphi \gamma \xi (r,t)$ | $\xi \varphi \gamma (r,t)(31)$ | Compound LFP signal from subset of neurons $i\u2208X\u2032\u2282X$ with $NX\u2032=\gamma NX$, $\gamma \u2208[0,1]$, rescaled by factor ξ |

We ran our simulations for a total duration of *T*=5200 ms using a temporal resolution of d*t*=0.1 ms (cf. Tables 5 and 11). However, LFP and CSD signals were resampled prior to file storage to a temporal resolution of $dt\psi =1ms$ by 1) applying a fourth-order Chebyshev type I low-pass filter with critical frequency $fc=400Hz$ and 0.05 dB ripple in the passband using a forward–backward linear filter operation, and 2) then selecting every 10th time sample.

#### Post-processing and Data Analysis

As the contributions to the CSD and LFP of the different cells sum linearly (cf. Eq. 12), we compute population-resolved signals as the sum over contributions from all cells in a population, and the full compound signals as the sum over all population signals (see Table 9). In Results, we derive a rescaled “low-density predictor” $\varphi \gamma \xi (r,t)$ of the LFP from random subsets of neurons in all populations. Thereby, we make a downscaled LFP-generating model setup with the same column volume, but with neuron density reduced to a factor $\gamma \u2208(0,1)$ of the original density, while preserving the in-degrees, that is, the number of synaptic connections onto individual neurons. The LFP from the downscaled setup is multiplied by an overall scaling factor *ξ* chosen to roughly preserve the LFP from the full-scale model. For analysis and plotting, the initial 200 ms of results after simulation onset was removed, and the signal mean was subtracted from LFP and CSD traces emulating DC filtering during experimental data acquisition. The contribution from each population to the overall LFP signal was assessed by computing and comparing the signal variances (see Table 10).

Cross-correlations between single-cell LFPs $\varphi i(r,t)$ were quantified by analyzing the power spectrum of the compound extracellular signals. Given that the compound LFP/CSD is a linear superposition of single-cell LFP/CSD contributions, the power spectrum of the compound signal $P\varphi (r,f)$ can be obtained as the sum of all single-cell power spectra $P\varphi i(r,f)$ and all pairwise cross-spectra $C\varphi i\varphi j(r,f)$, or equivalently from the average single-cell power spectrum $P\varphi \xaf(r,f)$, the average pairwise cross-spectrum $C\varphi \xaf(r,f)$, and the total cell count *N*, as shown in Table 10. Note that $C\varphi \xaf(r,f)$ and hence the average pairwise single-cell LFP coherence $\kappa \varphi \xaf(r,f)$ (Eq. 35) are real, while $C\varphi i\varphi j(r,f)$ is complex. Note that this definition allows for negative values of $\kappa \varphi \xaf(r,f)$. In the sum over all *i* and *j* in Equation (34), the imaginary parts of $C\varphi i\varphi j(r,f)$ and $C\varphi j\varphi i(r,f)$ cancel because $C\varphi i\varphi j(r,f)=C\varphi j\varphi i(r,f)*$ (∗ denotes the complex conjugate). The power spectrum $P\varphi \gamma \xi (r,f)$ of the compound signals of the “downscaled” network (i.e., low-density LFP predictor, see above) is given by the reduced cell count $\gamma N$, the average single-cell power spectrum $P\varphi \gamma \xi \xaf(r,f)$, and the average pairwise single-cell cross-spectrum $C\varphi \gamma \xi \xaf(r,f)$ calculated from that subset of neurons (see Table 10). These single-cell averages are the same as the respective single-cell averages of the full-scale model setup, apart from variability due to subsampling.

Throughout this paper, signal power spectra are estimated using Welch's average periodogram method (Welch 1967) (with the

Post-processing | ||
---|---|---|

Parameters | ||

Symbol | Value | Description |

$Ttrans$ | 200 ms | Start-up transient |

$dt\psi $ | 1 ms | Signal resolution |

Power spectral density settings | ||

Method | plt.mlab.psd * | Welch's average periodogram |

Symbol | Value | Description |

T_{ψ} | 5000 ms | Signal length |

NFFT | 256 | Number of data points used in each block for the FFT |

Fs | 1 kHz | Sampling rate |

noverlap | 128 | Number of overlapping data points between blocks |

window | plt.mlab.window _hanning * | Window filter (* plt denotes matplotlib ) |

Post-processing | ||
---|---|---|

Parameters | ||

Symbol | Value | Description |

$Ttrans$ | 200 ms | Start-up transient |

$dt\psi $ | 1 ms | Signal resolution |

Power spectral density settings | ||

Method | plt.mlab.psd * | Welch's average periodogram |

Symbol | Value | Description |

T_{ψ} | 5000 ms | Signal length |

NFFT | 256 | Number of data points used in each block for the FFT |

Fs | 1 kHz | Sampling rate |

noverlap | 128 | Number of overlapping data points between blocks |

window | plt.mlab.window _hanning * | Window filter (* plt denotes matplotlib ) |

#### The hybridLFPy Python Package

To facilitate usage of the hybrid scheme by other users, a novel

The source code release of

The present implementation of

## Results: LFP Generated by a Cortical Microcircuit

To illustrate the application of the hybrid scheme for predictions of LFPs from point-neuron networks, we here present results for a modified version of the point-neuron network model of Potjans and Diesmann (2014) with ∼78, 000 neurons mimicking a 1 mm^{2} patch of cat primary visual cortex (see Methods). The microcircuit model has realistic cell density and deliberately neglects many biological details on the single-cell level, focusing on the effect of the connectivity on the local network dynamics in such circuits. Despite this simplicity, the model displayed firing rates across populations in agreement with experimental observations (Potjans and Diesmann 2014), as well as propagation of spiking activity across layers. Likewise, the microcircuit model in conjunction with simplified, passive multicompartment populations is used here to study the effect of the (spatial) connectivity on the laminar pattern of spontaneous and stimulus-evoked CSD and LFP signals. For this, CSD and LFP signals for a laminar multi-electrode recording at different cortical depths are computed. The large network size of the model is further used to illustrate the effect of correlations and neuron density on CSD and LFP predictions.

### Spontaneous vs. Stimulus-Evoked LFP

We first consider the LFP generated by spontaneous network activity. The output of our hybrid scheme covers various scales and measurement modalities, from spikes of each neuron (Fig. 6*A*), population-averaged firing rates (Fig. 6*B*), excitatory, inhibitory, and total synaptic-input currents (Fig. 6*C*), and membrane voltages (Fig. 6*D*), to the compound CSD and LFP stemming from all populations of different cell types (Fig. 6*E*–*G*). For spontaneous activity, that is, no modulated thalamic input (cf. Table 1), we observe asynchronous irregular spiking in all populations (Fig. 6*A*) and firing rates similar to the original model (Potjans and Diesmann 2014) (Fig. 6*B*). In particular, layer 2/3 exhibits low firing rates, and in general inhibitory neurons fire with higher rates than excitatory neurons of the respective layers. The network is in a balanced regime (Brunel 2000), reflected by the substantial cancellation of population-averaged excitatory and inhibitory input currents (Fig. 6*C*). The population-averaged membrane potential fluctuates below the fixed firing threshold $\theta =\u221250mV$ down to $\u223c\u221280mV$ (Fig. 6*D*).

The corresponding compound CSD and LFP signals with contributions from all cortical populations are shown in Figure 6*F* and *G*, respectively. As expected for spontaneous cortical network activity, the LFP signal amplitudes are small, $\u22430.1mV$ (intriguingly close to what has been seen in experiments, see, e.g., Fig. 1 in Maier et al. 2010 for macaque visual cortex, Fig. 7 in Hagen et al. 2015 for mouse visual cortex), and exhibit strong cross-channel covariance and coherence in line with experimental observations, for example Einevoll et al. 2007; Maier et al. 2010; Riehle et al. 2013; Hagen et al. 2015. In the model, these correlations stem from dendritic cable properties and volume conduction effects (Pettersen et al. 2008; Lindén et al. 2011; Łe¸ski et al. 2013). The correlations across channels are, as expected, generally less visible in the more localized CSD signal where volume conduction effects are absent (Nicholson and Freeman 1975; Pettersen et al. 2006; 2008).

LFP and associated CSD studies have commonly been used to investigate stimulus-evoked responses in sensory cortices, see, for example, Mitzdorf and Singer (1979); Mitzdorf (1985); Di et al. (1990); Swadlow et al. (2002); Einevoll et al. (2007); Sakata and Harris (2009); Szymanski et al. (2009); Maier et al. (2010); Jin et al. (2011); Szymanski et al. (2011), as well as Einevoll et al. (2013) and references therein. To model the situation with a sharp onset of a visual stimulus (or direct electrical stimulation of the TC pathway (Mitzdorf and Singer 1979)), we drive the network with a short thalamic pulse mimicking a volley of incoming spikes onto primary visual cortex from the visual thalamus (lateral geniculate nucleus). The activation targets populations in layers 4 and 6 (see Fig. 1*A* or Fig. 5*A*,*B*) and propagates in the network to populations in layers 2/3 and 5 (Fig. 7*A*,*B*). At the level of spiking activity, the results match the behavior of the original model (Potjans and Diesmann 2014) and even agree qualitatively with experimental findings in rodents from stimulus-evoked activation of auditory cortex (Sakata and Harris 2009) and somatosensory cortex (Armstrong-James et al. 1992; Einevoll et al. 2007; Reyes-Puerta et al. 2015). This points to a general biological plausibility of this generic network model, based largely on data from cat V1.

The corresponding CSD and LFP profiles across depth associated with this spiking activity are determined by the synapse locations and dendritic filtering by cell-type specific morphologies (see Figs 3 and 5*D*). For the CSD (Fig. 7*C*), a complex alternating spatiotemporal pattern of current sinks (negative CSD) and sources (positive CSD) is observed. Due to volume conduction, this detailed spatial pattern is largely smeared out in the LFP profile (Fig. 7*D*), which displays a strong positivity across the middle layer around *t* = 910 ms. We also note the different spatiotemporal profiles of the spiking activity (Fig. 7*A*,*B*) compared with laminar CSD and LFP profiles. Not only do the CSD and LFP signals typically fade out 5–10 ms later than the spiking, the spatial profiles are also very different. For example, the LFP signal is very weak between channels 11 and 12 (i.e., between 1100 and $1200\mu m$ depth), even if the firing rate of layer 5 positioned between these channels is very high. We also note that the predicted LFP magnitudes are an order of magnitude larger compared with the LFP predicted for spontaneous activity, that is, ∼ 1 mV for stimulus-evoked versus ∼ 0.1 mV for spontaneous activity.

Although the present example has not been tuned to address specific experiments, we nevertheless observe that the model predictions display several features seen in experiments. For example, stimulus-evoked LFP amplitudes on the order of 1 mV are similar to the maximal amplitudes (∼ 1–3 mV) observed in cat V1 following electric stimulation of TC axons (optical radiation) (Mitzdorf 1985) and ventroposterior lateral pathway axons projecting to rat somatosensory cortex (Castro-Alamancos and Connors 1996, Fig. 3), and in rat somatosensory (barrel) cortex following whisker flicks (Di et al. 1990; Einevoll et al. 2007; Reyes-Puerta et al. 2015). Also some qualitative features of the spatiotemporal LFP and CSD patterns from the cat visual cortex experiments of Mitzdorf (1985) can be recognized. One example is the early CSD sink around layer 4 (i.e., at channel 6 close to the boundary between layers 2/3 and 4), another is the large positivity in LFP extending through most of the layers following the initial response to the thalamic input volley.

The observed time courses of the LFP and CSD in our simulations with *δ*-pulse thalamic activations are not as directly comparable with experimental stimulus-evoked activity, where the input is temporally filtered by several cell populations before reaching thalamus on the way to cortex. However, we note that very swift stimulus-evoked responses lasting not much longer than the ∼ 10 ms response volleys seen in our simulations also are observed in the somatosensory system (Di et al. 1990; Einevoll et al. 2007; Reyes-Puerta et al. 2015).

### Effect of Network Dynamics on LFP

The spiking activity in the microcircuit model is highly sensitive to modification of intrinsic model parameters and external input (Bos et al. 2016). In general LFPs reflect synaptic input both from local and distant neurons (Herreras et al. 2015), and also depend on network state (see, e.g., Kelly et al. 2010; Gawne 2010). With the hybrid scheme, we first illustrate as an example the dependence of the spontaneous LFP, that is, the LFP without thalamic input, on local network dynamics as determined by intrinsic network parameters. In particular, we compare our reference network model with the original model proposed by Potjans and Diesmann (2014) where the strength of connections from L4I to L4E neurons is weaker and the synaptic weights are drawn from a Gaussian distribution. We next investigate the LFP when the network is stimulated by sinusoidally modulated thalamic input.

For spontaneous activity (Fig. 8*A*–*E*), the spiking (Fig. 8*A*) is asynchronous irregular in all populations. The firing-rate power spectra (Fig. 8*B*) vary from relatively flat to more band-pass-like with a maximum power at ∼ 80 Hz. The suppressed power at lower frequencies arises from active decorrelation due to inhibitory feedback (Tetzlaff et al. 2012). In combination with the low-pass filtering involved in the generation of LFPs from spiking activity (Lindén et al. 2010; Łe¸ski et al. 2013), the firing-rate spectra translate into an LFP power spectrum that, depending on recording depth, has either low-pass or band-pass filter characteristics (Fig. 8*D*). The notably sharper attenuation of the LFP power spectra compared with the firing-rate power spectra above $\u2273100Hz$ is expectedly due to the intrinsic dendritic filtering effect (Lindén et al. 2010). As the effect of this filtering depends on the position of the electrode compared with the neuronal morphology (Lindén et al. 2010), the result is a variable LFP power spectral density (PSD) profile across cortical depths, even in the absence of any structured external input (Fig. 8*E*).

There is experimental evidence that cortical microcircuits receive oscillating input at various frequencies from remote areas and subcortical structures (Bastos et al. 2015, van Kerkoerle et al. 2014; Ito et al. 2014). To illustrate the effect of an oscillatory external input in our model, we model thalamic input as independent realizations of non-stationary Poisson processes with a sinusoidal rate profile (see Table 5, Fig. 8*F*–*J*). The spiking activity of all populations as well as the CSD and LFP across depth are sensitive to thalamic input, showing that the network response goes beyond the layers receiving thalamic input, that is, layers 4 and 6. The stimulus-evoked spiking activity follows the 15 Hz modulated rate of the thalamic input (Fig. 8*F*). The 15 Hz oscillation is reflected in the firing-rate spectrum as a peak at ∼15 Hz (Fig. 8*G*) and is robustly transferred to the LFP (Fig. 8*H*–*J*). The LFP oscillation strength varies with depth and is greatest in channels 1 and 2 and channels 8–14, while the oscillation is barely seen in channels 3–6. Populations L4E/I and L6E/I receive thalamic inputs around the depths of channels 8–10 and channels 14 and 15, and these channels are strongly affected by the stimulus. However, recurrent connections between populations, dendritic propagation of currents, and volume conduction produce strong LFP oscillations also in other channels.

The LFP amplitudes are not only influenced by the temporal structure of the external input, but also by synaptic weights in two ways: first via their influence on the spiking dynamics in the point-network simulation and second via the influence on the size of the synaptic currents setting up the transmembrane currents in the multicompartment models in the LFP-computing step. In order to illustrate the weight dependence, we compare the LFP under spontaneous activity in our reference network model (Fig. 8*A*–*E*) with the corresponding spontaneous LFP in the original model by Potjans and Diesmann (2014) (Fig. 8*K–O*). The lower inhibition from population L4I onto L4E and the narrow weight distribution compared with our model gives a higher degree of spike synchrony in all populations. Although our model exhibits asynchronous irregular activity (Fig. 8*A*), the original model is closer to an SI regime (Brunel 2000) (Fig. 8*K*), resulting in high-frequency oscillations ∼80 Hz and in the 300–400 Hz band (Fig. 8*L*), both associated with delay loops in the multilayered network (Bos et al. 2016). The 80 Hz oscillation also appears in the LFP and its corresponding power spectra (Fig. 8*M–O*), but the magnitude of the peak is not constant across depth. The lowest magnitudes are located in the vicinity of layers 2/3 and 5 (in channels 3 and 4 and channels 11–13).

### Contributions from Individual Populations to CSD and LFP

The direct interpretation of CSD and LFP signals in terms of the underlying activity of different populations or input pathways is inherently ambiguous and thus difficult: for example, a CSD sink observed in cortical layer 2/3 can alternatively stem from excitatory synaptic inputs to the basal dendrites of layer 2/3 cells, similar inputs into the apical dendrites of layer 5 cells, or even return currents from appropriately placed inhibitory inputs onto the same cells (Lindén et al. 2010, Einevoll et al. 2013). Several schemes for decomposition of CSD and LFP data into contributions from cortical populations have thus been proposed: principal component analysis (Di et al. 1990), laminar population analysis (LPA, Einevoll et al. 2007), and independent component analysis (Łe¸ski et al. 2010; Makarov et al. 2010; Głąbska et al. 2014; Herreras et al. 2015). In our model, we have the benefit of having the contributions from the various populations, connections, or different synapse types to the CSD and LFP signals directly accessible.

Here, we focus on the LFP and CSD contributions from individual populations and different synapse types. To quantify the contributions from each postsynaptic population, that is, the CSD and LFP stemming from the transmembrane currents of a population, we simply summed all single-cell CSD and LFP contributions from all neurons in the population. Results for spontaneous network activity are shown in Figure 9. As different cell types are assigned appropriate morphologies and cortical depths based on available anatomical data, a trivial consequence is that the neurons in each population make their main contributions to the CSD and LFP at depths spanned by their dendrites (Fig. 9*B*,*C* for post-synaptic populations L23E and L6E, respectively). For example, L6E neurons have a high density of afferent synapses on apical dendrites in layer 4, and as a consequence, the amplitude of CSDs and LFPs generated by population L6E (Fig. 9*C*) is large in the vicinity of layer 4 and not just in layer 6.

To quantify the relative contributions from the various populations, we show in Figure 9*D*,*E* the CSD and LFP variances across time (corresponding to PSDs summed over all nonzero frequencies) for all depths (Lindén et al. 2011, Łe¸ski et al. 2013). For the compound CSD, only the L23E neurons contribute substantially in the superficial channels (channels 1–4) (Fig. 9*D*). At deeper contacts, the main contributing population is L6E. L4E and L5E populations make sizable contributions only in the channels closest to their somatic location reflecting that the net associated return currents of their distributed synaptic inputs are largely restricted to somatic regions (Lindén et al. 2010). The bulk of the variance of the CSD in layer 5 arises in about equal parts from the L6E and relatively sparse L5E populations. The magnitude of the depth-resolved CSD variances of each layer's inhibitory population (L23I, L4I, ...) is consistently one order of magnitude or more smaller than that of the corresponding excitatory cell populations. Moreover, they span comparatively small depth ranges as determined by the maximum extent of the dendrites.

The depth-resolved LFP variance (Fig. 9*E*) has similar features as the CSD variance, as expected from their common biophysical origin. However, volume conduction has some qualitative effects, such as the reduced relative contribution from the L6E population in layer 2/3. As correlated synaptic inputs are known to amplify and increase the spread of the LFP generated from a cortical population (Lindén et al. 2011, Łe¸ski et al. 2013), this may reflect that the synaptic inputs to the L23E population are more correlated than to the L6E population.

Overall, we conclude that for the spontaneous activity in our network, the L23E and L6E populations dominate the compound LFP and CSD with only smaller contributions from the other excitatory populations. The signal variances from the inhibitory populations are typically much smaller than the contributions from the excitatory populations, suggesting that they can be safely neglected, in line with recent findings of Mazzoni et al. (2015).

Although transmembrane currents in inhibitory neurons provide little of the observed CSD and LFP, the inhibitory synaptic inputs onto excitatory neurons provide a substantial contribution. In the present network, inhibitory synaptic currents have a 4-fold larger amplitude compared with most excitatory synapses (Table 5), inhibitory neurons have higher overall firing rates compared with excitatory neurons (Potjans and Diesmann 2014), and inhibition specifically targets soma-proximal sections (Markram et al. 2004). Since our LFP-generating model is linear (passive cable formalism and linear synapse model), we can decompose the compound signal into contributions from each synapse type. For example, selective removal of either inhibitory or excitatory synaptic currents in the CSD and LFP modeling (Fig. 9*F*,*G*) shows that the CSD and LFP signals (Fig. 9*H*) are dominated by inhibitory synaptic currents and their associated return currents (Fig. 9*I*,*J*) when the network operates in the spontaneous asynchronous irregular firing regime. Furthermore, at most depths the variance of the signal arising from inhibitory input exceeds the compound variance, implying that the inhibitory component is generally negatively correlated with the excitatory component (Fig. 9*I*,*J*). Visual inspection of Figure 9*F*–*H* also reveals that the inhibitory dominance appears particularly strong at high frequencies, in accordance with the firing-rate PSDs in Figure 8*B* showing less power of inhibitory spiking, and thus inhibitory synaptic-input currents, at low frequencies.

The relative contribution from excitation and inhibition to the CSD and LFP depends on thalamic input and network state, however. For oscillatory thalamic input, the CSDs and LFPs from excitatory (Fig. 10*A*) or inhibitory synapses (Fig. 10*B*) alone show much stronger oscillations than the compound signals (Fig. 10*C*). This follows from the observations that the contributions to the CSD and LFP from excitatory and inhibitory synapses are anticorrelated, which, in turn, is a consequence of the dynamical balance between excitation and inhibition in asynchronous-irregular states of balanced random networks (Hertz 2010; Renart et al. 2010; Tetzlaff et al. 2012).

For transient thalamic activation, the same cancellation of excitation and inhibition can be observed (Fig. 10*F*–*H*), although it is more pronounced at some depths (layers 5 and 6, i.e., channels 11–13) than at others, depending on whether excitation and inhibition are in phase or not (Fig. 10*I*,*J*). For such strong and transient thalamic input, the network activity is briefly imbalanced as inhibition cannot keep up with thalamic excitation on the short-time scales. Since thalamic input is most prominent in layer 4 (channels 7–10), fewer cancellation effects are present here: in channel 9 in Figure 10*I* the total CSD variance is, for example, seen to be larger than the individual contributions from inhibitory and excitatory synapses.

### Effect of Input Correlations

Synaptic inputs to two neighboring cells are typically correlated because 1) they receive, to some extent, inputs from the same presynaptic sources (“shared-input correlation”), and 2) the spike trains of the presynaptic neurons may be correlated (“spike-train correlation”). The net synaptic-input correlation is determined by the interplay between these two contributions, shared-input correlations and spike-train correlations (Renart et al. 2010; Tetzlaff et al. 2012). As the LFP is largely generated by synaptic inputs, synaptic-input correlations result in correlated single-cell LFP contributions $\varphi i(r,t)$ (for details, see Lindén et al. 2011; Łe¸ski et al. 2013). As outlined in the following, these single-cell LFP correlations play a dominating role for the spectrum of the compound LFP.

The power spectrum $P\varphi (r,f)$ (Eq. 36) of the compound LFP $\varphi (r,t)=\u2211i=1N\varphi i(r,t)$ of a population of *N* neurons is given by

*N*, the second term is proportional to $N(N\u22121)\u2248N2$ for large

*N*. Hence, for large

*N*, even small cross-correlations may dominate the spectrum of the compound LFP. Here, we investigate this situation by calculating the power spectrum $P\varphi 0(r,f)$ of the compound LFP under the assumption of zero cross-correlation (where it simply reduces to a sum over single-cell spectra $P\varphi i(r,f)$), and compare to the true spectrum $P\varphi (r,f)$. The ratio between these quantities is given by

*N*. Note also that anti-correlated signals ($\kappa \varphi \xaf(r,f)<0$) may lead to a ratio $P\varphi (r,f)/P\varphi 0(r,f)<1$.

In the example application, both for spontaneous (Fig. 11*A*,*B*) and for evoked activity (Fig. 11*C*,*D*) the compound power spectra $P\varphi (r,f)$ are systematically (across channels and frequencies) larger than $P\varphi 0(r,f)$, demonstrating the importance of cross-correlations in the present network. Depending on the recording depth and frequency, the ratio varies from ∼ 1 to 10^{3} (see Fig. 11*B*,*D*). For spontaneous activity (see Fig. 8*A*–*E*), the largest effects of cross-correlations are typically found at higher frequencies (Fig. 11*A*,*B*). At low frequencies, cross-correlations are suppressed by inhibitory feedback (cf. Fig. 8*B*, Tetzlaff et al. (2012)). The thalamic sinusoidally modulated input to the network (Fig. 8*F*–*J*) synchronizes single-cell CSDs and LFPs at the stimulus frequency and gives a large boost of the LFP power at this frequency (see peak in Fig. 11*C*,*D*). Close inspection reveals that there is also in fact a boost of the power at ∼ 80 Hz, but much less so than at ∼ 15 Hz. Note that the external activation hardly affects the single-cell spectra (see red curves in Fig. 11*A*,*C*). LFP synchronization is thus mainly encoded in the phase of $\varphi i(r,t)$.

In conclusion, cross-correlations between single-cell LFP contributions play a pivotal role in shaping the compound LFP spectra (similar for CSD spectra, results not shown). To account for the dominant features of the LFP (CSD) in such models, it is therefore essential to include the main factors determining the synaptic-input correlations, that is, realistic correlations in presynaptic spike trains and shared-input structure. The findings presented here for the cortical microcircuit model hold in general in the presence of correlated activity. Only the details of the spectra depend on the specific underlying network dynamics.

### Network Downscaling

Due to the computational cost associated with modeling LFPs, it would be desirable to downsize the postsynaptic populations of multicompartment model neurons to a fraction $\gamma N$ ($\gamma \u2208(0,1)$) while leaving the point-neuron network at full size (*N*) and at the same time preserve the in-degrees, that is, the number of inputs, of each postsynaptic cell. The power spectrum of the full-scale LFP can indeed be estimated from the population-averaged single-cell power spectra $P\varphi \gamma \xaf(r,f)\u2248P\varphi \xaf(r,f)$ and coherences $\kappa \varphi \gamma \xaf(r,f)\u2248\kappa \varphi \xaf(r,f)$ computed for downsized networks by means of Equation (15). These quantities are preserved except for deviations due to smaller sampling size $\gamma N$ (Eq. 36 in Table 10). However, due to lack of phase information in the power spectra, one cannot estimate the LFP time course.

One could attempt to obtain a time-course estimate $\varphi \gamma \xi (r,t)$, that is, a “low-density LFP prediction”, of the full-scale signal $\varphi (r,t)$ by upscaling single-cell LFPs $\varphi i(r,t)$ computed in the downsized setup by a scalar factor *ξ* (cf. Eq. 31). Such a naive upscaling can grossly recover the amplitude of the full-scale LFP $\varphi (r,t)$, but it still only partially reconstructs its detailed time course (Fig. 12*A*,*E*). Also, this approach does not generally give accurate power spectra as the two terms in Equation (15) scale differently with *ξ*: The rescaling introduces a prefactor $\xi 2$ in the population-averaged single-cell power spectra $P\varphi \gamma \xi \xaf(r,f)\u2248\xi 2P\varphi \xaf(r,f)$, while the coherences $\kappa \varphi \gamma \xi \xaf(r,f)\u2248\kappa \varphi \xaf(r,f)$ are unchanged. Thus, the compound spectra $P\varphi (r,f)$ and $P\varphi \gamma \xi (r,t)$ of the full-size LFP and the low-density LFP predictor, respectively, differ. Their ratio

*ξ*which allows for the recovery of the full-size compound LFP power, that is, makes the ratio in Equation (17) equal to one for all spatial positions

**r**and frequencies

*f*. This can only be done in the special case where $\kappa \varphi \xaf(r,f$) is a constant

*c*. Here the two extreme cases correspond to no correlation ($\kappa \varphi \xaf(r,f)=0$ with $\xi =1/\gamma $) and full correlation between all single-cell signals ($\kappa \varphi \xaf(r,f)=1$ with $\xi =1/\gamma $).

The substantial scaling effects observed for our microcircuit model in the asynchronous state (Fig. 12*A*–*D*) suggest that correlations cannot be neglected even when modeling the LFP for spontaneous network activity. Choosing the scaling factor $\xi =1/\gamma $ corresponding to $\kappa \varphi \xaf(r,f)=0$ (red lines in Fig. 12*A*,*C*) leads to a severe underestimation of the full-size compound power spectrum (Fig. 12*C*,*D*). Even though the correlation (i.e., Pearson's correlation coefficients) between the full-size LFP signals and low-density LFP predictions are quite high (Fig. 12*B*), the power ratios (Fig. 12*D*) reveal that the rescaled signals are systematically wrong in frequency bands where single-cell LFPs are most strongly correlated (i.e., the frequencies for which the compound spectra are much larger than the predictions when omitting cross-correlations (Fig. 11*A*,*B*)). Assuming the full-correlation scaling factor $\xi =1/\gamma $, however, typically overestimates the full-size compound power spectrum (cf. gray spectra in Fig. 12*C*), particularly at low frequencies.

The results for the sinusoidally stimulated network (Fig. 12*E*–*H*) are quite similar to the spontaneous-activity results, except around the stimulation frequency 15 Hz where the modulated input leads to strongly correlated single-cell LFP contributions and a strong boost of the compound LFP. The approximate downscaling procedure assuming the full-correlation scaling factor $\xi =1/\gamma $ thus essentially agrees with the full-size compound spectrum for 15 Hz (while giving a strong overestimation for frequencies other than ∼ 15 and ∼ 80 Hz).

We note in passing that in contrast to the power spectra, the computation of the spike-trigger-averaged LFP (stLFP) (Swadlow et al. 2002; Nauhaus et al. 2009; Jin et al. 2011; Denker et al. 2011) in downsized networks do not have a principled problem due to cross-correlations between single-cell LFPs. As stLFPs are linearly dependent on the single-neuron LFP contributions, the only principled problem with downsizing is increased noise in the estimates due to sampling over fewer postsynaptic neurons.

### LFP Prediction from Population Firing Rates

An important question in systems neuroscience is to what extent the dynamics of networks of thousands or millions of neurons can be described by much simpler mathematical descriptions in terms of neural *populations* (Deco et al. 2008; Blomquist et al. 2009). Likewise, we here ask the question of whether LFPs can be predicted from knowledge of the population firing rate (Einevoll et al. 2007; Moran 2008; Einevoll et al. 2013). The hybrid scheme is excellently suited for testing and development of simplified numerical schemes for LFP prediction as the ground truth, that is, the LFP from the full network, is available as benchmarking data.

The use of current-based synapses and passive dendrites in the present application of the hybrid scheme renders synaptic events independent of each other in the LFP prediction. This inherent linearity results in a unique spatio-temporal relation $HXi(r,\tau )$ for $\tau \u2208[\u2212\u221e,\u221e]$ between a spike event of a point neuron *i* in population *X* and its contribution to the compound LFP $\varphi (r,t)$ from all its postsynaptic multicompartment model neurons. In this scheme, the link is causal, that is, the spikes drive the LFP, so that $HXi(r,\tau )=0$ for $\tau <0$ (as in LPA (Einevoll et al. 2007)). $HXi(r,\tau )$ encompasses connectivity, spike transmission delays and all postsynaptic responses including effects of synaptic-input currents and passive return currents. With a linear, current-based model such as our example cortical column, it is in principle possible by linear superposition to fully reconstruct the compound LFP if $HXi(r,\tau )$ and the spike times *t*^{i}_{l} are known for all neurons *i* in each population of the network.

It is, however, in the case of large networks impractical to assess each $HXi(r,\tau )$, as the LFP response needs to be determined for every neuron separately. In contrast, a large reduction in dimensionality can be achieved by determining the population-averaged LFP responses $H\xafX(r,\tau )$ of a spike within each population *X*. We thereby ignore heterogeneity in kernels $HXi(r,\tau )$ due to the variability in the connections from neurons in population *X*. An approximate compound LFP $\varphi *(r,t)$ based on population firing rates (Einevoll et al. 2007) can be computed from these extracted population kernels by means of the convolution $\varphi *(r,t)=\u2211X(\nu X*H\xafX)(r,t)$, where $\nu X(t)$ are the instantaneous population firing rates.

Here we estimate the population LFP kernels by computing the response to synchronous activation of all neurons in a population (Fig. 13). The spatio-temporal kernels $H\xafX(r,\tau )$ are extracted from time slices $[tX\u221220ms,tX+20ms]$ of the compound LFP response $\varphi (r,t)/NX$, where *N*_{X} is the number of neurons in a presynaptic population. The procedure results in unique kernels $H\xafX(r,\tau )$ for each excitatory and inhibitory population in the network (Fig. 13*A*).

In the example application, the population kernels $H\xafX(r,\tau )$ differ significantly between populations. Excitatory spike events result in prominent LFP negativities at depths where most connections are made, such as in layer 4 (channels 7–10) for TC connections ($H\xafTC(r,\tau )$, column 1 in Fig. 13*A*, cf. Fig. 5*C*). In contrast, spikes of inhibitory point neurons on average produce prominent LFP positivities in their corresponding layer, such as in layer 2/3 (channels 3–6) for population L23I (column 2 in Fig. 13*A*). In all cases, the signatures of opposite-sign return currents and also other, weakly connected populations are seen across depth.

As seen in Figure 13*C*,*F* the population-rate predictions are in good qualitative agreement with the ground-truth LFP for both spontaneous and sinusoidally modulated network activity with correlation coefficients (*cc*) between 0.48 and 0.94 for spontaneous activity (Fig. 13*D*) and 0.52 and 0.98 for thalamically evoked oscillations (Fig. 13*G*). Overall, the population-rate predictions appear to be best for the lower frequencies (Fig. 13*E*,*H*), while the inherent variability in the individual point-neuron kernels $HXi(r,\tau )$ (which is not accounted for in the population approximation) has a larger effect on the higher frequencies. This can be understood on biophysical grounds, as the lower frequencies are expected to mainly reflect the gross anatomical features of the postsynaptic populations and their presynaptic connections patterns where the individual variability plays a lesser role (Lindén et al. 2010; Pettersen et al. 2012). The correlation coefficients and power spectra of the population-rate prediction thus show that the population-rate LFP predictor is more accurate than the low-density LFP predictors (Fig. 12*C*,*F*) in case of substantial downscaling.

Although the spike-trigger-averaged LFP (stLFP) (Swadlow et al. 2002; Nauhaus et al. 2009; Jin et al. 2011; Denker et al. 2011), calculated as the cross-covariance between the population spike rate $\nu X(t)$ and the compound LFP $\varphi (r,t)$ divided by the total number of spikes (i.e., $\u222b0T\nu X(t)dt$), is related to our LFP population kernels, it measures very different aspects of cortical dynamics. The population kernels $H\xafX(r,\tau )$ are causal and independent of effects of spike-train correlations. The stLFP, however, is non-causal and strongly depends on spike-train correlations, and thus also network state (Einevoll et al. 2013). The stLFP is thus not only very different from $H\xafX(r,\tau )$, it also varies strongly between the spontaneous and sinusoidally modulated network state (see example for L5E neurons in Fig. 13*B*).

## Discussion

We have here described a hybrid modeling scheme for computing the LFP incorporating both large-scale neural-network dynamics and the biophysics underlying LFP generation on the single-neuron level. The hybrid modeling scheme was illustrated with a full-scale network model of a cortical column in early sensory cortex (Potjans and Diesmann 2014), and the impact of individual populations, network dynamics, and cell density on the mesoscopic LFP signal was investigated.

### The Hybrid LFP Modeling Scheme

The hybrid scheme combines the simplicity and efficiency of point-neuron network models with the biophysics-based modeling of LFP by means of multicompartment model neurons with detailed dendritic morphologies. The neuronal network dynamics are governed by the point-neuron network model independent of LFP predictions. The spikes of the point-neuron network are distributed to the synapses of the multicompartment model neurons with realistic cell-type and layer-specific connectivity. Synapse activation results in spatially distributed transmembrane currents, which are mapped to an LFP signal according to well-established volume-conduction theory.

A main motivation for developing the hybrid LFP modeling scheme was to obtain the ability to compute LFPs for a key class of network models that are amenable to mathematical analysis and can provide intuitive understanding of emerging network dynamics, namely point-neuron models. Similar to networks of anatomically and biophysically detailed neuron models, point-neuron networks can generate realistic spiking activity. In addition, the hybrid modeling scheme brings a substantial computational advantage: with present-day computing and software technologies, point-neuron networks with ∼ 100, 000 neurons can be modeled with laptop computers, and networks comprising millions of point neurons can be routinely simulated on high-performance compute facilities (Helias et al. 2012; Kunkel et al. 2014). Until now, the largest simulation of LFPs based on networks of multicompartmental neuron models with reconstructed morphologies, in contrast, comprised 12, 000 neurons and was done on a Blue Gene/P supercomputer with 4096 CPUs (Reimann et al. 2013). The linearity of electromagnetic theory allows for the implementation of the LFP hybrid modeling scheme as an “embarrassingly” parallel operation (Foster 1995). Therefore, the results for the cortical microcircuit application with ∼ 78,000 neurons were obtained with only 256 CPUs, and could even be acquired with much smaller computing architectures.

A full implementation of the hybrid scheme is provided by the freely available

The present hybrid LFP scheme involves several assumptions with respect to 1) the generation of realistic spiking activity, 2) forward modeling of extracellular potentials, and 3) the combined use of point-neuron networks and multicompartment modeling. In the following, we review the main assumptions and discuss potential extensions:

1. *Spike-train generation by point-neuron networks:* Although highly simplified, single-compartment models of individual neurons (point-neuron models) can mimic realistic spiking for a variety of cell types (Izhikevich 2003; Kobayashi et al. 2009; Yamauchi et al. 2011) and can make accurate predictions of single-cell firing responses under in-vivo like conditions (Jolivet et al. 2008; Gerstner and Naud 2009). Moreover, networks of point neurons can reproduce a number of activity features observed in vivo, such as spike-train irregularity (Softky and Koch 1993; van Vreeswijk and Sompolinsky 1996; Amit and Brunel 1997; Shadlen and Newsome 1998), membrane-potential fluctuations (Destexhe and Paré 1999), asynchronous firing (Ecker et al. 2010; Renart et al. 2010; Ostojic 2014), correlations in neural activity (Gentet et al. 2010; Okun and Lampl 2008; Helias et al. 2013), self-sustained activity (Ohbayashi et al. 2003; Kriener et al. 2014), and realistic firing rates across laminar cortical populations (Potjans and Diesmann 2014). Recently Rössert et al. (2016) demonstrated that the multicompartment-neuron network of Markram et al. (2015) could be reduced to a point-neuron network using an automated step-by-step workflow without significantly affecting the spiking statistics in the network. Nevertheless, under circumstances where simple point neurons may fail to produce plausible spike patterns for a particular input, for example, in the presence of highly nonlinear shaping of postsynaptic responses or dendritic computation (reviewed, e.g., by London and Häusser 2005), more sophisticated point neurons (Breuer et al. 2014) may be required. However, the present hybrid scheme makes no assumptions on generators of spiking activity. In principle, network spikes could, as an example, be replaced by statistical models of spike generation (Lindén et al. 2011; Łe¸ski et al. 2013), multicompartment neuron spikes or even experimentally measured spiking activity.

2. *Biophysical forward modeling of LFPs:* The biophysical forward model described by Equation (10), implemented in

*infinite*,

*isotropic*(same in all directions),

*homogeneous*(same in all positions), and

*ohmic*(frequency-independent) extracellular medium represented by a scalar conductivity $\sigma e$. However, one could generalize the forward model in a straightforward manner to account for anisotropy (Nicholson and Freeman 1975; Logothetis et al. 2007; Goto et al. 2010), or jumps in conductivities at tissue interfaces, for example, using the so-called method of images (Pettersen et al. 2006; Gold et al. 2006; Hagen et al. 2015; Ness et al. 2015). For even more complicated geometrical spatial variations of the conductivity, the forward modeling problem can always be solved by means of finite element modeling (Ness et al. 2015). Recent experiments have only found a small frequency dependence of the extracellular conductivity $\sigma e$ at LFP frequencies ($f\u2272500Hz$, Logothetis et al. 2007; Wagner et al. 2014), but see Gabriel et al. (1996; 2009), Bédard and Destexhe (2009). In any case the forward model could still be applied with frequency-dependent conductivity by means of Fourier decomposition where each frequency component of the LFP signal is considered separately. For more information on possible generalizations of the biophysical forward-modeling scheme, see Pettersen et al. (2012). Finally, we assumed the so-called disc-electrode approximation and averaged the computed LFP signal across the electrode surface (Moulin et al. 2008; Lindén et al. 2014; Ness et al. 2015). Although electrode impedance will affect the measurement, it appears that confounding effects from this can easily be avoided with present-day LFP recording techniques (Nelson and Pouget 2010), hence few compelling reasons exist to incorporate additional temporal filters.

3. *Combined use of point-neuron and multicompartmental models:* The key approximation in the hybrid LFP scheme comes from the combined use of point-neuron (single-compartment) and multicompartmental neuron models. The multicompartment neurons are mutually unconnected, have no outgoing (efferent) connections, and are solely used to compute the LFP. Furthermore, due to dendritic filtering, the somatic postsynaptic potentials in the multicompartment model neurons are not identical to those of their point-neuron counterparts. This inconsistency could, at least partially, be resolved by adjusting the amplitudes and temporal shapes of the synaptic currents in either the multicompartment neurons or the point neurons (Koch and Poggio 1985; Wybo et al. 2013; 2015).

### Applications of the Hybrid LFP Scheme

The hybrid scheme is not limited to the example point-neuron network model and the particular multicompartment neuron models chosen here. It can be applied to networks 1) of arbitrary topology (graph structure, distance dependencies, and dimensions), 2) with any number of populations, 3) with arbitrarily complex point-neuron (e.g., LIF, Izhikevich, MAT, and Hodgkin-Huxley) and synapse dynamics (e.g., current-based, conductance-based, static, and plastic), and 4) any level of biophysical detail in multicompartment neuron models (e.g., morphologies and active channels).

For illustration, we used the hybrid scheme to compute LFPs along a virtual laminar multi-electrode from activity in a multilayered spiking point-neuron network, modeling signal processing in a patch of primary visual cortex. The network consisted of ∼ 78, 000 neurons organized in 4 layers, each with an excitatory and an inhibitory population, representing a cortical patch of ∼ 1 mm^{2} (Potjans and Diesmann 2014). Altogether 16 different cell types and 10 different morphologically reconstructed neurons were used in the LFP calculation. This cortical microcircuit model is well suited for the illustration of the hybrid scheme due to its 1) minimum level of detail in single-neuron dynamics of both point neurons (LIF) and multicompartment neurons (passive membranes), 2) realistic neuron density allowing investigation of the effects of correlations and scaling of network size, and 3) its spatial organization of multiple populations across cortical layers which yields cancellation effects not captured by LFP proxies such as in Mazzoni et al. (2015).

Even though the example application was based on a generic network model biased toward cat visual cortex and not tuned to address specific experiments, its spiking activity nevertheless matched experimental findings on distributions of firing rates, asynchronous irregular activity and propagation of evoked activity (Potjans and Diesmann 2014). We even observed the predicted LFP to be in qualitative accordance with LFP measurements in primary sensory cortices from a variety of animal species and sensory modalities in terms of 1) LFP amplitude, for both spontaneous ($\u22430.1\u20131mV$, see Maier et al. 2010; Hagen et al. 2015; Reyes-Puerta et al. 2016), electrically stimulated activity ($\u223c1\u20133mV$, see Mitzdorf 1985; Castro-Alamancos and Connors 1996) and 2) stimulus-evoked spatiotemporal LFP and CSD patterns (Di et al. 1990; Einevoll et al. 2007; Reyes-Puerta et al. 2015). The fact that our spontaneous and evoked LFP amplitudes are comparable to the observed amplitudes supports the overall biological plausibility of the hybrid LFP scheme.

For the present example, the LFP was dominated by synaptic inputs, and their associated return currents, on excitatory neurons, in particular onto pyramidal cells in layers 2/3 and 6. Furthermore, contributions from inhibitory synaptic inputs typically dominated the contributions from the excitatory inputs, particularly for LFPs stemming from spontaneous network activity. Although the main point of employing the present example was to illustrate the use of the hybrid LFP scheme and not to make predictions for specific neural systems, we note in passing that a dominance of inhibitory synaptic inputs appears to be in agreement with LFPs generated in the CA3 region of hippocampus as observed in an in vitro setting (Bazelot et al. 2010) and in cortex of both humans and macaque in vivo (Telenczuk et al. 2016). In accordance with a previous study (Lindén et al. 2011), we found that correlations in synaptic input play a major role in determining the CSD and LFP stemming from the network activity, for both spontaneous and stimulus-evoked activity. We further showed that due to inevitable correlations between synaptic-input currents the main features of the LFP can only be correctly predicted by a full-scale model. As our ambition is to compute LFPs also for extended point-neuron networks with millions of neurons covering, for example, entire cortical areas, we finally demonstrated how the hybrid modeling scheme can predict LFPs from population firing rates rather than from spikes of individual neurons (Einevoll et al. 2007).

The present microcircuit model application involving simplified, passive multicompartment populations was here used to study the effect of the spatial connectivity on the laminar pattern of spontaneous and stimulus-evoked CSD and LFP signals. The application thus represented a minimal approach incorporating spatial features in LFP predictions of multilaminar point-neuron networks. However, several of the simplifying model assumptions made in the present example application can straightforwardly be generalized. In particular, such generalizations concern 1) the synaptic connectivity between point neurons and the equivalent multicompartment neurons, 2) the absence of active conductances, 3) the positioning of the cells, 4) the reconstructed morphologies, 5) the representation of external inputs, and 6) the fact that the model only encompassed the local circuitry of a ∼ 1 mm^{2} patch of cortex.

(1) *Synaptic connectivity:* While the population-specific connection probabilities, delay distributions, synapse time constant and mean synaptic weights were identical for the connections in the point-neuron network and those between point neurons and LFP-generating multicompartment neurons, the exact realizations of the two types of connectivities were different. In contrast to the point-neuron network, each cell of a particular type had a fixed in-degree, that is, a fixed number of synaptic inputs, a fixed synaptic current amplitude and a layer-specific synapse positioning that depended on compartment surface area, in the LFP modeling step. The use of the hybrid LFP scheme, however, is not restricted to these or any other specific assumptions about the synaptic connectivity patterns. One could, for example, gather all point-neuron network connections and corresponding weights and delays for use with the LFP-generating multicompartment neuron populations. However, for large networks this would require additional computing and memory resources as the number of recorded connection weights and delays grow proportionally to *N*^{2} in a network of *N* neurons with fixed connection probabilities. The present layer-specific synapse positioning could also be specialized further, for example, to account for dendrite-specific connections that may be relevant for LFP predictions in particular if nonlinear synapse and dendrite models are used.

(2) *Active conductances:* In the present study, we neither included the active channels underlying spike generation, nor active dendritic conductances in the multicompartment neuron models (Remme and Rinzel 2011). Experiments suggest that the contribution to the LFP from the former is small in stimulus-evoked recordings from sensory cortex, at least for the low frequencies of the LFP (Pettersen et al. 2008), but see Ray and Maunsell (2011). In any case the contribution of the spikes to the extracellular potentials, including the LFP, could be included in the present scheme by giving each spike produced in the point-neuron network simulations a cell-type-specific spatiotemporal signature in the computation of the extracellular potential (e.g., as calculated in Holt and Koch 1999; Hagen et al. 2015). Substantial effects of active dendritic conductances on the LFP were observed in a two-layered model by Reimann et al. (2013), but this should be further explored. Recent modeling results (Ness et al. 2016) suggest that for the purposes of LFP prediction, active dendritic conductances can, at least for subthreshold potentials, be effectively described by means of “quasi-active linearized” theory (Sabah and Leibovic 1969; Koch 1984; Remme and Rinzel 2011; Ness et al. 2016). This simplifies the LFP modeling substantially and makes the computation similar in complexity to the present case with purely passive membranes.

(3) *Soma positioning:* For model conciseness, the somas of all neurons belonging to the same specific cortical network populations were set to have the same range across cortical depth, see Figure 6*E*. By instead assuming a biologically more plausible distribution of soma depths, the CSD profiles are expected to be spatially smoothed compared with the profiles seen for the present distributions, for example, in Figure 6*F*. Also the LFP profiles will be affected, but to a lesser degree since the LFP profiles already are spatially smoothed due to volume conduction effects.

(4) *Reconstructed morphologies:* We further chose to rely on a small number of highly detailed reconstructed dendritic morphologies from experimental preparations, and partly reused morphologies across neural populations. Obviously, larger sets of distinct reconstructed dendritic morphologies can be used in future studies as they become available. The effect of dendritic morphologies on generated LFP can be further assessed by the use of stylized (Tomsett et al. 2014; Głąbska et al. 2014) or artificially grown morphologies resembling real neurons (Cuntz et al. 2010; 2011; Torben-Nielsen and De Schutter 2014; Mazzoni et al. 2015). While our morphologies were used as is, effects on the predicted LFP from varying morphologies could be quantified using our present measures of LFP power spectra and variance across depth. Such measures could also be extended to include effects of lateral displacement of the recording electrode (Lindén et al. 2011; Einevoll et al. 2013; Łe¸ski et al. 2013; Tomsett et al. 2014), but also cross-channel coherences and stLFPs should be considered as well.

(5) *External input:* In our point-neuron network modified from Potjans and Diesmann (2014), we assumed the depolarizing input from surrounding cortex and remote areas to be represented by deterministic, fixed-amplitude input currents (DC inputs) rather than independent Poisson spike trains. We thus avoided the generation and storage of $O(105)$ high-rate uncorrelated Poisson spike trains and distributing the corresponding spike events onto the multicompartment model neurons, but this can be introduced in future applications.

(6) *Scale:* Our network model represents only an isolated cortical column under $\u223c1$ mm^{2} of pial surface (Potjans and Diesmann 2014), but work is underway to extend the network to a larger scale, for example, by incorporating additional cortical areas (Schmidt et al. 2016) and extending the network size in the lateral directions (Senk et al. 2015). In addition to allowing for predictions of LFPs at several lateral positions as measured by multishank electrodes, the anticipated outcome is also an altered network dynamics and consequently altered LFPs. The present columnar model lacks structured input from other parts of cortex which corresponds to ∼ 50% of all excitatory synapses (see Potjans and Diesmann (2014) and references therein). For example, accounting for different cortical areas and their interactions would allow for the detailed investigation of pathway-specific LFP contributions as recently reviewed in hippocampus by Herreras et al. (2015).

### Effects of Correlations and Network Size

Synaptic inputs to neurons are typically correlated, and as shown in this article and in earlier studies (Lindén et al. 2011; Łe¸ski et al. 2013), these correlations have a major impact on the properties of the generated LFP (and corresponding CSD). There are different contributions to these input correlations, shared presynaptic neurons and correlations in the presynaptic spiking activity, see, for example, Renart et al. (2010), Tetzlaff et al. (2012), and correct LFP predictions require that both effects are properly taken into account. The generation of presynaptic spike trains with realistic correlation structure requires networks of realistic size (van Albada et al. 2015). The contribution from shared presynaptic input can only be accounted for by using realistic statistics of inputs to the LFP-generating multicompartment model neurons, that is, realistic synaptic connection probabilities resulting in realistic statistics of shared inputs. If these two requirements are fulfilled, the properties of single-cell LFPs as well as the correlations between pairs of single-cell LFPs, will be correctly accounted for.

Given synaptic inputs with realistic statistics from sufficiently large point-neuron networks, do we actually need to represent the full population of LFP generating neurons in order to predict a realistic compound LFP? Or can we alternatively get a good estimate of the compound population-LFP signal from a downscaled population of neuronal LFP generators if we know the correct single cell and pairwise LFP statistics, thereby reducing computational costs? As shown in this article, the answer is negative: in the presence of (even tiny) synaptic-input correlations, a realistic compound LFP can only be generated by multicompartmental neuron populations with realistic size and cell densities.

The hybrid modeling scheme allows one to account for both, that is, realistic sizes of networks to generate spike trains with correct correlation structure and realistic sizes and cell densities of the LFP-generating multicompartmental neuron populations. This can be achieved since point-neuron networks can be simulated very efficiently, and the multicompartmental neurons are independent and can be simulated serially (or in an embarrassingly parallel manner).

### Outlook

While we here have focused on the computation of the LFP based on output from spiking point-neuron networks, a similar hybrid approach could be used when the network dynamics is rather modeled in terms of firing rates or even neural fields (Deco et al. 2008). For our example case of a 4-layered cortical network with an excitatory and an inhibitory population in each layer, the present scheme could be adapted directly by replacing the set of spike trains for each population with the corresponding population firing rates (Schuecker et al. 2015) in the LFP-generating step. However, the feasibility and prediction accuracy of such a scheme would have to be investigated in detail.

Another natural development would be to consider other measurement modalities. The present LFP scheme already incorporates the prediction of ECoG (electrocorticography) signals, that is, the electrical signals recorded at the cortical surface, although the LFP forward-modeling scheme may have to be adjusted to account for the discontinuity in electrical conductivity at the cortical surface (Nunez and Srinivasan 2006). One straightforward method to account for discontinuities in our predictions is the so-called method of images, as previously used with our present forward modeling scheme (Hagen et al. 2015; Ness et al. 2015). An extension to EEG (electroencephalography) and MEG (magnetoencephalography) would in principle also be straightforward as the key variable linking single-neuron activity and the measured signal is the single-neuron current dipole moment (Hämäläinen et al. 1993; Nunez and Srinivasan 2006). This dipole moment can be computed from multicompartmental neuron models when the transmembrane or axial currents are known (Lindén et al. 2010; Ahlfors and Wreh 2015). Given the magnitude and orientation of the current dipole moments for all contributing neurons, the EEG and MEG signal can be computed by a linear superposition of single-cell contributions given an appropriate extracellular volume conductor model for the EEG signal. EEG dipole-source localization techniques and EEG forward models have sometimes assumed simplified spherical head models composed of *N* concentric shells, each with different conductivities (Nunez and Srinivasan 2006). Such models could be used also in the present hybrid scheme as well as realistic head models (Kaiboriboon et al. 2012). Note that MEG signal prediction from the current dipole moment is simplified due to its independence on conductivity, and that the variation in magnetic permeability is insignificant across soft tissue, bone, and air (Hämäläinen et al. 1993). Another measurement that could be modeled is voltage-sensitive dye imaging (VSDi) where the signal largely reflects average membrane potentials of dendrites close to the cortical surface. The spatial profile of the weights in the averaging procedure of the VSDi forward model will be determined by the spatial distribution of dye and the propagation of light in the neural tissue (Chemla and Chavane 2010a, 2010b; Tian et al. 2011).

Even though the LFP has been measured for more than half a century, the interpretation of the recorded data has so far largely been qualitative (Einevoll et al. 2013). The biophysical origin of the signal on the single-cell level appears well understood (Rall and Shepherd 1968; Holt and Koch 1999; Gold et al. 2006; Linden 2010; Buzsáki et al. 2012; Pettersen 2012; Einevoll et al. 2013), and several modeling studies have explored the link between neuron and network activity (Pettersen et al. 2008; Lindén et al. 2010; 2011; Łe¸ski et al. 2013; Reimann et al. 2013; Tomsett et al. 2014; Głąbska et al. 2014). However, the computation of LFPs from network activity has until now been too cumbersome and computer intensive to allow for practical exploration of the links between different types of network dynamics and the resulting LFP. Thus, a validation of network models against measured LFP data has essentially been absent. With the present hybrid LFP scheme, accompanied by the release of the simulation tool

## Funding

European Union Seventh Framework Programme (FP7/2007–2013); (604102) Human Brain Project, HBP; (269912) BrainScaleS, the Helmholtz Association through the Helmholtz Portfolio Theme “Supercomputing and Modeling for the Human Brain” (SMHB), Juelich Aachen Research Alliance (JARA), the German Research Foundation (DFG; grant DI 1721/3-1 [KFO219-TP9]), the Danish Council for Independent Research and FP7 Marie Curie Actions COFUND (grant id: DFF-1330–00226); the Research Council of Norway (NFR, through ISP, NOTUR -NN4661 K).

## Notes

We thank Zoltán F. Kisvárday and Armen Stepanyants for useful discussions in the early phase of this study and providing us with experimentally reconstructed morphologies of cat visual cortex neurons. *Conflict of Interest*: None declared.