An efficient analytical reduction of nonlinear detailed neuron models

Detailed conductance-based neuron models, consisting of nonlinear branched dendrites and thousands of synapses, are essential for understanding the integrative and computational properties of single neurons and large neuronal networks, and for interpreting experimental results. Simulations of such models are computationally expensive, severely limiting their utility. We introduce a novel analytic approach to simplify complex nonlinear neuron models while preserving the identity of individual dendrites and synapses. Neuron_Reduce represents each stem dendrite by a unique cylindrical cable, keeping its specific membrane and axial properties. Neuron Reduce maps synapses and active membrane ion channels to the respective cylinder while preserving their transfer impedance to- and from- the soma as in the detailed model. The reduced model accelerates the simulation speed by up to 200-fold while closely replicating the sub- and supra- threshold voltage dynamics for a variety of cell types and inputs, including the nonlinear “ping pong” interaction between somatic Na+- and dendritic Ca2+- spikes, found in L5 neocortical pyramidal cells. Neuron_Reduce also replicates dendritic computations discriminating spatiotemporal input sequences. The reduced neuron models will enable realistic simulations of neural networks at unprecedented scale, including of biologically-inspired “deep networks” and facilitate the construction of neuromorphic-based systems. Neuron_Reduce is publicly available (https://github.com/orena1/neuron_reduce) and is straightforward to implement.


Introduction
Compartmental models (CMs) were first employed by Wilfrid Rall to study the integrative properties of neurons 1 . They enabled him to explore the impact of spatiotemporal activation of conductance-based dendritic synapses on the neuron's output and the effect of the dendritic location of a synapse on the time course of the somatic excitatory post-synaptic potential 2 . By simulating electrically-distributed neuron models, Rall demonstrated how the cable properties of dendrites explain the variety of somatic EPSP shapes that were recorded at the soma of α-motoneurons, thus negating the dominant explanation at that time that the differences in shapes of the somatic EPSPs in these cells result from differences in the kinetics of the respective synapses. This was an impressive example that faithful models of the neuron (as a distributed rather than a "point" electrical unit) are essential for the correct interpretation of experimental results. Since Rall's 1964 and 1967 studies using CMs, EPSP "shape indices" measured at the soma are routinely used for estimating the electrotonic distance of dendritic synapses from the soma.
Over the years, detailed CMs of neurons have provided key insights into hundreds of experimental findings, both at the single cell and the network levels. A notable example at the single cell level is the explanation as to why the somatic Na + action potential propagates backward in the soma-to-dendrites direction and (typically) not vice versa 3 .
CMs have also pinpointed the conditions for the generation of local dendritic Ca 2+ spikes [4][5][6] and provided an explanation for the spatial restriction of the active spread of dendritic spikes from distal dendrites to the soma 7 and see also [8][9][10][11][12][13][14] . Today, detailed CMs models are even been used for simulating signal processing in human pyramidal neurons, including their large numbers of dendritic spines/synapses 15 .
At the network level, detailed CMs are utilized for such noteworthy projects as largescale simulations of densely in silico reconstructed cortical circuits 16,17 and the overarching goal of the Allen Institute to simulate large parts of the visual system of the mouse 18,19 . Because detailed compartmental modeling is increasingly becoming an essential tool for understanding diverse neuronal phenomena, major efforts have been invested in developing user-friendly computer software that implements detailed CMs, the best known of which are NEURON 20 , GENESIS 21  Modern personal computers can simulate tens of seconds of electrical activity of single neurons comprising thousands of nonlinear compartments and synapses. However, they handle poorly cases where many model configurations need to be evaluated such as in large-scale parameter fitting for single-neuron models 5,28 , or when the dendritic tree is morphologically and electrically highly intricate and consists of tens of thousands of dendritic synapses, as with cortical pyramidal neurons 15 . When the aim is to simulate a neuronal network consisting of hundreds of thousands of such neurons, only very powerful computers can cope. For example, the simulation of a cortical network consisting of 200,000 detailed neuron models on the BlueGene/Q supercomputer takes several hours to simulate 30 seconds of biological time 17 .
To overcome this obstacle, two approaches have been pursued. The first involves developing alternative, cheaper and more efficient computing architectures (e.g., neuromorphic-based computers 29,30 ). These have not yet reached the stage where they can simulate large-scale network models with neurons consisting of branched nonlinear dendrites having a realistic number of synapses. The other approach is to simplify neuron models while preserving their input/output relationship as faithfully as possible.
Rall 31 was the first to suggest a reduction scheme in his "equivalent cylinder" model, which showed that, for certain idealized passive dendritic trees, the whole tree could be collapsed into a single cylinder that was analytically identical to the detailed tree. The "equivalent cylinder" preserves the total dendritic membrane area, the electrotonic length of the dendrites, and most importantly, the post-synaptic potential (amplitude and time-course) at the soma for a dendritic synapse when mapped to its respective electrotonic location on the "equivalent cylinder" 32,33 . However, this method is not applicable for dendritic trees with large variability in their cable lengths (e.g., pyramidal neurons with a long apical tree and short basal trees), conductance-based synapses, or for dendrites with nonlinear membrane properties.
Over the years, several different reduction schemes have been proposed; e.g., a recent work mapped all the synapses to a single compartment, taking the filtering effect of the dendrites into account 34 . Other methods reduce the detailed morphology to a simplified geometric model while preserving the total membrane area [35][36][37] or the axial resistivity 38 ; see also 12,39,40 . However, these methods have a variety of drawbacks; in particular, they are either "hand fitted" and thus lack a clear analytical underpinning or are complicated to implement, and in some cases, their computational advantage for realistic numbers (thousands) of synapses is not quantified. Most of these methods do not support dendrites with active conductances 35,38,39,41,42 and they have not been tested on a broad range of neuron types. Importantly, none of the previous methods provided an easy-touse open access implementation. Thus, today there is no simple, publicly available reduction method for neuron models that can be used by the extensive neuroscience and machine-learning communities.
To respond to this need, the present study provides an analytic method for reducing the complexity of detailed neuron models while faithfully preserving the essential input/output properties of these models. Neuron_Reduce is based on key theoretical insights from Rall's cable theory, and its implementation for any neuron type is straightforward without requiring hand tuning. Depending on the neuron modeled and the number of synapses, Neuron_Reduce accelerates the simulation run-time by a factor of up to 250 while preserving the identity of individual synapses and their respective dendrites. It also preserves specific membrane properties and dendritic nonlinearities, hence maintaining specific dendritic computations. Neuron_Reduce is easy to use, fully documented, and publicly available on GitHub (https://github.com/orena1/neuron_reduce).

Analytic mapping of a detailed dendritic tree to a reduced multi-cylinder tree
The thrust of our analytical reduction method (Neuron_Reduce) is described in Fig. 1a-c. This method is based on representing each of the original stem dendrites by a single cylindrical cable which has the same specific membrane resistivity (Rm, in Ωcm 2 ), capacitance (Cm, in F/cm 2 ) and axial resistivity (Ra, in Ωcm) as in the detailed tree (Fig.   1a). Also, each cylindrical cable satisfies two constraints: (i) the magnitude of the transfer impedance, "# $,& (()" = |, $ (()/. & (()|, from its distal sealed end (X = L) to its origin at the soma-end (X = 0) is identical to the magnitude of the transfer impedance from the electrotonically most distal dendritic tip to the soma in the respective original dendrite; (ii) at its proximal end (X = 0), the magnitude of the input impedance, , $ ((), to an input current, . / ((), injected at a dendritic location, x, will be identical in both the detailed and the reduced cylinder models (see Methods). Consequently, synapses and nonlinear ion channels are mapped to their respective loci in the reduced cylinder while preserving the respective transfer resistance to the soma (see Fig. 1 Step B and Methods). Based on Eqs. (1) - (11), Neuron_Reduce generates a reduced multicylindrical tree for any ω value (different reduced models for different ω values).
Conveniently, we found a close match between the detailed and the reduced models for ω = 0 (the steady-state case); Therefore, all figures in this work are based on reduced models with ω = 0 (see Discussion) Figure 1. An analytic method for reducing neuron model complexity (Neuron_Reduce). a Detailed passive model of 3D reconstructed L5 thick tufted pyramidal cell from rat neocortex. Its nine stem dendrites (one apical and 8 basal) are depicted in different colors. b Each original stem dendrite is reduced to a single cylinder that retains the specific passive cable properties (Rm, Cm, and Ra) of the original tree. The diameter and length of the respective cylinders are computed analytically using Eqs.
(1) - (11), such that each cylinder preserves both the transfer resistance from the most electrotonically distal dendritic tip to the soma as well as the input resistance at the soma end of the corresponding stem dendrite. This generates a unique cylindrical cable for each of the original stem dendrites. Scale bars in a and b are 100 µm. c Synapses with similar transfer resistance to the soma (exemplar synapses are marked as 1-4 at top right) are all mapped to the respective locus in the reduced cylinder so that their transfer resistance is similar in the two models. In the reduced model, these synapses are merged into one "NEURON" process (red synapse in b), but they retain their individual activation time (see Methods and Supplementary Fig. 1). The same mapping also holds for active membrane conductances (schematic yellow region, denoting the Ca 2+ "hot spot" in the apical tree). d Transfer impedance (# 0,$ = # $,0 ) between point d on the apical tree (shown in a and b) and the soma (X = 0) as a function of the input frequency in both the detailed (black trace) and the reduced (red trace) models. e Composite somatic EPSPs resulting from sequential activation of the four distal apical synapses shown in c in the detailed model (black trace) and in the reduced model (red trace). In this simulation the dendritic tree was passive. The synapses were activated in temporal order 1, 2, 3, 4 as shown by the vertical lines below the composite EPSP. The respective peak conductances of these AMPA-based synapses were 0.6, 0.3, 0.4 and 0.4 nS (see details in Supplementary Table 2 and see Supplementary Fig. 1 for the active case).

Neuron_Reduce implemented on layer 5 pyramidal cell model with synapses
In Figure 1, Neuron_Reduce is implemented on a detailed compartmental model of a 3D reconstructed layer 5 pyramidal neuron from the rat somatosensory cortex (same model as in 5 ). This neuron consists of eight basal dendrites and one apical dendrite (shown in different colors) stemming from the soma. This neuron model has active membrane ion channels at both the soma and dendrites (see below). However, Neuron_Reduce first treats the modeled tree as passive by abolishing all voltagedependent membrane conductances, and only retaining the leak conductance.
Implementing Eqs. (1) - (11) for this cell produced a reduced, multi-cylindrical, passive model ( Fig. 1b, Step A) consisting of only 50 compartments rather than the 642 compartments in the detailed model.

Figure 1c
shows an example of four synapses located at different apical branches.
These synapses all have the same transfer resistance to the soma in the detailed tree.
Therefore, Neuron_Reduce maps these synapses to a single respective locus in the respective cylinder, such that their transfer resistance is identical in both models. In the reduced model, these synapses are merged into one "NEURON" process (red synapse in Fig. 1b). However, they retain their individual activation times (see Methods).

Figure 1d
compares the transfer impedance between a specific point in the apical tree (marked by 'd' in Fig. 1a and b) and the soma. By construction, for the passive case, the transfer resistance (for ω = 0) is equivalent for the respective loci in the detailed and the reduced model. This is indeed the case in Fig. 1d (left-most point on the x-axis), thus validating the implementation of the Neuron_Reduce analytic method. Note that although constructed using ω = 0, the similarity between the detailed and reduced model also holds for higher input frequencies. However, for ( around 10 -100 Hz, the transfer impedance from d-to soma (and vice versa, due to the reciprocity theorem for passive systems 43 ) is somewhat larger in the reduced model (compare the red and black lines).
To test the performance of Neuron_Reduce on transient synaptic inputs (composed of mixed input frequencies), we sequentially activated the four synapses shown in Fig. 1c in both the detailed and the reduced models (see Methods and Supplementary Table   2). Figure 1e shows the close similarity in the composite somatic EPSPs between the two models, further validating that the mapping of the detailed model to the reduced model using ω = 0 provides satisfactory results for the passive case (see also Supplementary Fig. 2).

Accuracy and speed-up of Neuron_Reduce for nonlinear neuron models
To measure the accuracy of Neuron_Reduce for a fully-active nonlinear neuron model we ran a comprehensive set of simulations using the well-established case of the L5 pyramidal cell model 5 shown in Fig. 2a (same cell as in Fig. 1). This neuron model includes a variety of nonlinear dendritic channels including a voltage-dependent Ca 2+ "hot spot" in the apical tuft (schematic yellow region in Fig. 1c) and an Na + -based spiking mechanism in the cell body. We randomly distributed 8,000 excitatory and 2,000 inhibitory synapses on the modeled dendritic tree (the synaptic parameters are listed in Supplementary Table 2 between the two spike trains peaked around zero (Fig. 2c), and the inter-spike interval distributions of the two models were similar (Fig. 2d).
The full range of responses to a random synaptic input for the two models was explored by varying the firing rate of the excitatory (AMPA-and NMDA-based) synapses and measuring the degree of similarity between the firing rates of the two models, which indicated a good fit between the two (Fig. 2e). We used the SPIKE-synchronization was high for all output frequencies but was low for low output frequencies when the synapses consisted of an NMDA component, although improving significantly for output frequencies above ~7 Hz (see Discussion). Figure 2g shows the SPIKEsynchronization as a function of the firing rate of the detailed model, for active and passive dendrites and with/without NMDA-based synaptic conductance, demonstrating again that when NMDA synapses are involved, the performance of the reduced model is low for low output rates. We also tested other spike trains similarity metrics 46,47 ( Supplementary Fig. 3) and found comparable results to those shown in Fig. 2. We have also analyzed the performance of Neuron_Reduce on two additional patterns of synaptic input. In one case, the synaptic input was activated in an oscillatory manner at different frequencies (see Methods). In these cases, the spike synchronization measure ranged between 0.75 -1 (Supplementary Fig. 4a, b). In the other case the synaptic input was taken from a spontaneously active Blue Brain circuit 17 (see Methods). In this case the spike synchronization measure was 0.708 ( Supplementary Fig. 4c).  Table 2 for synaptic parameters). Excitatory synapses were activated randomly at 5Hz and the inhibitory synapses at 10 Hz. This detailed model consists of a dendritic Ca 2+ " hot spot" (as in Fig. 1c) and an Na + spiking mechanism at the cell body. Scale bar 100 µm. b An example of the voltage dynamics at the soma of the detailed model (black trace) and the reduced model (red trace); spike times are represented by the black and red dots above the respective spikes. c Cross-correlation between spikes in the reduced versus the detailed models. d Inter-Spike Interval (ISI) distributions for the two models. e Output firing rate of the reduced (red) versus the detailed (black) models as a function of the firing rate of the excitatory synapses. Grey dots represent the case shown in b. f SPIKE-synchronization measure between the two models as a function of the firing rate of the detailed model for the case of only AMPA (blue) and AMPA + NMDA synapses (orange). The performance of the reduced model with NMDA-synapses was lower for low output frequency but improved significantly for output frequencies above ~7 Hz (see Discussion). g SPIKE-synchronization between the detailed and the reduced models as a function of the firing rate of the detailed model, for active and passive dendrites, and with/without NMDA-based synaptic conductance.
We compared the performance of our reduction method to two other reduction approaches, one of which was Rall's "equivalent cable" reduction method 31,48 . The other method maps all the dendritic synapses to the somatic compartment, after computing the filtering effect of the dendritic cable for each synapse 34 (see Methods).
Neuron_Reduce outperformed both these reduction methods (Supplementary Fig. 5).  In Figure 4 we compared the dendritic voltage in the detailed and in the respective location in the reduced model. We found that: (i) the voltage transients could differ significantly in dendritic branches that are all mapped to the same compartment in the reduced model (e.g., compare the grey traces in the yellow compartments in Fig. 4b).
(ii) the average voltage trace of these different dendritic branches (black trace in Fig.   4b) is similar to the voltage in the respective compartment in the reduced model (red trace in Fig. 4b).
The implications of the latter finding for capturing highly nonlinear local dendritic events is elaborated in the Discussion.

Neuron_Reduce replicates key dendritic nonlinearities and computations
To determine the capabilities of the reduced models to support nonlinear dendritic phenomena and dendritic computations, we repeated two classical experiments in both the detailed and the reduced model of the L5 pyramidal cell shown in Fig. 1. The first simulated experiment started by injecting a brief depolarizing step current to the soma of the detailed model to generate a somatic Na + action potential (AP, black trace in Fig.   5a). This AP propagated backward to the apical dendrite, the BPAP (red trace in Fig.   5a1). Repeating the same current injection in the reduced model led to a similar phenomenon, but with a larger BPAP (Fig. 5b1). The detailed model also included a "hot region" with voltage-dependent calcium conductances in its apical dendrite (see also Fig. 1). Combining somatic current injection with synaptic-like transient depolarizing current injected to the apical nexus evoked a prolonged Ca 2+ spike in the distal apical dendrite (red trace at the apical tree) which, in turn, generated a burst of somatic Na + spikes (the BAC firing 4,5,50 , Fig. 5a2). Neuron_Reduce maps the nonlinear dendritic "hot" Ca 2+ region to its respective location in the reduced model (see Fig. 1 and Methods). Fig. 5b1-b2 shows that the exact same combination of somatic and dendritic input currents also produced the BAC firing phenomenon in the reduced model. However, the reduced model was somewhat more excitable than the detailed model; this resulted in a burst of three spikes with a higher frequency (and sometimes with an additional spike) in the reduced model (compare Fig. 5a2 to b2). Left, detailed L5 pyramidal cell model with nonlinear Ca 2+ " hot spot" (same model as in Fig. 2). a1 Injecting a depolarizing step current to the soma (0.95 nA for 8.5 ms) in the detailed model evoked a somatic action potential, AP (black trace) that propagated backward semi-actively into the apical tree (red trace). a2 Combining the somatic input with a transient synaptic-like current injection (0.95 nA peak value with 0.5 ms and 5 ms rise-time and decay-time respectively) to the "hot region" in the apical dendrite evoked a prolonged local Ca 2+ spike which, in turn, triggered a burst of two extra somatic Na + spikes (the BAC firing phenomenon 50 ). b1-b2 Same as in a1-a2, but for the reduced model. Scale bars for the respective morphologies are 100 µm.
The second simulated experiment attempted to replicate theoretical and experimental results reported in previous studies 1,51,52 . In these studies, several excitatory synapses were activated sequentially in time, on a stretch of a basal dendrite, either in the somato-dendrites (OUT) direction or vice versa (the IN direction). Rall showed that the shape and size of the resultant composite somatic EPSP depended strongly on the spatiotemporal order of synaptic activation; it was always larger and more delayed for the centripetal (dendrites-to-soma) than for the centrifugal (soma-to-dendrites) sequence of synaptic activation (this difference can serve to compute the direction of motion 51 ). It was shown that the difference in the resulting somatic voltage peak between these two spatio-temporal sequences of synaptic activation was enhanced when nonlinear NMDA-dependent synapses were involved and that it made it possible to discriminate between complex patterns of dendritic activation 52 .
To simulate these phenomena, twelve excitatory synapses were placed along one basal branch in the detailed model (red dots on the green basal tree, Fig. 6a1 Fig.   6a2). Neuron_Reduce merged these twelve synapses into five point processes along the respective cylinder (Fig. 6b1). We repeated the same experiment in the reduced model and found that the EPSP resulting from the IN direction was larger and delayed, with a similar EPSP waveform to that of the detailed model (see also Supplementary Fig. 2 and Discussion). Next, an NMDA component was added to the 12 simulated synapses; this resulted in larger somatic EPSP amplitudes in both directions (and both models) and a smaller difference in the peak timing between the different directions in both the detailed and the reduced model (compare Fig. 6a3 to b3). To generalize the impact of the spatio-temporal order of synaptic activation, we used a directionality index suggested in a previous study 52 . This measure estimates how different a given synaptic sequence is from the IN sequence by calculating the number of synaptic swaps needed to convert this given pattern into the IN pattern (using the bubble-sort algorithm, see Methods). We tested the EPSPs that resulted from different temporal combinations of synaptic activation (each having a different directionality index), both without (Fig. 6c1) and with an NMDA component (Fig. 6d1) the behavior of the two models was similar when the somatic voltage in the two models was subtracted by the peak value obtained in the OUT direction ( Fig. 6c2 and 6d2).
Then, the difference between the reduced and the detailed models was, on average, only 0.11 ± 0.43 mV for the AMPA-only case and 0.35 ± 0.43 mV for the AMPA + NMDA case. Thus, although the detailed and the reduced models differ to a certain extent (see Discussion), the capability of the reduced model to discriminate between spatiotemporal patterns of synaptic activation is similar to that of the detailed model.

Neuron-Reduce can be implemented successfully on a variety of neuron models
We next tested the utility of Neuron_Reduce on 13 different neuron models from different brain regions (Fig. 7). Four models were obtained from the Blue Brain database 17,53 : L6 Tufted Pyramidal Cell, L4 Double Bouquet Cell, L4 Spiny Stellate Cell and L5 Martinotti cell from the rat somatosensory cortex. Two additional models were obtained from the Allen Institute cell type database 11 : an L4 Spiny cell and an L1 Aspiny Cell from the mouse visual cortex. Medium spiny neuron from the mouse basal ganglia 54 ; two rat thalamocortical neurons 55 ; Golgi cell from mouse cerebellar cortex and one inhibitory hippocampal neuron from the rat 56 . We also took two additional neuron models from our laboratory: rat L2/3 large basket cell 57 and a model of a human L2/3 pyramidal cell from the temporal cortex 58 . All these models were based on 3D reconstructions and were constrained by experimental recordings (see Supplementary   Table 2 for details on the various neuron models and input parameters).
Neuron_Reduce successfully generated a reduced model for all these different cell types, with highly faithful response properties in all cases (Fig. 7). Three examples with their respective morphologies for the detailed and reduced models are shown in Fig.   7a-c. For a given input we measured the spiking activity of the detailed and reduced models ( Fig. 7d-f) and calculated the corresponding SPIKE-synchronization values.
For the L6 tufted PC model (Fig. 7a, d), the L2/3 large basket cell model (Fig. 7b, e) and the L6 double bouquet model (Fig. 7c, f) Fig. 7d-f). The SPIKE-synchronization values for additional inputs, and for the other 10 neuron models and their corresponding reduced models, are shown in Fig.   7g. We have also tested the performance of Neuron_Reduce and the variability of the SPIKE-synchronization measure using eight neocortical neuron types, with 11 cell models per type taken from the Blue Brain cells dataset 17,53 . Supplementary Figure 6 shows that, for all cells the SPIKE-synchronization measure remains similar to that found in Figure 7, mean values per cell type ranging between 0.43-0.86, increasing with output frequency of the modelled cell.

Neuron_Reduce: theoretical underpinning and advantages
Neuron_Reduce is a new tool for simplifying complex neuron models while enhancing their simulation run-time. It analytically maps the detailed tree into a reduced multicylindrical tree, based on Rall's cable theory and linear circuit theory (Fig. 1). The underpinning of the reduction algorithm is that it preserves the magnitude of the transfer impedance "# $,1 (()" from each dendritic location, j, to the soma (the dendro-somatic direction, Eqs. (1) - (11) in Methods). Since in linear systems, it holds that "# $,1 (()" = "# 1,$ (()", for passive dendritic trees it also preserves the transfer impedance in the somato-dendritic direction (e.g., current injection at the soma will result in the same voltage response at the respective sites in the detailed and reduced models 59 ).
Note that dendritic voltage transients (e.g., synaptic potentials) contain a range of frequencies, (. We however had to select one frequency to use for the mapping of the detailed-to-the-reduced tree. Consequently, we examined a whole range of possible ( values for this mapping. Conveniently we found that ω = 0 is the preferred frequency for generating the reduced model (namely, when the mapping from detailed to reduced model is performed based on the transfer resistance "# $,1 (( = 0)" = "3 $,1 "), see voltage. In other words, when using the transfer resistance for our mapping procedure, we preserved the total charge transfer (which in our case, was proportional to the voltage time integral) from the synapse to the soma (and vice versa), but not e.g., the EPSP peak value.
Neuron_Reduce was proven to be accurate in replicating voltage dynamics and spike timing for a large regime of input parameters and a variety of neuron types (Fig. 7,   Supplementary Figure 6 and Supplementary Table 2). This claim is based on using several metrics for assessing the quality of the performance of the reduced model ( Supplementary Fig. 3). Moreover, we implemented STDP rule on synaptic weights in both the detailed and reduced model and found great similarity in the weight trajectories of the respective synapses (not shown). Neuron_Reduce is straightforward to use, it is fast and generally applicable, thus enabling its implementation on any neuron morphology with any number (even tens of thousands) of synapses. One key advantage of Neuron_Reduce is that it retains the identity of individual dendrites and synapses and that it maps dendritic nonlinearities to their respective loci in the reduced model, hence preserving local excitable dendritic phenomena and therefore maintaining nonlinear dendritic computations. Neuron_Reduce also preserves the passive cable properties (Rm, Ra, and Cm) of the detailed model, thus preserving synaptic integration and other temporal aspects of the detailed model. We have also examined Neuron_Reduce for the case of a gap junction connecting dendrites of two L2/3 large basket cells. As Neuron_Reduce preserves the transfer resistance from the location of the synapses (in this case the gap junction) to the soma and vice versa, one expects that the coupling coefficient between the two connected cells will be preserved in the reduced models, after mapping the gap junction to its appropriate location in the reduced model. This is indeed what we have found (not shown).
Neuron_Reduce enhances the computational speed by a factor of up to several hundred folds, depending on the simulated morphology and the number of simulated synapses ( Fig. 3 and Supplementary Table 1). This combination of capabilities, together with its user-friendly documentation and its public availability, make Neuron_Reduce a promising method for the community of neuronal modelers and computational neuroscientists, and for the growing community interested in "biophysical deep learning".
For a large number of synapses and complex morphologies, the run-time of Neuron_Reduce models can be accelerated by up to 250 fold as compared to their respective detailed models ( Fig. 3 and Supplementary Table 1). This is achieved in two associated steps. First, the algorithm reduces the number of compartments of the neuron model; e.g., for the reconstructed tree in Fig. 1, it reduced the number of compartments from 642 to 50. Then, synapses (and ion channels) that are mapped to the same electrical compartment in the reduced tree (because they have similar transfer resistance to the soma) are merged into one point process in NEURON. Each of these steps on its own has a relatively small effect on the run-time. However, when combined, a large (supra-linear) improvement in the computational speed is achieved (Supplementary Table 1). This is because at each time step, NEURON computes both the voltage in each electrical compartment as well as the currents and states of each point process and membrane mechanism (synapses and conductances). Reducing the number of compartments in a model decreases the number of equations to be solved and the number of synapses to be simulated (due to the reduced number of compartments, a larger number of synapses are merged together). Importantly, merging synapses preserves the activation time of each synapse. Note, however, that in its present state, Neuron_Reduce cannot merge synapses with different kinetics.

Other reduction methods
Several other reduction methods for single neurons have been proposed over the years 12,[34][35][36][37][38][39] . Most are not based on an analytic underpinning and thus require handtuning of the respective biophysical and morphological parameters. In addition, most of these methods have not been examined using realistic numbers of dendritic synapses and are incapable of systematic incorporation of dendritic nonlinearities. In most cases, their accuracy has not been assessed for a range of neuron types (but see 41 ). Many of these methods are not well-documented, thus making it hard to compare them directly with Neuron_Reduce. Nevertheless, we did compare the performance of Neuron_Reduce to two other reduction methods and showed that it outperformed them ( Supplementary Fig. 5).

Limitations of Neuron_Reduce and future work
It should be noted that although the transfer impedance from a given dendritic locus to the soma is preserved in the reduced model, the input impedance at that locus is not preserved (is lower) in the reduced model. Consequently, the conditions for evoking local dendritic events, and the fine details of these events are not identical in the detailed and the reduced models (e.g., compare Fig. 5a1 and a2 to Fig. 5b1 and b2 and see Fig.   4). Indeed, if there were highly local dendritic Na + spikes (as in 60 ) then Neuron_Reduce will not capture them, as this local dendritic spike will be averaged-out in the respective lumped cable. Similarly, because the local voltage response to a current injection in the dendrite depends on the dendritic impedance, the local synaptic responses are somewhat different in the detailed versus the reduced cases, especially when voltagegated ion channels (such as NMDA-dependent synaptic channels) are involved. In fact, when large dendritic NMDA signals are involved, the resultant somatic EPSPs are expected to be different in the detailed as compared to the reduced model, as is the case in Figures 2 and 6. Indeed, if one insists on preserving highly local nonlinear dendritic events, then the full dendritic tree should be modeled.
Despite these local differences, the reduced model for L5PC did generate a local dendritic Ca 2+ spike in the cylinder representing the apical dendrite and was able to perform an input classification task (enhanced by NMDA-conductance), as in the detailed tree (Figures 5 and 6). Moreover, when embedded in large circuits, individual neurons are likely to receive semi-random dendritic input, rather than a clustered input on specific dendrites. For such inputs, the reduced models generated by Neuron_Reduce capture most of the statistics of the membrane voltage dynamics as in the detailed model (Figures 2 and 7 and Supplementary Figures 6 and 4). where and 6 = √1 + >(K where K is the membrane time constant, 3 C L C . From Eq. (1), the input impedance at X = 0 is .
From Eqs. (2) and (4) we get and thus "# $,$ (()" = f 2 A6 from which we compute the diameter, d, for that cylinder Equations (1)  This distal dendritic locus, x, was mapped to the distal end, X = L, of the corresponding cylinder. We then used Eqs. (1) - (11) to calculate the unique respective cylinder for each stem dendrite.
In order to map synapses from the detailed model to the reduced one we computed, for each synapse at location j in the detailed model, "# $,1 (()", and then mapped this synapse to the respective location y in the reduced model, such that "# $,k (()" = "# $,1 (()". This reduced model is then compartmentalized into segments (typically with spatial resolution of 0.1l, see Fig. 3b). We then merged all synapses with identical kinetics and reversal potential, that are mapped to a particular segment, onto a single point process object in NEURON (Fig. 1, step B). These synapses retain their original activation time and biophysical properties through the connection of each of their respective original NetStim objects to the single point process that represents them all (each of these connections was mediated by the synapse's original NetCon object). As shown in Supplementary Table 1, this step dramatically reduced the running time of the model. We note that all the results presented in this study were obtained using ( = 0 in Eqs. (1) - (11), since running the same simulations with ( = 0 provided the best performance (see Supplementary Fig. 7). However, ( is a parameter in the algorithm code and can be modified by the user. Note also that "# $,$ (()", "# $,1 (()" and "# $,& (()" were analytically computed for each original stem dendrite using the NEURON impedance tool 64 .

Neuron models used in the present study
To estimate the accuracy of the reduction method, we ran 50-second simulations on cell morphologies of different types, in both the reduced and detailed models. Models of 13 neurons were used in this study; their details are available in Supplementary Table 2.
For each of the models, we distributed 1,250 -10,000 synapses on their dendritic trees.
Eighty percent of the synapses were excitatory, and the rest were inhibitory. The synaptic conductances were modeled using known two-state kinetic synaptic models 17 .
For simplicity, we did not include synaptic facilitation or depression. All models had one type of GABAA -based inhibitory synapses and either AMPA -or AMPA+ NMDA -based excitatory synapses. The synaptic rise and decay time constants were taken from various works cited in Supplementary Table 2. When no data were available, we used the default parameters of the Blue Brain Project synaptic models 17,53 . Inhibitory synapses were activated at 10 Hz, whereas the activation rate of the excitatory synapses was varied to generate different output firing rates in the range of 1 to 20 Hz (Figures   2, 3, 4, 7 and Supplementary Figures 3-7); the values used for each model are listed in Supplementary Table 2. In all models except Supplementary Fig. 4, synaptic activation time was randomly sampled from a homogenous Poisson process. In Supplementary Fig. 4 a, b the activation time was sampled from an inhomogeneous Poisson process with a time-dependent intensity l(m) = n * sin(m * q * 2A) + 1, where t is time in seconds, r is the firing rate of the synapse and f is the oscillation frequency.
In Supplementary Fig. 4c we extracted a single Layer 5 Thick-tufted Pyramidal Cell with an early bifurcating apical tuft (L5_TTC2; gid 75586) from active blue brain microcircuit 17 . Calcium and potassium concentration were set to 1.23 mM and 5.0 mM, respectively. The synaptic activation from the microcircuit was replayed to this detailed model and also to its respective reduced model. Synaptic depression and facilitation were disabled, and the synapse time constants, which varied in the microcircuit, were set to their mean value (the decay time constant was set to 1.74 ms and 8.68 ms for AMPA and GABAA respectively; the rise time constant for GABAA was set to 4.58 ms) all other variables were as in the blue brain simulations.

Estimating the accuracy of the reduced models
Cross-correlations were calculated between the spike trains of the detailed and the reduced models. The window size was 500 ms, and the bin size was 1 ms. The resulting cross-correlations were normalized by the number of spikes in the detailed model (Fig.   2c). Inter spike intervals (ISI) were binned in windows of 21 ms to create the ISI distribution in Fig. 2d.
SPIKE-synchronization measure is a parameter-and scale-free method that quantifies the degree of synchrony between two spike trains 44 . SPIKE-synchronization uses the relative number of quasi-simultaneous appearances of spikes in the spike trains. In this study, we used the Python implementation of this method 65 . To allow comparison to the literature, Supplementary Figure 3 depicts three additional metrics against which to compare the performance of the detailed and the reduced models: Trace accuracy 39 , ISI distance 44 , and Γ coincidence factor 66 .

Comparison to other reduction algorithms
We compared Neuron_Reduce to two classical reduction algorithms ( Supplementary   Fig. 5): 1. Equivalent cable using the d3/2 rule for reduction. Rall and Rinzel 32,33 showed that for idealized passive dendritic trees, the entire dendritic tree can be collapsed to a single equivalent cylinder that is analytically identical (from the point of view of the soma) to the detailed tree. However, neurons do not have ideal dendritic trees, mostly because dendritic terminations typically occur at different electrotonic distances from the soma.
Nevertheless, it is possible to collapse any dendritic tree using a similar mapping (Rall's "d3/2 rule") as in the idealized tree; this will provide an "equivalent cable" (rather than an "equivalent cylinder") with a varying diameter for the whole dendritic tree (see details in Rall et al. 48 ). The electrotonic distances to the soma of synapses and nonlinear dendritic mechanisms were computed in the original model and then each synapse and mechanism was mapped to the corresponding segment in the "equivalent cable" preserving the electrotonic distance to the soma as in the original tree.

Mapping all synapses to the soma. Another recent reduction scheme was proposed
where all dendritic synapses are mapped, after implementing cable filtering for each synapse, to the somatic compartment 34 . Here we used a modified version of this method. We used Neuron_Reduce to generate a multi-cylindrical model of the cell as in Fig. 1b. Then, all the synapses in the original tree were mapped to the model soma.
To account for dendritic filtering, for each synapse, we multiplied the original synaptic conductance, r skt , by the steady-state voltage attenuation factor from the original dendritic location, j, of the synapse to the soma. Specifically, r skt * = r skt * "# $,1 " "# $,$ " = r skt * , $,1 , $,$ where r skt * is the new synaptic weight for synapse j when placed at the soma of the reduced model.

Discriminating between spatio-temporal patterns of synaptic activation
In Figure 6, twelve synapses, placed at 25 µm intervals, were distributed on a stretch of one basal dendrite. The peak AMPA conductance per synapse was 5 nS. In cases where the synapses also had an NMDA component, the NMDA-based peak conductance was 3.55 nS. The synapses were activated in a specific temporal order with a time delay of 3.66 ms between them. This resulted in an input velocity of 7 µm/s for the sequential IN and OUT patterns in Fig. 6. In addition, the temporal order of synaptic activation was randomized and scored according to the directionality index 52  The Neuron_Reduce algorithm is publicly available on GitHub (http://github.com/orena1/neuron_reduce).