Data-driven reduction of dendritic morphologies with preserved dendro-somatic responses

Abstract
Introduction
Results
Discussion
Materials and methods
Data availability
References
Article and author information
Metrics

Abstract

Dendrites shape information flow in neurons. Yet, there is little consensus on the level of spatial complexity at which they operate. Through carefully chosen parameter fits, solvable in the least-squares sense, we obtain accurate reduced compartmental models at any level of complexity. We show that (back-propagating) action potentials, Ca²⁺ spikes, and N-methyl-D-aspartate spikes can all be reproduced with few compartments. We also investigate whether afferent spatial connectivity motifs admit simplification by ablating targeted branches and grouping affected synapses onto the next proximal dendrite. We find that voltage in the remaining branches is reproduced if temporal conductance fluctuations stay below a limit that depends on the average difference in input resistance between the ablated branches and the next proximal dendrite. Furthermore, our methodology fits reduced models directly from experimental data, without requiring morphological reconstructions. We provide software that automatizes the simplification, eliminating a common hurdle toward including dendritic computations in network models.

Introduction

Morphological neuron models have been instrumental in neuroscience (Segev and London, 2000). Major experimental discoveries, for instance that N-methyl-D-aspartate (NMDA) channels (MacDonald and Wojtowicz, 1982) produce local dendritic all or none responses (Schiller et al., 2000; Major et al., 2008), or that dendritic Ca²⁺ spikes mediate coincidence detection between distal inputs and somatic action potentials (APs) (Larkum et al., 1999), have been combined in morphological models to arrive at a consistent picture of dendritic integration: the dendrite is an intricate system of semi-independent subunits (Mel, 1993; Poirazi et al., 2003b; Poirazi et al., 2003a), amenable to dynamic regulation (Poleg-Polsky et al., 2018; Wybo et al., 2019), and able to distinguish specific input patterns (Branco and Häusser, 2010; Laudanski et al., 2014).

Nevertheless, morphological models are not without shortcomings. They are highly complex and consist of thousands of coupled compartments, each receiving multiple non-linear currents. The parameters of the models, typically fitted with evolutionary algorithms to electro-physiological recordings (Hay et al., 2011; Almog and Korngreen, 2014; Van Geit et al., 2016), number in the tens of thousands. Since recordings can only be obtained at a few dendritic sites, these fits are under-constrained, and thus susceptible to over-fitting. Accordingly, the single-neuron fitting challenge, where model performance was measured on unseen spike trains, was not won by a biophysical model, but by an abstract spiking model (Gerstner and Naud, 2009; DiLorenzo and Victor, 2013). Finally, many network-level observations can be explained without morphological models (Gerstner et al., 2012).

These shortcomings underline the need to find the essential computational repertoire of a neuron: the set of computations needed to understand brain functions such as learning and memory. Model simplification is crucial in this endeavor, as it elucidates the lowest level of complexity at which computational features are preserved. Conceptually, the simplification thus extracts the essential elements required for the computation from the underlying biophysics. Simulation-wise, the reduced model requires few resources and can thus be integrated in large-scale networks. Experimentally, the reduced model likely leads to a well-constrained fit.

Past simplification efforts can be grouped in two categories: approaches that use traditional compartments, but with adapted parameters (Pinsky and Rinzel, 1994; Destexhe, 2001; Tobin et al., 2006; Bahl et al., 2012; Marasco et al., 2013; Amsalem et al., 2020), and approaches that rely on more advanced mathematical techniques (Kellems et al., 2009; Kellems et al., 2010; Wybo et al., 2013; Wybo et al., 2015). Many authors apply ad hoc morphological simplifications and then adjust model parameters to conserve geometrical quantities, such as surface area (Davison et al., 2000; Hendrickson et al., 2011; Marasco et al., 2013), or electrical quantities, such as attenuation (Destexhe, 2001) and transfer impedance (Amsalem et al., 2020). Other authors propose two-compartment models whose parameters are adjusted to reproduce certain response properties (Pinsky and Rinzel, 1994; Naud et al., 2014). Nevertheless, these approaches often lack the flexibility to be adapted to a wide range of dendritic computations. Furthermore, particular afferent spatial connectivity motifs might be lost, so that essential computations are not captured by such a reduction. Hence, in a network model, the computational relevance of these motifs might unknowingly be ignored. Advanced mathematical techniques on the other hand are very flexible, as they can incorporate the response properties of the morphology implicitly. (Wybo et al., 2013 and Wybo et al., 2015 propose a system of convolutions whose kernels capture the passive properties of the compartmental model, whereas Kellems et al., 2009 and Kellems et al., 2010 linearly transform the compartmental model into a low-dimensional basis.) However, these techniques are not supported by standard simulation software, such as NEURON (Carnevale and Hines, 2004), a considerable hurdle toward their integration in the canonical neuroscience toolset.

Here, we introduce a simplification method with the flexibility of advanced mathematical techniques, but using traditional compartments. The approach accommodates any morphological model in public repositories (McDougal et al., 2017) and commonly used ion channels (Podlaski et al., 2017). The reduced models represent the optimal approximation to the dendritic resistance matrix, evaluated at any spatial resolution of choice. We fit ion channels by approximating the quasi-active resistance matrices (Koch and Poggio, 1985). All our fits are uniquely solvable linear least-squares problems. The obtained models extrapolate well to non-linear dynamics, and reproduce back-propagating APs (bAPs), Ca²⁺ spikes (Larkum et al., 1999), and NMDA spikes (Schiller et al., 2000; Major et al., 2008) with few compartments. Additionally, we investigate whether a dendritic tree with given afferent spatial connectivity motifs can be simplified by ablating subtrees or branches and grouping synapses in the next proximal compartment. We find that effective weight-rescale factors for synapses can be computed if temporal conductance fluctuations stay below a limit that depends on the difference in input resistance between the ablated branch and the next proximal dendrite. Under application of these factors, voltage responses in the simplified tree are preserved. Finally, we demonstrate that our approach can fit reduced models to experimental recordings, without the need to reconstruct full morphologies. We have created a Python toolbox (NEural Analysis Toolbox – NEAT – https://github.com/unibe-cns/NEAT) (copy archived at swh:1:rev:1cb15f36aa0a764105348541d046c85ef38e21ee) that implements this method together with an extension to NEURON to simulate the obtained models.

Results

A systematic simplification of complex neuron morphologies

Our simplification strategy fits compartments with voltage-gated channels to a reduced set of dendritic locations of interest (Figure 1A, left-middle). The locations, together with the morphology, provide a tree structure for the reduced model (Figure 1A, middle-right) where each node corresponds to a traditional compartment and each edge to a coupling conductance. The fit requires that the reduced model responds similarly to perturbative current steps as the full neuron, at the set of chosen locations (Figure 1B). We consider perturbations around any spatially uniform holding potential $v_{h}$ . Experimentally, $v_{h}$ may be reached by injecting a constant current $i_{h}$ , and thus $v_{h}$ is the effective equilibrium potential under this current injection. In models, we do not need to know $i_{h}$ ; we simply assume that $v_{h}$ is the equilibrium potential around which to linearize the model. Note that if $i_{h} = 0$ , $v_{h}$ is the resting membrane potential.

Figure 1 with 1 supplement see all

Download asset Open asset

Flexible and accurate reduction methodology.

(A) For any set of locations on a given morphology (left, here an L2/3 pyramidal cell [Branco and Häusser, 2010]), a reduced compartmental model can be derived (middle), with an associated schematic representation (right). (B) Steps of our approach: (1) choice of locations at which the reduced model should reproduce the full model’s voltage, (2) coupling, leak and channel conductances are fitted to resistance matrices derived from the full model at different holding potentials, and (3) capacitances are fitted to mimic the largest eigenmode of the full model. (C) The resistance matrix of the passive full model (top) restricted to the five locations in A is approximated accurately by the inverse of the conductance matrix of the passive reduced model (bottom). Labels correspond to locations in A. (D) Example components of the quasi-active resistance matrix of the full model, equipped with a Na⁺-channel, as a function of the holding potential $v_{h}$ . Red lines show the four holding potentials at which our methodology evaluates the resistance matrix. Singularities correspond to holding potentials where the linearization is invalid and should be avoided in the fit. Labels correspond to locations in A. (E) Temporal shape of exemplar input impedance kernels of the full model (gray) and their reduced counterparts (blue, dashed). (F) Same as in E, but for transfer impedance kernels.

Suppose we have chosen $N$ locations. Let $δ 𝐢$ be a vector describing perturbative input currents to each of those $N$ locations. The linearized voltage response of the full neuron, at those locations, is given in the Fourier domain for an input of a particular frequency $ω$ by

δ 𝐯 (ω) = Z_{v_{h}} (ω) δ 𝐢 (ω),

with $Z_{v_{h}}$ the quasi-active (Koch and Poggio, 1985) $N \times N$ impedance matrix. In experiments, $Z_{v_{h}}$ is extracted from $δ 𝐯$ measured in response to $δ 𝐢$ . In models, $Z_{v_{h}}$ is algorithmically computed based on the full morphology, the parameters of the passive currents, and the dynamics of the voltage-gated ion channels (Koch and Poggio, 1985). For our reduction method, we consider perturbative current steps, so that after an equilibration period, we obtain a steady-state voltage response $δ 𝐯$ . Note that this corresponds to evaluating (Equation 1) at zero frequency, and that the corresponding matrix $Z_{v_{h}} (ω = 0)$ contains the input and transfer resistances. Throughout the manuscript, we maintain the notation that $z$ without argument signifies $z (ω = 0)$ , representing an individual input or transfer resistance, and that $Z$ signifies $Z (ω = 0)$ , representing the resistance matrix.

Linearizing the reduced compartmental model around a holding potential (see Materials and methods — The conductance matrix) yields an $N \times N$ conductance matrix $G_{v_{h}}$ , with on its diagonal

g_{L i} + \sum_{n \in 𝒩_{i}} g_{C, i n} + \sum_{d \in ℐ_{i}} {\bar{g}}_{d i} l_{d} (v_{h}),

where $g_{L i}$ is the unknown leak conductance of the $i$ -th compartment, $𝒩_{i}$ the set of nearest neighbor compartments, and $g_{C, i n}$ the unknown coupling conductance between compartment $i$ and compartment $n \in 𝒩_{i}$ . The second sum runs over the set $ℐ_{i}$ of all linearized currents of ion channels present in compartment $i$ , where ${\bar{g}}_{d i}$ is the unknown maximal conductance of ion channel $d$ and $l_{d} (v_{h})$ is a factor that follows from linearizing the dynamics of channel $d$ around the holding potential $v_{h}$ (see Materials and methods – Quasi-active channels). Thus, the unknowns $g_{L i}$ , $g_{C, i n}$ for $n \in 𝒩_{i}$ , and ${\bar{g}}_{d i}$ for $d \in ℐ_{i}$ are to be fitted for each compartment, and the factors $l_{d} (v_{h})$ are known and determined by the ion-channel dynamics. The $i j$ ’th off-diagonal element of $G_{v_{h}}$ is $- g_{C, i j}$ – the negative of the coupling conductance between compartments $i$ and $j$ – if $i$ and $j$ are nearest neighbor compartments on the reduced tree structure, and zero otherwise. The conductance matrix relates the reduced model’s voltage response to the perturbative input current steps

G_{v_{h}} δ 𝐯 = δ 𝐢 .

The full neuron and reduced model will behave similarly for all possible perturbative input steps $δ 𝐢$ if their responses $δ 𝐯$ match. From (Equation 1) and (Equation 3), it follows that $Z_{v_{h}}$ should be the inverse of $G_{v_{h}}$ . Consequently, we require that multiplying the known $Z_{v_{h}}$ (measured or calculated) by the parametric $G_{v_{h}}$ yields the identity

Z_{v_{h}} G_{v_{h}} \approx I .

From (Equation 2), it can be seen that (Equation 4) is linear in the parameters that have to be fitted (leak, coupling, and maximal ion-channel conductances of the reduced compartments). By consequence, (Equation 4) can be cast into a least-squares problem and solved accurately (Figure 1C, Figure 1—figure supplement 1A).

Since the linearized ion channel activation $l_{d} (v_{h})$ depends on the holding potential, $Z_{v_{h}}$ changes with $v_{h}$ (Figure 1D). The fit must disentangle the changes in $Z_{v_{h}}$ induced by the various channels. In models, we first block all voltage-gated ion channels in the full and reduced models and fit $G_{pas}$ to $Z_{pas}$ according to (Equation 4), and thus obtain leak and coupling conductances for each compartment. Then, we unblock one ion channel at a time and decompose the conductance matrix into $G_{v_{h}} = G_{pas} + G_{v_{h}, chan}$ , with $G_{v_{h}, chan}$ a diagonal matrix containing the conductances of the unblocked channel, linearized around $v_{h}$ , at each compartment. Thus $G_{v_{h}, chan}$ depends linearly on the unknown maximal ion channel conductance parameters. With $Z_{v_{h}}$ and $G_{pas}$ known, we optimize these maximal conductances, so that the left-hand side of

Z_{v_{h}} G_{v_{h}, chan} \approx I - Z_{v_{h}} G_{pas} .

matches the right-hand side, and that it does so for multiple holding potentials (we chose $v_{h} = - 75, - 55, - 35$ , and $- 15$ mV, Figure 1D, Figure 1—figure supplement 1D). Thus, we obtain an overdetermined system of equations in the unknown maximal channel conductances at each compartment and compute its solution in the least mean squares sense.

Having fixed the steady-state behavior of the reduced model, we are left with the capacitance of each compartment of the reduced model to match the temporal dynamics of the full model as closely as possible. In our reductions, the capacitance of each compartment is thus a parameter to be fitted. To do so, we found it sufficient to consider passive membrane dynamics. Blocking all voltage-gated ion channels in the reduced model ( $G_{v_{h}, chan} = 0$ ), temporal voltage responses follow

diag (𝐜) \dot{𝐯} = - G_{pas} 𝐯 + 𝐢,

with $𝐜$ a vector containing the unknown capacitance of each compartment. We require that voltage-decay back to rest, as described by (Equation 6), matches voltage decay in the full neuron with all active channels blocked, at the $N$ chosen locations. According to linear dynamical systems theory (Strogatz, 2000), this decay can be decomposed as a sum of exponentially decaying eigenmodes, each with an associated eigenvalue $α_{k}$ and eigenvector $ϕ_{k}$ . The former gives the exponential time scale $τ_{k} = 1 / α_{k}$ of the decay and the latter the spatial profile of the mode. The most prominent of these modes has the largest time scale $τ_{0} = 1 / α_{0}$ , and primarily models the voltage decay back to rest through trans-membrane currents (Holmes et al., 1992). In experiments, we extract this mode by fitting an exponential to the voltage decay. In full models, we compute the eigenvalue $α_{0}$ and eigenvector $ϕ_{0}$ – restricted to the $N$ chosen locations – with the separation-of-variables method (Major et al., 1993). In the reduced model, eigenmodes are found as the eigenvalues and eigenvectors of the matrix $S = - diag {(𝐜)}^{- 1} G_{pas}$ . To fit the capacitances, we require that $S$ has $α_{0}$ as eigenvalue and $ϕ_{0}$ as corresponding eigenvector:

S ϕ_{0} = α_{0} ϕ_{0} .

This system of equations is linear in the reciprocals of the capacitances entering in $S$ and can thus be solved efficiently.

To accurately reproduce the spatio-temporal voltage responses of the full model, the reduced model must approximate the ‘impedance kernels’ $z (t)$ (Wybo et al., 2015), the inverse Fourier transforms of the frequency-dependent elements of the impedance matrix. Despite only fitting the largest time scale $1 / α_{0}$ , we accurately reproduce kernels at all times (Figure 1E,F, Figure 1—figure supplement 1B).

Unless a large number of closely spaced compartment sites is retained on the morphology, our reductions replace the true axial currents flowing along the dendritic tree, which are shaped by the distribution of ion channels in the dendritic membrane, with a strong abstraction: a single coupling conductance between compartments that may be far apart. Hence, it is not guaranteed that the resting voltage of reduced model will be equal to the resting voltage of the full model, at the $N$ chosen locations. We ensure that this will be the case by fitting the leak reversals for each compartment. The product of leak conductance and reversal represents a constant current for each compartment, and thus does not influence model dynamics beyond fixing it’s resting voltage. To fit the leak reversals, we evaluate all ion-channel currents in the reduced model at the full model’s resting voltage and require that temporal voltage derivatives are zero. The obtained equations are linear in the reduced model’s leak reversals, and hence can be solved efficiently. Note that our methodology does not guarantee that the obtained reversals are within physiological limits. Thus, they should purely be seen as fit parameters that determine the resting membrane potential.

Reduced models match the voltage response of their full counterparts

We demonstrated the reduction on two computations that require accurate spatio-temporal interaction. We reproduced sequence discrimination in an L2/3 pyramidal cell model (Branco and Häusser, 2010), where a neuron responds more strongly to a centripetal sequence of inputs than to a centrifugal one, by only retaining compartments on a single branch (Figure 2A). We also reproduced input-order detection (Torben-Nielsen and Stiefel, 2009), where a neuron responds more strongly when one input arrives before the other, by only retaining compartments at the soma and the two input sites (Figure 2B).

Figure 2

Download asset Open asset

Voltage-match between full and reduced models for spatio-temporal dendritic computations.

(A) Reduction of full model (right) to a single branch (middle) reproduces sequence discrimination (right), full model in gray, and reduced model colored for different time-steps between inputs, centripetal (dashed), and centrifugal (dash-dotted). (B) Full model (left) and three-compartment reduction (middle, bottom) discriminate temporal order of inputs, where the response to inputs (middle, top) ordered $1 \to 2$ is stronger than $2 \to 1$ . Voltage responses (right) in full model (gray) and reduced model (colored). (**C, D**) Reductions of resp. L5 pyramidal and Purkinje cells with active ion channels at the soma, and excitatory (AMPA+NMDA) and inhibitory (GABA) synapses at dendritic compartment sites. From left to right: full model with compartment sites (soma S and a selected dendrite site D are labeled), reduced model, somatic voltage with zoom on a single AP (full model in gray and reduced model in purple, input spikes at the bottom), dendritic voltage at D (full model in gray and reduced model in green, input spikes at the bottom), and the relative root mean square voltage errors at each compartment (root mean squared error normalized by the standard deviation $R R M S E (v) = \sqrt{Avg {(v_{full} - v_{reduced})}^{2}} / σ_{v_{full}}$ ). Spike coincidence factor $Γ$ is also shown.

We then tested our approach with non-linear currents concentrated at discrete points (hot-spots) along the morphology. In two biophysical models – an L5 neocortical pyramid (Hay et al., 2011; Figure 2C) and a Purkinje cell (Chen et al., 2013; Miyasho et al., 2001; Figure 2D) – we removed active dendritic channels, but added their total conductance at rest to the passive leak (we term this the 'passified’ model), while retaining all active channels at the soma. The soma, with its AP-generating channels, naturally forms such a non-linear hot-spot. We implemented dendritic hot-spots by clustering sets of excitatory (α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid [AMPA] + NMDA) and inhibitory (γ-aminobutyric acid [GABA]) synapses at randomly selected points on the morphology. Both the reduced pyramidal cell model and the reduced Purkinje cell model accurately reproduced the somatic and dendritic membrane voltage traces of their full counterparts (Figure 2C,D, middle panels). Furthermore, 97% and 89% of APs coincided within a 6 ms window, respectively (Jolivet et al., 2008; Figure 2C,D, right panels).

Models with passive dendritic membranes but non-linear synapse clusters can already reproduce much of the canonical computational repertoire of dendrites. Clusters of AMPA+NMDA synapses at the distal tips of basal dendrites can make neurons function as two-layer networks (Mel, 1993; Poirazi et al., 2003b; Poirazi et al., 2003a). Such reduced models also exhibit the same robustness to input noise as the full models (Poleg-Polsky, 2019), and can include modulation of the co-operativity between AMPA+NMDA clusters through compartments with shunting conductances (Wybo et al., 2019). Furthermore, because our reduced models reproduce dendritic input impedance properties, long time-scale NMDA currents activate strongly in distal compartments as high distal input impedance helps to overcome their voltage-dependent Mg²⁺ block (Jahr and Stevens, 1990a; Jahr and Stevens, 1990b), whereas short time-scale AMPA currents dominate in more proximal compartments. Thus passive models with non-linear synapse clusters also capture the shift in computational strategy from a proximal temporal code to a distal rate code (Branco and Häusser, 2011).

What level of spatial complexity reproduces characteristic responses generated by somatic and dendritic ion channels? These responses arise through the interplay of integrative properties of the morphology – captured in the impedance kernels $z (t)$ – and channel dynamics. If a full model (here the L2/3 pyramid) is reduced to a single somatic compartment (Figure 3A), the impedance kernel has only one time scale – the membrane time-scale (Figure 3B, top – purple). The kernel of the full model contains additional shorter time scales (Figure 3B, top – gray), leading to a faster AP onset than in the reduced model (Figure 3B, bottom). This effect is thus a fundamental limitation of the single-compartment reduction. Extending the compartmental model with four nearby dendritic sites (Figure 3A) adds fast components to the somatic impedance kernel (Figure 3B, top – green) that model the spread of charge to the other compartments. These components increase AP amplitude and decrease AP delay (Figure 3C), leading to a better match with the full model.

Figure 3 with 1 supplement see all

Download asset Open asset

Dendritic computations with active channels are captured by our reduced models.

(**A–C**) Effect of compartment distribution on AP dynamics in reduced models. (A) One- and five-compartment reductions of the L2/3 pyramid, equipped with somatic and dendritic ion channels. (B) Differences in short time-scale behavior in somatic input impedance kernels (top) between full model (gray) and one- (purple) and five-compartment (green) reductions result in different AP delays (bottom). (C) AP amplitude (top) and AP delay (bottom) for the three models (colors as in B). (**D–F**) Effect of compartment distribution on AP back-propagation in basal dendrites. (D) Four- and two-compartment reductions of a basal branch. (E) APs at soma (top) and most distal compartment site (bottom) for four models (full in gray, four compartments in blue, two compartments in orange, and full but with a passive dendrite in green). (F) Amplitude (top) and delay (bottom) for bAPs at different distances from soma (if compartment is present in model), averaged over all basal branches longer than 150 μm (error bars indicate standard deviation). (**G–I**) AP back-propagation in the apical dendrite of the L5 pyramid. (G) Reductions of the apical dendrite with increasing inter-compartment spacing. (H) Voltage waveform at soma (top, full in gray, $Δ x = 100$ μm in blue, $Δ x = 200$ μm in orange) and two dendritic sites (middle, bottom). (I) Waveform amplitude as a function of distance to soma for various inter-compartment spacings. (**J, K**) Ca²⁺-spike-mediated coincidence detection. (J) Reduction of the L5-pyramid’s apical dendrite to 11 compartments. (K) Response to a somatic current pulse (top, left), a dendritic synaptic current waveform (top, right), and the coincident arrival of both inputs (bottom).

Dendritic ion channels support the back-propagation of APs (Stuart et al., 1997). For each basal branch of at least 150 μm in the L2/3 pyramid, we derived reduced models with three compartments (at 50, 100, and 150 μm from the soma) and with one compartment. We compared bAP amplitudes and delays relative to the somatic AP (Figure 3E) and found that even models with a single distal compartment support active back-propagation (Figure 3E,F). In the apical dendrite of the L5 pyramid (Figure 3G), we found that a distance step of 100 μm between compartments was required to support bAPs to the same degree as in the full model (Figure 3H,I).

Finally, we considered the pairing of a somatic current pulse with a dendritic post-synaptic potential waveform. We found that an 11-compartment reduction (Figure 3J) – with compartments spaced at 100 μm through the apical trunk to support bAPs until the Ca²⁺ hot-zone – reproduced the Ca²⁺-mediated coincidence detection mechanism (Larkum et al., 1999; Pérez-Garci et al., 2013).

Conditions under which afferent spatial connectivity motifs can be simplified

Our method allows us to reduce morphological complexity by removing any branch or subtree without affecting the integrative properties of other branches or subtrees. Up until now, we have only considered reductions where the to-be removed branch or subtree does not receive direct synaptic input. However, if this branch or subtree does receive synaptic input, the morphological reduction also simplifies the afferent spatial connectivity motifs targeting it, as synapses have to be grouped at the nearest proximal intact compartment. We assess whether the computational repertoire of the neuron remains intact under such a simplification, or which elements of the repertoire might be lost. We study two proxies for this repertoire: the voltage difference between full and reduced models at the intact compartments – in order to assess whether the output of local computations in the ablated branch or subtree can be recovered, and the difference between full and reduced models in NMDA-spike threshold at the nearest proximal intact compartment (quantified as the smallest number of AMPA+NMDA synapses with a weight of 1 nS that need to be activated simultaneously to elicit an NMDA spike) – in order to assess whether local integrative properties at the intact compartments can also be retained under the shift of synapses. We allow ourselves the liberty of rescaling the weights of the shifted synapses by a temporally constant factor $β$ (Figure 4A). Which spatial synapse distributions admit simplification in this sense, and under which input conditions?

Figure 4

Download asset Open asset

Simplification of afferent spatial connectivity motifs.

(A) Removal of a branch with a synapse (red triangle) is considered possible if the correct voltages at the compartment sites (here, dendritic – green square and somatic – purple square) can be obtained by shifting the synapse to the compartment site and rescaling its weight with a fixed factor $β$ for the all input conditions under consideration. (B) Comparison between weight-rescale factors for current-based (left) and conductance-based (right) input. Voltage trace for dendritic compartment without rescaling (gray) and with the current-based scale factor $β_{c u r r}$ compensating for attenuation (blue). Bottom right panel shows voltage trace with conductance-based scale factor $β_{c o n d}$ . (C) Current- and conductance-based scale factors in the green dendritic branch, for a shift of a synapse at a given distance from the soma to the dendritic compartment site (green square). (D) Spatial peak voltage without (left, blue) and with (right, orange) application of $β_{c o n d}$ . (E) For a cluster of AMPA+NMDA synapses in isolation, scale factors for physiological constants can be obtained (see Materials and methods — Synaptic weight-rescale factors) that reproduce the correct voltage waveform. Colors as in D. (F) Maximal amplitude of NMDA-spike waveform upon activation of increasing numbers of synapses – NMDA-spike thresholds indicated with vertical lines. Colors as in D. (**G–H**) Removing a whole subtree and shifting multiple synapses (red triangles) to the next proximal compartment site (green square). NMDA-spike generation (gray voltage trace) at the compartment site through burst inputs to local AMPA+NMDA synapses (green inputs), with AMPA (blue) and GABA (red) background inputs spread throughout the subtree. Reductions shown without rescaling (blue, dash-dotted), with the analytical single-site rescaling rule (red, full) and the numerical multi-site rule (orange, dashed). (I) Error in NMDA-spike threshold for the three cases in H. (J) Dependence of the error in NMDA-spike threshold on the factor $Δ z g$ , with $Δ z$ the average input resistance difference between synapse sites and the compartment and $g$ the average synaptic conductance.

For a single, current-based synapse, shifted from its original location on the ablated branch to the next proximal compartment, we analytically compute the weight-rescale factor and find $β_{curr} = z_{c s} / z_{c c}$ , with $z_{c s}$ the transfer resistance from synapse site $s$ to compartment site $c$ and $z_{c c}$ the input resistance at $c$ . Since $z_{c s} < z_{c c}$ (Koch, 1998), the weight-rescale factor weakens the shifted synapse, so that it matches the attenuation of its distal counterpart. Nevertheless, when $s$ is located more distally on the same branch or subtree as $c$ , the transfer resistance $z_{c s}$ is often close to the input resistance $z_{c c}$ (Figure 1—figure supplement 1E). Thus, $β_{curr}$ is often close to one, and hence negligible (Figure 4B,C). For a conductance-based synapse, rescaling weights by $β_{curr}$ does not accurately fit the full model (Figure 4B). We decompose the synaptic conductance $g (t)$ as $g + δ g (t)$ , with $g$ the temporal average and $δ g (t)$ the fluctuations around $g$ . We then analytically compute (see Materials and methods — Synaptic weight-rescale factors) that the weight-rescale factor for a conductance-based synapse

β_{cond} = \frac{1}{1 + Δ z g}

recovers the full model’s voltage if

Δ z δ g (t) ≪ 1 + Δ z g,

where $Δ z = z_{s s} - z_{c c}$ is the difference between input resistance $z_{s s}$ at the original synapse site $s$ and input resistance $z_{c c}$ at the nearest proximal compartment site $c$ . Multiplying the synaptic weights by $β_{cond}$ weakens the synapse so that it matches the shunt effect of its distal counterpart, and so that the voltage at all compartments more proximal than $c$ and in their respective side branches is reproduced (Figure 4D).

If $s$ and $c$ are close, $Δ z$ is small and (Equation 9) is satisfied for plausible conductance fluctuations. Thus, shifting the synapse to the next proximal compartment keeps the voltage response intact. If the separation between $s$ and $c$ increases, $Δ z$ increases accordingly, and the magnitude of tolerated conductance fluctuations shrinks (Equation 9). Thus, for a given magnitude of $δ g (t)$ , it is possible that a thick branch can be removed, but not a thin branch, as the input resistance increase in the former is much more gradual than in the latter. Nevertheless, (Equation 9) suggests that even for large $Δ z$ the full model’s voltage can be recovered by downscaling synaptic weight by $β_{cond}$ , but only if $δ g (t) ≪ g$ . A tonic level of AMPA and/or GABA activation, as in the high-conductance state (Destexhe et al., 2003), can implement this input condition. For NMDA synapses, this input condition does not hold, as NMDA spikes arise through a strong, transient conductance increase. However, when a cluster of such synapses – to be moved to a proximal site – is considered in isolation, a workaround can be found by applying constant rescale factors not just to synaptic weight, but also to threshold and width of the Mg²⁺ block, as well as to the reversal (see Materials and methods — Synaptic weight-rescale factors). These factors recover NMDA-spike shape (Figure 4E) and threshold (Figure 4F) when the cluster is moved to from $s$ to $c$ .

Our analytical weight-rescale factor $β_{cond}$ treated reductions where a single site $s$ was moved to $c$ . We consider now the ablation of a whole subtree and the shift of all its synapses to $c$ (Figure 4G). In our simulations, we distributed 100 AMPA and 100 GABA synapses on the to be ablated subtree, and activated them with a Poisson process with fixed rate. We studied how the reduction influenced the local integrative properties at the intact compartment $c$ . As a proxy for those integrative properties, we checked the shape of the waveform elicited by activating AMPA+NMDA synapses at $c$ in a short burst (Figure 4H), and we quantified the threshold for NMDA-spike generation. We found that neither shape nor threshold were preserved by the analytical weight-rescale factors $β_{cond}$ , derived for the case where only a single input site is present on the subtree. We numerically extend the derivation of the weight-rescale factors to multi-site inputs (see Materials and methods – Synaptic weight-rescale factors). These extended weight-rescale factors $β_{cond}^{ext}$ reproduced local voltage (Figure 4H) and NMDA-spike threshold (Figure 4I). We also investigate whether $β_{cond}^{ext}$ can still be conceptualized as depending on $Δ z_{avg} g_{avg}$ – the product between (1) a difference in average input resistance between synapses on the to be ablated subtree and the compartment site and (2) the average conductance load exerted by those synapses. To do so, we activated the 100 AMPA and 100 GABA synapses on the to be ablated subtree at a wide range of input conditions (see Materials and methods – Simulation-specific parameters), and measured the error in threshold for NMDA-spike generation at $c$ . We found that this error depended strongly on $Δ z_{avg} g_{avg}$ and remained small for $Δ z_{avg} g_{avg} ≪ 1$ (Figure 4J).

In conclusion, we can now assess whether given spatial synapse motifs on branches or subtrees admit simplification by shifting the synapses to the nearest proximal compartment site and ablating the branch or subtree. We find that reduction in this sense is possible for any motif if fluctuations are small, so that $δ g (t) ≪ g$ , by applying temporally constant weight-rescale factors. Otherwise, the synapses must be located sufficiently closely to the compartment site, so that $Δ z$ is small and that $Δ z δ g$ remains below $1 + Δ z g$ for all relevant activation levels.

Active and passive reduced dendrites under synaptic bombardment

We next investigated how many compartments are needed to reproduce somatic and dendritic responses under a synaptic bombardment and quantified the contribution of active dendritic Na⁺, Ca²⁺, and K⁺ channels. We mimicked the in vivo state in the basal dendrites of the L2/3 pyramid by distributing 200 AMPA and 200 GABA synapses receiving Poisson inputs and 20 clusters of AMPA+NMDA synapses receiving bursts of inputs (Figure 5A). We sampled 10 different sets of meta-parameters governing synaptic activation (Poisson and burst rates and synaptic weights, see Materials and methods — Simulation-specific parameters), resulting in output spike rates between 0 and 8 Hz (Figure 5B). We derived reduced models, once based on the full model with all active channels (‘active reduced dendrite’ in Figure 5C) and once based on the passified full model (‘passive reduced dendrite’ in Figure 5C). We distributed increasing numbers of compartment sites on the basal branches and measured somatic and dendritic voltage responses for the active full model and for the reduced models with active and passive dendrites (Figure 5C, Figure 5—figure supplement 1 for dendritic traces). Spike coincidence (Figure 5D), subthreshold somatic voltage error (Figure 5E), and dendritic voltage error (Figure 5F) showed a significant improvement when compartment numbers increased from 0 (point-neuron – top row in Figure 5C) to 2 per 100 μm. In the presence of active channels, the error decreases until up to three compartments per 100 μm (Figure 5E,F). Spike coincidence factors were consistently ∼10% higher than in the passive case (Figure 5D).

Figure 5 with 1 supplement see all

Download asset Open asset

Reductions with active and passive dendritic compartments under in vivo like conditions.

(A) AMPA and GABA synapses, and AMPA+NMDA synapse clusters are spread randomly throughout the basal dendrites of the L2/3 pyramid (inset). AMPA and GABA synapses receive Poisson inputs and AMPA+NMDA synapses receive bursting input (right). (B) Output spike rates of the full model for 10 different input configurations. (C) Reductions with increasing numbers of compartment sites along the basal dendrites (left). Somatic voltage traces and spike times for the full model (gray), for a reduction with passive dendrites (middle, blue) and a reduction with active dendrites (right, orange). (**D–F**) Spike coincidence factors (D), relative somatic voltage errors (E) and relative dendritic voltage errors (F) for reductions with increasing numbers of compartments, and with passive (blue) and active dendrites (orange), averaged over all 10 input configurations. Relative voltage errors are computed as in Figure 2C,D.

Fitting reduced models directly to experimental data

Traditionally, creating a morphological neuron model involves a reconstruction of the morphology and an evolutionary fit of the electrical parameters (Hay et al., 2011; Almog and Korngreen, 2014; Van Geit et al., 2016). Our method allows skipping these resource-intensive steps (Figure 6A). We adapt our method to a common experimental paradigm where hyper- and depolarizing current steps are injected under applications of various ion-channel blockers (Marti Mengual et al., 2020). We extract resistance matrices and holding potentials from voltage step amplitudes (Figure 6B), and the time scale of the largest eigenmode from the average voltage decay back to rest (Figure 6C). If the modeling goal is a reduced model, the traditional reconstruction and subsequent reduction both introduce errors (Figure 6D). Thus, we find that while our reduction of the full model is faithful within an expected accuracy margin (Figure 6E), fitting a reduced model directly to the experimental traces is more accurate (Figure 6F).

Figure 6

Download asset Open asset

Fitting reduced models directly to experimental data.

(A) Comparison between the traditional neuron model creation paradigm (morphological reconstruction and evolutionary fit, possibly followed by a reduction), and the proposed direct experiment to reduced dendrite model paradigm. (B) Resistance matrices and holding potentials are extracted from the response amplitude to hyper- and depolarizing current step inputs, here measured once under application of an h-channel antagonist and once under control conditions. (C) Time scale of largest eigenmode is extracted from average decay back to the resting membrane potential. (D) Combined root mean square voltage error $R M S E (v) = \sqrt{Avg {(v_{full} - v_{reduced})}^{2}}$ of full model (gray) – fitted to current step data through an evolutionary algorithm, of the reduction of the full model (red), and of the direct fit of the reduced model to the data (brown). (E) Experimentally recorded traces (gray), traces from the full reconstruction (black), and its reduction (green and purple, dashed). (F) Experimentally recorded traces (gray) and traces from the directly fitted reduced model (green and purple, dashed).

Discussion

We have presented a flexible yet accurate method to construct reduced neuron models from experimental data and morphological reconstructions. The method consists of linear parameter fits at holding potentials that are informative for the full non-linear dynamics. First, we derive leak and coupling conductances from the passified version of the full model. Second, we fit the maximal conductances of the active ion-channels, so that the reduced model, at various holding potentials, responds similarly to input perturbations as the full model. Third, we fit the capacitances to reproduce spatio-temporal integrative properties. Finally, we ensure that resting membrane potentials in the full and reduced models match. The resulting method can be adapted to a wide variety of use cases.

A fundamental use case is deriving reduced models that retain specific elements of the dendritic computational repertoire, to explore how those dendritic computations may improve network function. To that purpose, one can apply our method to compartments placed at input sites needed for the computation of interest. For instance, if a computation requires independent computational subunits (Poirazi et al., 2003b; Kastellakis et al., 2016), one would place compartments at the dendritic tips. If one wants to modulate co-operativity of excitatory inputs at distal sites, one would place compartments with shunting conductances at the bifurcations between the distal sites (Wybo et al., 2019). Detecting input sequences (Branco and Häusser, 2010) or implementing a dendritic integration gradient (Branco and Häusser, 2011) would require placing two or more compartments along the branches of interest. For many of these computations, a passive morphology with non-linear (NMDA) synapses suffices. Furthermore, our method can be combined with abstract models of AP generation (Pozzorini et al., 2013; Pozzorini et al., 2015).

A second use case is constructing models directly from electro-physiological recordings. Our method avoids the labor-intensive detour of reconstructing the morphology and optimizing model parameters with an evolutionary algorithm (Marti Mengual et al., 2020). In combination with advances in voltage-sensitive dye imaging, our method may thus form the basis of a high-throughput experiment-to-model paradigm.

A third use case is elucidating the effective complexity of a given dendritic tree. We assessed whether afferent spatial input motifs could be simplified by removing a branch or subtree and grouping the affected synapses at the next proximal compartment. We found that such simplification is possible when the input resistance of the affected synapse sites is close to the input resistance of the compartment site. If this condition does not hold, the simplification can still proceed if conductance fluctuations are small compared to the average conductance magnitude, but not otherwise. What is the simplest model for a given dendritic tree that captures its full computational repertoire? We find that spike coincidence and subthreshold voltage metrics reach satisfying accuracy at, but do not improve much beyond ∼3 compartments per 100 μm. Further research is required to tease apart mere numerical errors from potentially missed computational features.

The process of simplifying complex dynamical systems, such as morphological neuron models, is inherently approximate. Due to the vast variety of morphologies, voltage-gated ion-channel configurations and input patterns, it is beyond the scope of any single paper to give an exhaustive overview of possible inaccuracies. Rather, we propose three lines of inquiry if a reduction fails to reproduce a particular response pattern. First, we suggest checking the discrepancy between the impedance kernels $z (t)$ of the full and reduced models. Often, response mismatches arise from slight inaccuracies in how ion-channel dynamics are integrated by the passive system (Figure 3A–C). Second, we suggest assessing to what degree afferent spatial input motifs, present in the full model, are simplified, and whether that simplification still admits a sufficient accuracy in reproducing local voltage-dependent responses, such as NMDA spikes (Figure 4). Finally, we suggest checking whether the distribution of compartments on the morphology admits ion channel dynamics that are sufficiently rich to reproduce the response characteristic in question. To illustrate these lines of inquiry, we have further simplified the L5 pyramid to five compartments, by omitting the compartments in the apical trunk (Figure 3—figure supplement 1A vs Figure 3J). This reduction results a change in AP shape, and in a failure to elicit AP bursts through the pairing of a somatic and a dendritic input, even though a Ca²⁺ spike can still be elicited (Figure 3—figure supplement 1B). Likely, the change in AP shape can be traced back to the difference in the somatic input impedance kernels (Figure 3—figure supplement 1C), while the failure to elicit an AP burst may be due to insufficient Na⁺-driven cross-talk between soma and apical tuft.

Computational tools have played a key role in neuroscience. Take NEURON (Carnevale and Hines, 2004), which accelerated our understanding of morphological neurons, or NEST (Gewaltig and Diesmann, 2007), which enabled simulation of large-scale spiking point-neuron networks. We have implemented a Python toolbox that automatizes the simplification process (NEural Analysis Toolkit – NEAT – https://github.com/unibe-cns/NEAT; copy archived at swh:1:rev:1cb15f36aa0a764105348541d046c85ef38e21ee). The toolbox reads any morphology in the standard ‘.swc’ format (Ascoli, 2006) and returns the parameters of the reduced models, while also providing tools to export the models to NEURON (Carnevale and Hines, 2004).

Our method and toolbox fill a void in between the extremes of modeling large-scale networks with abstract models and modeling single cells in all their detail. By enabling the efficient systematic derivation of simplified neurons, amenable to simulation at the network level, our work bridges the gap between two branches of neuroscience that historically have remained separate.

Share this article

Cite this article

Flexible and accurate reduction methodology.

Voltage-match between full and reduced models for spatio-temporal dendritic computations.

Dendritic computations with active channels are captured by our reduced models.

Simplification of afferent spatial connectivity motifs.

Reductions with active and passive dendritic compartments under in vivo like conditions.

Fitting reduced models directly to experimental data.

Author details

Willem AM Wybo

Contribution

For correspondence

Competing interests

Jakob Jordan

Contribution

Competing interests

Benjamin Ellenberger

Contribution

Competing interests

Ulisses Marti Mengual

Contribution

Competing interests

Thomas Nevian

Contribution

Competing interests

Walter Senn

Contribution

Competing interests

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

Categories and tags

Further reading