Inverse correlation processing by neurons with active dendrites

In many neuron types, the dendrites contain a significant density of sodium channels and are capable of generating action potentials, but the significance and role of dendritic sodium spikes are unclear. Here, we use simplified computational models to investigate the functional effect of dendritic spikes. We found that one of the main features of neurons equipped with excitable dendrites is that the firing rate of the neuron measured at soma decreases with increasing input correlations, which is an inverse relation compared to passive dendrite and single-compartment models. We first show that in biophysical models the collision and annihilation of dendritic spikes causes an inverse dependence of firing rate on correlations. We then explore this in more detail using excitable dendrites modeled with integrate-and-fire type mechanisms. Finally, we show that the inverse correlation dependence can also be found in very simple models, where the dendrite is modeled as a discrete-state cellular automaton. We conclude that the cancellation of dendritic spikes is a generic mechanism that allows neurons to process correlations inversely compared to single-compartment models. This qualitative effect due to the presence of dendrites should have strong consequences at the network level, where networks of neurons 1 . CC-BY 4.0 International license peer-reviewed) is the author/funder. It is made available under a The copyright holder for this preprint (which was not . http://dx.doi.org/10.1101/137984 doi: bioRxiv preprint first posted online May. 14, 2017;


Introduction
Increasing evidence shows that nonlinear integration of synaptic inputs in dendrites is crucial for the computational properties of neurons. A major role in the integration is played by dendritic spikes: regenerative currents through Na + , Ca 2+ or NMDAr channels. The first evidence of dendritic spikes came from field recordings [1,2,3,4,5], corroborated by the intracellular recordings [6,7,8]. The repertoire of techniques was further enlarged by patch clamp [9,10,11,12,13] and optical methods. Calcium imaging allowed for the direct observation of calcium spikes [14,15,16,17,18], and glutamate uncaging and voltage sensitive dyes led to the discovery of NMDA spikes [19,20].
Dendritic spikes allow for more subtle integration of synaptic input than in a passive dendrite. A single dendritic branch can act as a coincidence detector, generating a spike when exposed to synchronized input [21,22,20]. The propagation of dendritic spikes generated in the distal part of dendritic tree can be gated by synaptic input in the proximal region, as was shown for hippocampal CA1 pyramidal neurons [23], and L5 pyramidal neurons [24]. After the initiation of a dendritic spike, sodium channels inactivate and the branch switches into a refractory state which crucially affects integration [25]. Backpropagating action potentials also play an essential role in spike time-dependent plasticity [26,27,28], and the participation of local dendritic spikes has been implicated in long-term potentiation [29,30].
Sodium spikes can propagate in many cell types: neocortical pyramidal cells [31,9,32,33,34], hippocampal CA1 and CA3 pyramidal cells [35,36,37,38,39], interneurons [40] or thalamic neurons [41]. Models predicted that, in vivo, the presence of synaptic background activity should greatly enhance the initiation and propagation of dendritic spikes [42,43]. This suggests that there can be a heavy traffic of dendritic spikes, as indeed found in recent dendritic recordings in awake animals [44]. Therefore interactions between dendritic spikes likely play an important role in dendritic integration.
Here, using different types of computational models, we show that interactions between dendritic spikes can change the way the correlation of synaptic input affects the firing rate of the cell. For single-compartmental neurons, an increase of correlation of synaptic input is known to cause an increase of fir-ing rate [45,46]. Because the propagation of a dendritic spike can be stopped by the refractory periods caused by the other dendritic spikes, the relation between correlation and firing rate can be different for multi-compartmental neurons.

Results
We first show, using biophysical models of neurons with active dendrites, the phenomenon of inverse correlation processing. Next, we investigate this phenomenon using active dendrites modeled by the integrate-and-fire model. Finally, we show that this phenomenon is also present in simplified discretestate models.

Biophysical models of inverse correlation processing
In all models studied here, we considered neurons subject to in vivo-like synaptic activity. In particular, we aimed at investigating the effects of synaptic noise correlations on the firing of the cell. The synapses were located on the dendritic compartments and were of two types: excitatory AMPA synapses (E e = 0 mV) and inhibitory GABA A synapses (E i = −75 mV), with four times more excitatory synapses than inhibitory ones.
We considered the case where the pre-synpatic spikes triggering the excitatory and inhibitory inputs are correlated. There might be two biological sources of such correlation. First of all, in cortical networks neurons are known to share some of their synaptic inputs which inevitably makes their spikes correlated. Since this type of correlation is common to larger population of neurons, we call it the global correlation [47]. In addition, a single axon can create several synapses on close-by segments of the dendrite activating them at the same time [48], so that the correlation of the input spikes within the same dendritic segments (local correlation) can be higher than the global correlation.
To generate such locally and globally correlated synaptic inputs, we used a previously proposed algorithm [49]. From a global spike train G obtained from Poisson process, we draw random spikes with probability r G . The spike train L i thus obtained corresponds to all synapses on a single dendritic compartment. Finally, the spike train S j i for a synapse situated on compartment i is built by drawing spikes from L i with local probability r L [ Fig. 1]. In this model, the ratio of shared spikes between synapses situated on the same compartment is equal to c L = r L and the ratio of spikes shared by synapses situated on different compartments is equal to c G = r L r G . Note that c L is always greater than or equal to c G . At the end, we added a jitter to each spike time to obtain desynchronized input. The jitter is drawn from exponential distribution with a standard deviation τ j with equal probability for positive and negative values. To obtain the time-dependent synaptic conductances, we convolve the resulting spike trains with an exponential function reflecting the change of synaptic conductance due to a single spike.
We first apply this model of correlations to simulate in vivo-like activity in a biophysical model. We used the well-known FitzHugh-Nagumo (FN) model, which is a simplified two-variable system derived from the Hodgkin-Huxley model [50,51]. See Methods for equations of the model and its parameters. We compared a point neuron with a multi-compartmental model [ Fig. 2]. The firing rate increases with correlation for a point neuron, but for a multi-compartmental neuron, the relationship is inversed. In the remaining of the paper, we will refer to this property as inverse correlation processing. This phenomenon is due to the fact that, with active dendrites, the correlation increases the local dendritic spike initiation, but these spikes cancel (due to the refractory period following a dendritic spike), and thus fewer spikes reach the soma. In the following, we investigate this phenomenon in more detail.

Integrate and fire models of excitable dendrites
To investigate the details of the local dendritic integration mechanism with dendritic spikes, we consider an integrate-and-fire type model. In particular, we used the Adaptive Exponential (AdEx) model [52], which can simulate various phenomena like adaptation or bursting mechanisms. In this model, we can freely adjust the refractory period and the adaptation parameters. We extended the model with axial currents and spike waveform to adapt it for modeling multi-compartmental dendrites. See Methods for the equations of the multi-compartmental AdEx model. After this modification, the model can produce a spike propagating across the dendrite, but the propagation can be stopped by the tunable refractory period caused by other spikes.
The velocity of spikes depends on the passive properties of the dendrite, as well as on the specific AdEx parameters. We obtained an approximate analytical formula for the speed of the spike (Supplementary materials). The For each run, we measured the number of dendritic spikes that reach the soma and we asked how this quantity is affected by the ratio of shared spikes. We performed the scan over local correlation c L and global correlation c G for different weights of synapses. While increasing synaptic weights, we decreased the frequency of inputs spikes, such that the total current injected to the dendrite was constant. We found that for the point neuron the number of spikes in the soma increased with the correlation of synaptic input for all weights while for the multi-compartmental neuron we can see the opposite behavior [ Fig. 3]. We compared the behavior of a multi-compartmental neuron equipped with an active dendrite and a neuron with a passive dendrite [ Fig. 4 ]. For the neuron with a passive dendrite the frequency of somatic spikes monotonically increases with the correlation of synaptic input. For the neuron with an active dendrite the frequency of somatic spikes initially increases with the correlation, but starts to decrease for larger values of the correlation leading to an existence of a global maximum of the firing rate.
In the multi-compartmental neuron with an active dendrite, the dependence on correlation is strikingly different for different synaptic weights [ Fig. 5]. For low synaptic weights, we observed an increase in the number of spikes reaching the soma for an increased correlation. For larger weights, in the lowcorrelation regime we observed an increase in number of spikes, but further increase of correlation decreased the number of spikes reaching soma.
The reason for this effect is that for low weights, most synaptic potentials remain sub-threshold and summed synaptic inputs from single or multiple compartments are necessary to trigger a single spike. Hence, a higher input correlation leads to an increase in the probability of generating a spike. On the other hand, when we increase synaptic weights few synaptic inputs may trigger a spike, whose propagation to the soma is limited by the effects of the refractory period and collisions. In this regime, a higher correlation means higher probability of stopping spike propagation by the refractory periods of preceding spikes. Therefore, we observe a decrease in the number of spikes reaching the soma.
The length of refractory period plays crucial role in the number of spikes reaching soma as a function of input correlations. We performed the scan over local c L and global c G correlation for different lengths of refractory period [ Fig. 6]. For longer refractory periods, the probability of stopping spikes by refractory periods induced by other dendritic spikes is higher. This effect is  much more visible for correlated synaptic inputs than for uncorrelated ones [ Fig. 6].
The response for correlated inputs depends also on the intensity of synaptic bombardment. For low intensity, we observe an increase in the number of spikes reaching soma with correlation. This is due to low traffic of dendritic spikes and low probability of interference between propagations of spikes. For higher intensities, we can see the inverse relation [ Fig. 7].
Simplified discrete-state model Finally, we designed a mathematical model, called discrete-state model, aiming at grasping the core mechanism of spike propagation / annihilation. Compared to the more detailed models described so far, this one has much smaller simulation time and allows to easily modify parameters independently (e.g. propagation speed). The simulation is exact in that it is not based on a (time) discretization of the propagation of dendritic fronts. Finally, it allows to see whether the annihilation of spikes is enough to explain decrease of mean somatic spikes number. More precisely, suppose that we are given a set of spatio-temporal synaptic inputs S = {(t i , x i )} N i=1 and assume that each input (t i , x i ) produces two contra-propagating fronts (e.g. dendritic spikes). We also assume that S is generated with the procedure described above. In Figure 8, the synaptic inputs are the blue dots and the lines, the propagating fronts. Starting from S, we build the set of fronts and annihilation events recursively (see Methods). Figure 8 shows the networks of propagating fronts for different propagating speeds but same inputs S. One can see that increasing speed increases the dendritic spikes reaching the soma (x = 0). Intuitively, this occurs because a front has less chance to meet another front for higher propagation speeds. Finally, in Figure 9 we present the dependency of the mean number of spikes reaching the soma as a function of the ratio of shared spikes. In this scenario, each synaptic input triggers a dendritic spikes, so the results corresponds to the case of large synaptic weights of the multi-compartmental AdEx model [ Fig. 5].

Discussion
In this paper, we have investigated the behavior of excitable neuronal dendrites using different models. The main and remarkable property that we    explored is the ability of the dendrite to initiate and propagate dendritic spikes which gives a very particular property of inverse correlation processing, which itself cannot be obtained with a single-compartment model. We have shown this property in different models ranging from biophysical models, integrate-and-fire models, and simplified discrete models. We discuss below possible implications of this finding, and ways to test it experimentally.
Our model challenges the well-known property that the presence of correlations in the pre-synaptic activity augments the firing of neurons has been investigated extensively in point-neuron models [45,46]. We find here that neurons equipped with excitable dendrites can operate in a regime which is opposite to the point neuron: dendritic neurons can be set to decrease their firing when correlations are present in their afferent activity. Furthermore, this inverse correlation processing mode was dominant in models where dendrites generate fast spikes (such as Na + spikes) whenever the intensity of synaptic bombardment was high enough. Recent experiments [44] provided strong evidence that Na + dendritic spikes are much more frequent than so-matic spikes in awake animals (10 times more frequent during exploration), so there is a serious possibility that neurons can be set to this inverse correlation processing mode. The inverse responses to synaptic input correlation shown here can be present also in individual branches of dendritic tree exposed to strong synaptic bombardment. These responses are strongly affected by synaptic weights, which should be understood relatively to dendritic spike threshold.
We demonstrated the consequences of active dendritic conductances and the inverse correlation phenomenon on a simplified integrate-and-fire-type model (AdEx). Although this type of models is usually applied to point neurons without any spatial extent, we extended it over multiple compartments connected through axial currents such that propagating spikes were observed. Our motivation for using the AdEx model is, first, that the AdEx model is general approach to model single neurons and it can reproduce many neuronal types [53]. Simpler models such as leaky integrate-and-fire can be obtained as a special case of AdEx by specific choice of the model parameters. Secondly, the AdEx allows one to simulate different intrinsic properties such as bursting and adaptation. This allows us to later extend the present model by including such mechanisms in dendrites. Third, it is compatible with neuromorphic hardware, which have implemented neurons equipped with the AdEx mechanism [54], and which is presently extended to include neurons with excitable dendrites [55]. In addition, to allow for analytical treatment of the dendritic spikes, we further simplified the model to a discrete state model capturing the dynamics of the collisions of dendritic spikes. The full mathematical treatment of this model will be demonstrated in another report.
What are the consequences of such inverse correlation processing? One obvious effect will be to cancel correlations at the network level, because the neurons subject to correlated input will fire less, and neurons with decorrelated input will be the ones dominating the population, and will most likely shift the population towards decorrelated firing. Such dynamics should be examined by future network models with neurons equipped with excitable dendrites.
Finally, it is important to test the inverse correlation processing experimentally. Dynamic-clamp was used previously to investigate the fact that neurons increase their synchrony [56], but unfortunately, this technique cannot be used to test the present mechanism, because dynamic-clamp emulates inputs in the soma (the site of the recording electrode), while it is impor-tant that the inputs are dendritic. One possibility would be to use 2-photon imaging with glutamate uncaging [57] to simulate inputs with controlled frequency and synchrony. Another possibility would be to use voltage-sensitive dye imaging of dendrites in vitro [58] combined with controlled network activity in the slice, to directly monitor the genesis, propagation and possible collision of dendritic spikes.

FitzHugh-Nagumo model
The FN model is described by the following set of equations: where V is the membrane voltage, f (V ) = V − V 3 is a nonlinear function that allows a regenerative activity through a positive feedback, and W is a recovery variable which provides a slow negative feedback. a, b and c are constant parameters (for a recent overview, see [59]).
The dendrite was composed of connected FN compartments, such that axial current was flowing between compartments. The parameters were uniform along dendrite, a = 0.1, b = 0.25, c = 1 and axial resistance r i = 0.1.

m-AdEx model
Adaptive exponential integrate-and-fire model (AdEx) [52] is a synthesis of two improvements of leaky integrate-and-fire (LIF) model: 1. enrichment with spike triggering depolarizing currents [60] providing smooth spike initiation zone in place of fixed threshold of LIF model, 2. enrichment with neuronal adaptation represented by adaptation current governed by additional differential equation [61,62].
The model was introduced at first for a single-compartmental neuron. For a multi-compartmental model a membrane voltage V (t, x) and an adaptation current w(t, x) are governed by a set of partial differential equations: where c m is a specific capacitance, d is a diameter of a dendrite, r i is an intracellular specific resistance, g L is a leak conductance, E L is a leak reversal potential, V T is a spike threshold and ∆ T is a slope factor. g e (t) and g i (t) are synaptic excitatory and inhibitory conductances, E e and E i are synaptic reversal potentials. These equations were obtained by combining the classic cable equations [64] with the AdEx model.
The adaptation current w is governed by a second differential equation: where a represents subthreshold adaptation and τ w is the adaptation time constant. This set of equations is solved numerically by discretization of the equations over multiple dendritic compartments. When a membrane potential of any compartment is near threshold V T depolarizing exponential current surpass other currents and membrane voltage quickly tends to infinity. Whenever the membrane potential crosses peak value V p we detect a spike and compartment enters into refractory period in which voltage is repolarized according to equation: The exponential decay characteristic of this equation models the falling phase of an action potential. After emitting a spike the compartment stays in the refractory state during time t r . At the time of the spike the adaptation current is incremented according to the equation where the parameter b represents spike-triggered adaptation.

The refractoriness for detailed model and for m-Adex model
In a simulation of the detailed model of L5 pyramidal neuron [63] we can see how the initiation of the first spike can prevent the propagation of second spike. The active propagation can be impossible even 20 ms after the initiation of the first spike, [Fig. 10]. The time delay needed for propagation of the second spike varies from branch to branch and depends on the density of channels and on the diameter of the dendrite (the propagation of the second spike is facilitated in the retrograde direction because of decrease of dendrite's diameter).
The refractoriness of a detailed model can be mimicked in m-AdEx model. In Fig. 11, we can see how refractoriness affects propagation of dendritic spike. EPSPs which occur within refractory period can propagate only passively to soma.

Correlation measure for spike trains
To measure the correlation between two spike trains S i (t) and S j (t) we applied the cross-covariance function CCVF ij (s)

Simulation of the discrete-state model
Given a set S of somatic inputs, we describe our procedure to re-construct the propagating dendritic fronts leading to somatic spikes. Take the latest element s 1 = (t 1 , x 1 ) ∈ S. We look for possible annihilation events with its outgoing front, i.e. going away from the soma located at x = 0. Hence, we select all inputs S ↑ (s 1 ) in S \ {s 1 } which are located in the cone above s 1 defined as S ↑ (s 1 ) ≡ {(t, x) : x ≥ ±(t − t 1 ) + x 1 }. Consider the somatic input s 2 ∈ S ↑ (s 1 ) closest to the line x = −(t − t 1 ) + x 1 : it will produce an annihilation event with s 1 . Perform the same for s 2 recursively until possible. This gives a chain of annihilation events (s 1 , s 2 , ..., s k ). We do the same for possible annihilation events with the ongoing front from s 1 , i.e. going toward the soma located at x = 0. It gives another chain of events (s −p , ..., s 0 , s 1 ). Hence, we found a chain of annihilation events C(s 1 ) = (s −p , ..., s k ): it produces exactly one somatic spike. We continue the same procedure with S \ C(s 1 ) until S = ∅. The overall procedure runs in O(n 2 ) where n is the cardinal of S.
Supporting information S1 Estimation of the speed of spike To estimate velocity of spike we measure a time needed to rise voltage in the compartment in front of the spike. We take into account only exponential and leak currents and the current which flows from the compartment with a spike: where V p is an amplitude of a spike, V R is a resting membrane voltage, and l is an axial length of compartment. From this equation we can calculate the time needed to bring the membrane voltage from resting value to peak value. The velocity of the spike is approximately the length of the compartment divided by this time.