Deep supervised learning with local weight updates requires that each neuron receive signals that can be used to determine its ‘credit’ for the final behavioral output. We explored the idea that the cortico-cortical feedback signals to pyramidal cells could provide the required information for credit assignment. In particular, we were inspired by four observations from both machine learning and biology:

1. Current solutions to credit assignment without weight transport require segregated feedforward and feedback signals (Lee et al., 2015; Lillicrap et al., 2016).

2. In the neocortex, feedforward sensory information and higher-order cortico-cortical feedback are largely received by distinct dendritic compartments, namely the basal dendrites and distal apical dendrites, respectively (Spratling, 2002; Budd, 1998).

3. The distal apical dendrites of pyramidal neurons are electrotonically distant from the soma, and apical communication to the soma depends on active propagation through the apical dendritic shaft, which is predominantly driven by voltage-gated calcium channels. Due to the dynamics of voltage-gated calcium channels these non-linear, active events in the apical shaft generate prolonged upswings in the membrane potential, known as ‘plateau potentials’, which can drive burst firing at the soma (Larkum et al., 1999; 2009).

4. Plateau potentials driven by apical activity can guide plasticity in pyramidal neurons in vivo (Bittner et al., 2015; Bittner et al., 2017).

With these considerations in mind, we hypothesized that the computations required for credit assignment could be achieved without separate pathways for feedback signals. Instead, they could be achieved by having two distinct dendritic compartments in each hidden layer neuron: a ‘basal’ compartment, strongly coupled to the soma for integrating bottom-up sensory information, and an ‘apical’ compartment for integrating top-down feedback in order calculate credit assignment and drive synaptic plasticity via ‘plateau potentials’ (Bittner et al., 2015; Bittner et al., 2017) (Figure 3A).

Figure 3

Download asset Open asset

As an initial test of this concept we built a network with a single hidden layer. Although this network is not very ‘deep’, even a single hidden layer can improve performance over a one-layer architecture if the learning algorithm solves the credit assignment problem (Bengio and LeCun, 2007; Lillicrap et al., 2016). Hence, we wanted to initially determine whether our network could take advantage of a hidden layer to reduce error at the output layer.

The network architecture is illustrated in Figure 3B. An image from the MNIST data set is used to set the spike rates of $\mathrm{\ell }=784$ Poisson point-process neurons in the input layer (one neuron per image pixel, rates-of-fire determined by pixel intensity). These project to a hidden layer with $m=500$ neurons. The neurons in the hidden layer (which we index with a ‘0’) are composed of three distinct compartments with their own voltages: the apical compartments (with voltages described by the vector ${\displaystyle {\mathit{V}}^{0a}(t)=[V_1^{0a}(t),...,V_m^{0a}(t)]}$), the basal compartments (with voltages ${\displaystyle {\mathit{V}}^{0b}(t)=[V_1^{0b}(t),...,V_m^{0b}(t)]}$), and the somatic compartments (with voltages ${\displaystyle {\mathit{V}}^0(t)=[V_1^0(t),...,V_m^0(t)]}$). (Note: for notational clarity, all vectors and matrices in the paper are in boldface.) The voltage of the $i^{th}$ neuron in the hidden layer is updated according to:

where $g_l$, $g_b$ and $g_a$ represent the leak conductance, the conductance from the basal dendrites, and the conductance from the apical dendrites, respectively, and $\tau =C_m/g_l$ where $C_m$ is the membrance capacitance (see Materials and methods, Equation (16)). For mathematical simplicity we assume in our simulations a resting membrane potential of 0 mV (this value does not affect the results). We implement electrotonic segregation in the model by altering the $g_a$ value—low values for $g_a$ lead to electrotonically segregated apical dendrites. In the initial set of simulations we set $g_a=0$, which effectively makes it a feed-forward network, but we relax this condition in later simulations.

We treat the voltages in the dendritic compartments simply as weighted sums of the incoming spike trains. Hence, for the $i^{th}$ hidden layer neuron:

where $W_{ij}^0$ and $Y_{ij}$ are synaptic weights from the input layer and the output layer, respectively, $b_i^0$ is a bias term, and ${\displaystyle {\mathit{s}}^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}}$ and ${\displaystyle {\mathit{s}}^1}$ are the filtered spike trains of the input layer and output layer neurons, respectively. (Note: the spike trains are convolved with an exponential kernel to mimic postsynaptic potentials, see Materials and methods Equation (11).)

The somatic compartments generate spikes using Poisson processes. The instantaneous rates of these processes are described by the vector ${\displaystyle {\mathit{\varphi }}^0(t)=[{\varphi }_1^0(t),...,{\varphi }_m^0(t)]}$, which is in units of spikes/s or Hz. These rates-of-fire are determined by a non-linear sigmoid function, $\sigma (\cdot )$, applied to the somatic voltages, that is for the $i^{th}$ hidden layer neuron:

where ${\varphi }_{\max }$ is the maximum rate-of-fire for the neurons.

The output layer (which we index here with a ‘1’) contains $n=10$ two-compartment neurons (one for each image category), similar to those used in a previous model of dendritic prediction learning (Urbanczik and Senn, 2014). The output layer dendritic voltages (${\displaystyle {\mathit{V}}^{1b}(t)=[V_1^{1b}(t),...,V_n^{1b}(t)]}$) and somatic voltages (${\displaystyle {\mathit{V}}^1(t)=[V_1^1(t),...,V_n^1(t)]}$) are updated in a similar manner to the hidden layer basal compartment and soma:

where $W_{ij}^1$ are synaptic weights from the hidden layer, ${\displaystyle {\mathit{s}}^0}$ are the filtered spike trains of the hidden layer neurons (see Equation (11)), $g_l$ is the leak conductance, $g_d$ is the conductance from the dendrites, and $\tau$ is given by Equation (16). In addition to the absence of an apical compartment, the other salient difference between the output layer neurons and the hidden layer neurons is the presence of the term $I_i(t)$, which is a teaching signal that can be used to force the output layer to the correct answer. Whether any such teaching signals exist in the real brain is unknown, though there is evidence that animals can represent desired behavioral outputs with internal goal representations (Gadagkar et al., 2016). (See below, and Materials and methods, Equations (19) and (20) for more details on the teaching signal).

In our model, there are two different types of computation that occur in the hidden layer neurons: ‘transmit’ and ‘plateau’. The transmit computations are standard numerical integration of the simulation, with voltages evolving according to Equation (1), and with the apical compartment electrotonically segregated from the soma (depending on $g_a$) (Figure 3C, left). In contrast, the plateau computations do not involve numerical integration with Equation (1). Instead, the apical voltage is averaged over the most recent 20–30 ms period and the sigmoid non-linearity is applied to it, giving us ‘plateau potentials’ in the hidden layer neurons (we indicate plateau potentials with $\alpha$, see Equation (5) below, and Figure 3C, right). The intention behind this design was to mimic the non-linear transmission from the apical dendrites to the soma that occurs during a plateau potential driven by calcium spikes in the apical dendritic shaft (Larkum et al., 1999), but in the simplest, most abstract formulation possible.

Importantly, plateau potentials in our simulations are single numeric values (one per hidden layer neuron) that can be used for credit assignment. We do not use them to alter the network dynamics. When they occur, they are calculated, transmitted to the basal dendrite instantaneously, and then stored temporarily (0–60 ms) for calculating synaptic weight updates.

To train the network we alternate between two phases. First, during the ‘forward’ phase we present an image to the input layer without any teaching current at the output layer (${\displaystyle I(t{)}_i=0,\mathrm{\forall }i}$). The forward phase occurs between times $t_0$ to $t_1$. At $t_1$ a plateau potential is calculated in all the hidden layer neurons (${\displaystyle {\mathit{\alpha }}^f=[{\alpha }_1^f,...,{\alpha }_m^f]}$) and the ‘target’ phase begins. During this phase, which lasts until $t_2$, the image continues to drive the input layer, but now the output layer also receives teaching current. The teaching current forces the correct output neuron to its max firing rate and all the others to silence. For example, if an image of a ‘9’ is presented, then over the time period $t_1$-$t_2$ the ‘9’ neuron in the output layer fires at max, while the other neurons are silent (Figure 4A). At $t_2$ another set of plateau potentials (${\displaystyle {\mathit{\alpha }}^t=[{\alpha }_1^t,...,{\alpha }_m^t]}$) are calculated in the hidden layer neurons. The result is that we have plateau potentials in the hidden layer neurons for both the end of the forward phase (${\displaystyle {\mathit{\alpha }}^f}$) and the end of the target phase (${\displaystyle {\mathit{\alpha }}^t}$), which are calculated as:

Figure 4

Download asset Open asset

where $\mathrm{\Delta }{t}_s$ is a time delay used to allow the network dynamics to settle before integrating the plateau, and $\mathrm{\Delta }{t}_i=t_i-(t_{i-1}+\mathrm{\Delta }{t}_s)$ (see Materials and methods, Equation (22) and Figure 4A).

Similar to how targets are used in deep supervised learning (LeCun et al., 2015), the goal of learning in our network is to make the network dynamics during the forward phase converge to the same output activity pattern as exists in the target phase. Put another way, in the absence of the teaching signal, we want the activity at the output layer to be the same as that which would exist with the teaching signal, so that the network can give appropriate outputs without any guidance. To do this, we initialize all the weight matrices with random weights, then we train the weight matrices ${\displaystyle {\mathit{W}}^0}$ and ${\displaystyle {\mathit{W}}^1}$ using stochastic gradient descent on local loss functions for the hidden and output layers, respectively (see below). These weight updates occur at the end of every target phase, that is the synapses are not updated during transmission. Like Lillicrap et al. (2016), we leave the weight matrix ${\displaystyle \mathit{Y}}$ fixed in its initial random configuration. When we update the synapses in the network we use the plateau potential values ${\displaystyle {\mathit{\alpha }}^f}$ and ${\displaystyle {\mathit{\alpha }}^t}$ to determine appropriate credit assignment (see below).

The network is simulated in near continuous-time (except that each plateau is considered to be instantaneous), and the temporal intervals between plateaus are randomly sampled from an inverse Gaussian distribution (Figure 4B, top). As such, the specific amount of time that the network is presented with each image and teaching signal is stochastic, though usually somewhere between 50–60 ms of simulated time (Figure 4B, bottom). This stochasticity was not necessary, but it demonstrates that although the system operates in phases, the specific length of the phases is not important as long as they are sufficiently long to permit integration (see Lemma 1). In the data presented in this paper, all 60,000 images in the MNIST training set were presented to the network one at a time, and each exposure to the full set of images was considered an ‘epoch’ of training. At the end of each epoch, the network’s classification error rate on a separate set of 10,000 test images was assessed with a single forward phase for each image (see Materials and methods). The network’s classification was judged by which output neuron had the highest average firing rate during these test image forward phases.

It is important to note that there are many aspects of this design that are not physiologically accurate. Most notably, stochastic generation of plateau potentials across a population is not an accurate reflection of how real pyramidal neurons operate, since apical calcium spikes are determined by a number of concrete physiological factors in individual cells, including back-propagating action potentials, spike-timing and inhibitory inputs (Larkum et al., 1999, 2007, 2009). However, we note that calcium spikes in the apical dendrites can be prevented from occurring via the activity of distal dendrite targeting inhibitory interneurons (Murayama et al., 2009), which can synchronize pyramidal activity (Hilscher et al., 2017). Furthermore, distal dendrite targeting interneurons can themselves can be rapidly inhibited in response to temporally precise neuromodulatory inputs (Pi et al., 2013; Pfeffer et al., 2013; Karnani et al., 2016; Hangya et al., 2015; Brombas et al., 2014). Therefore, it is entirely plausible that neocortical micro-circuits would generate synchronized plateaus/bursts at punctuated periods of time in response to disinhibition of the apical dendrites governed by neuromodulatory signals that determine ‘phases’ of processing. Alternatively, oscillations in population activity could provide a mechanism for promoting alternating phases of processing and synaptic plasticity (Buzsáki and Draguhn, 2004). But, complete synchrony of plateaus in our hidden layer neurons is not actually critical to our algorithm—only the temporal relationship between the plateaus and the teaching signal is critical. This relationship itself is arguably plausible given the role of neuromodulatory inputs in dis-inhibiting the distal dendrites of pyramidal neurons (Karnani et al., 2016; Brombas et al., 2014). Of course, we are engaged in a great deal of speculation here. But, the point is that our model utilizes anatomical and functional motifs that are loosely analogous to what is observed in the neocortex. Importantly for the present study, the key issue is the use of segregated dendrites which permit an effective feed-forward dynamic, punctuated by feedback driven plateau potentials to solve the credit assignment problem.

To solve the credit assignment problem without using weight transport, we had to define local error signals, or ‘loss functions’, for the hidden layer and output layer that somehow took into account the impact that each hidden layer neuron has on the output of the network. In other words, we only want to update a hidden layer synapse in a manner that will help us make the forward phase activity at the output layer more similar to the target phase activity. To begin, we define the target firing rates for the output neurons, ${\displaystyle {\mathit{\varphi }}^{1\ast }=[{\varphi }_1^{1\ast },...,{\varphi }_n^{1\ast }]}$, to be their average firing rates during the target phase:

(Throughout the paper, we use ${\varphi }^{*}$ to denote a target firing rate and $\overline{\varphi }$ to denote a firing rate averaged over time.) We then define a loss function at the output layer using this target, by taking the difference between the average forward phase activity and the target:

(Note: the true loss function we use is slightly more complex than the one formulated here, hence the $\approx$ symbol in Equation (7), but this formulation is roughly correct and easier to interpret. See Materials and methods, Equation (23) for the exact formulation.) This loss function is zero only when the average firing rates of the output neurons during the forward phase equals their target, that is the average firing rates during the target phase. Thus, the closer $L^1$ is to zero, the more the network’s output for an image matches the output activity pattern imposed by the teaching signal, $I(t)$.

Effective credit assignment is achieved when changing the hidden layer synapses is guaranteed to reduce $L^1$. To obtain this guarantee, we defined a set of target firing rates for the hidden layer neurons that uses the information contained in the plateau potentials. Specifically, in a similar manner to Lee et al., 2015, we define the target firing rates for the hidden layer neurons, ${\displaystyle {\mathit{\varphi }}^{0\ast }=[{\varphi }_1^{0\ast },...,{\varphi }_m^{0\ast }]}$, to be:

where ${\alpha }_i^t$ and ${\alpha }_i^f$ are the plateaus defined in Equation (5). As with the output layer, we define the loss function for the hidden layer to be the difference between the target firing rate and the average firing rate during the forward phase:

(Again, note the use of the $\approx$ symbol, see Equation (30) for the exact formulation.) This loss function is zero only when the plateau at the end of the forward phase equals the plateau at the end of the target phase. Since the plateau potentials integrate the top-down feedback (see Equation (5)), we know that the hidden layer loss function, $L^0$, is zero if the output layer loss function, $L^1$, is zero. Moreover, we can show that these loss functions provide a broader guarantee that, under certain conditions, if $L^0$ is reduced, then on average, $L^1$ will also be reduced (see Theorem 1). This provides our assurance of credit assignment: we know that the ultimate goal of learning (reducing $L^1$) can be achieved by updating the synaptic weights at the hidden layer to reduce the local loss function ${\displaystyle L^0}$ (Figure 5A). We do this using stochastic gradient descent at the end of every target phase:

Figure 5 with 1 supplement see all

Download asset Open asset

Figure 5—source data 1

Fig_5B.csv.

The first two columns of the data file contain the hidden layer loss ($L^0$) and output layer loss ($L^1$) of a one hidden layer network in response to all ‘2’ images in the MNIST test set after one epoch of training. The last two columns contain the same data, except that the data in the third column (Shuffled data $L^0$) was generated by randomly shuffling the hidden layer activity vectors. Fig_5C.csv. The first 10 columns of the data file contain the mean Pearson correlation coefficient between the hidden layer loss ($L^0$) and output layer loss ($L^1$) of the one hidden layer network in response to each category of handwritten digits across training. Each row represents one epoch of training. The last 10 columns contain the mean Pearson correlation coefficients between the shuffled hidden layer loss and the output layer loss for each category, across training. Fig_5S1A.csv. This data file contains the maximum eigenvalue of ${(I-\overline{J_{\beta }}\overline{J_{\gamma }})}^T(I-\overline{J_{\beta }}\overline{J_{\gamma }})$ over 60,000 training examples for a one hidden layer network, where $\overline{J_{\beta }}$ and $\overline{J_{\gamma }}$ are the mean feedforward and feedback Jacobian matrices for the last 100 training examples.

https://doi.org/10.7554/eLife.22901.009

Download elife-22901-fig5-data1-v1.zip

 where ${\eta }_i$ and ${\displaystyle \mathrm{\Delta }{\mathit{W}}^i}$ refer to the learning rate and update term for weight matrix ${\displaystyle {\mathit{W}}^i}$ (see Materials and methods, Equations (28), (29), (33) and (35) for details of the weight update procedures). Performing gradient descent on $L^1$ results in a relatively straight-forward delta rule update for ${\displaystyle {\mathit{W}}^1}$ (see Equation (29)). The weight update for the hidden layer weights, ${\displaystyle {\mathit{W}}^0}$, is similar, except for the presence of the difference between the two plateau potentials ${\displaystyle {\mathit{\alpha }}^t-{\mathit{\alpha }}^f}$ (see Equation (35)). Importantly, given the way in which we defined the loss functions, as the hidden layer reduces ${\displaystyle L^0}$ by updating ${\displaystyle {\mathit{W}}^0}$, $L^1$ should also be reduced, that is hidden layer learning should imply output layer learning, thereby utilizing the multi-layer architecture.

To test that we were successful in credit assignment with this design, and to provide empirical support for the proof of Theorem 1, we compared the loss function at the hidden layer, $L^0$, to the output layer loss function, $L^1$, across all of the image presentations to the network. We observed that, generally, whenever the hidden layer loss was low, the output layer loss was also low. For example, when we consider the loss for the set of ‘2’ images presented to the network during the second epoch, there was a Pearson correlation coefficient between $L^0$ and $L^1$ of $r=0.61$, which was much higher than what was observed for shuffled data, wherein output and hidden activities were randomly paired (Figure 5B). Furthermore, these correlations were observed across all epochs of training, with most correlation coefficients for the hidden and output loss functions falling between $r=0.2$ - $0.6$, which was, again, much higher than the correlations observed for shuffled data (Figure 5C).

Interestingly, the correlations between $L^0$ and $L^1$ were smaller on the first epoch of training (see data in red oval Figure 5C) . This suggests that the guarantee of coordination between $L^0$ and $L^1$ only comes into full effect once the network has engaged in some learning. Therefore, we inspected whether the conditions on the synaptic matrices that are assumed in the proof of Theorem 1 were, in fact, being met. More precisely, the proof assumes that the feedforward and feedback synaptic matrices (${\displaystyle {\mathit{W}}^1}$ and ${\displaystyle \mathit{Y}}$, respectively) produce forward and backward transformations between the output and hidden layer whose Jacobians are approximate inverses of each other (see Proof of Theorem 1). Since we begin learning with random matrices, this condition is almost definitely not met at the start of training. But, we found that the network learned to meet this condition. Inspection of ${\displaystyle {\mathit{W}}^1}$ and ${\displaystyle \mathit{Y}}$ showed that during the first epoch, the Jacobians of the forward and backwards functions became approximate inverses of each other (Figure 5—figure supplement 1). Since ${\displaystyle \mathit{Y}}$ is frozen, this means that during the first few image presentations ${\displaystyle {\mathit{W}}^1}$ was being updated to have its Jacobian come closer to the inverse of ${\displaystyle \mathit{Y}}$'s Jacobian. Put another way, the network was learning to do credit assignment. We have yet to resolve exactly why this happens, though the result is very similar to the findings of Lillicrap et al. (2016), where a proof is provided for the linear case. Intuitively, though, the reason is likely the interaction between ${\displaystyle {\mathit{W}}^1}$ and ${\displaystyle {\mathit{W}}^0}$: as ${\displaystyle {\mathit{W}}^0}$ gets updated, the hidden layer learns to group stimuli based on the feedback sent through ${\displaystyle \mathit{Y}}$. So, for ${\displaystyle {\mathit{W}}^1}$ to transform the hidden layer activity into the correct output layer activity, ${\displaystyle {\mathit{W}}^1}$ must become more like the inverse of ${\displaystyle \mathit{Y}}$, which would also make the Jacobian of ${\displaystyle {\mathit{W}}^1}$ more like the inverse of ${\displaystyle \mathit{Y}}$’s Jacobian (due to the inverse function theorem). However, a complete, formal explanation for this phenomenon is still missing, and the the issue of weight alignment deserves additional investigation Lillicrap et al. (2016). From a biological perspective, it also suggests that very early development may involve a period of learning how to assign credit appropriately. Altogether, our model demonstrates that deep learning using random feedback weights is a general phenomenon, and one which can be implemented using segregated dendrites to keep forward information separate from feedback signals used for credit assignment.

Given our finding that the network was successfully assigning credit for the output error to the hidden layer neurons, we had reason to believe that our network with local weight-updates would exhibit deep learning, that is an ability to take advantage of a multi-layer structure (Bengio and LeCun, 2007). To test this, we examined the effects of including hidden layers. If deep learning is indeed operational in the network, then the inclusion of hidden layers should improve the ability of the network to classify images.

We built three different versions of the network (Figure 6A). The first was a network that had no hidden layer, that is the input neurons projected directly to the output neurons. The second was the network illustrated in Figure 3B, with a single hidden layer. The third contained two hidden layers, with the output layer projecting directly back to both hidden layers. This direct projection allowed us to build our local targets for each hidden layer using the plateaus driven by the output layer, thereby avoiding a ‘backward pass’ through the entire network as has been used in other models (Lillicrap et al., 2016; Lee et al., 2015; Liao et al., 2015). We trained each network on the 60,000 MNIST training images for 60 epochs, and recorded the percentage of images in the 10,000 image test set that were incorrectly classified. The network with no hidden layers rapidly learned to classify the images, but it also rapidly hit an asymptote at an average error rate of 8.3% (Figure 6B, gray line). In contrast, the network with one hidden layer did not exhibit a rapid convergence to an asymptote in its error rate. Instead, it continued to improve throughout all 60 epochs, achieving an average error rate of 4.1% by the 60^th epoch (Figure 6B, blue line). Similar results were obtained when we loosened the synchrony constraints and instead allowed each hidden layer neuron to engage in plateau potentials at different times (Figure 6—figure supplement 1). This demonstrates that strict synchrony in the plateau potentials is not required. But, our target definitions do require two different plateau potentials separated by the teaching signal input, which mandates some temporal control of plateau potentials in the system.

Figure 6 with 1 supplement see all

Download asset Open asset

Figure 6—source data 1

Fig_6B_errors.csv.

This data file contains the test error (measured on 10,000 MNIST images not used for training) across 60 epochs of training, for a network with no hidden layers, a network with one hidden layer, and a network with two hidden layers. Fig_6B_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 60 epochs for each of the three networks described above. Fig_6C.csv. The first column of this data file contains the categories of 10,000 MNIST images presented to a two hidden layer network (after 60 epochs of training). The next three pairs of columns contain the $x$ and $y$-coordinates of the t-SNE two-dimensional reduction of the activity patterns of the input layer, the first hidden layer, and the second hidden layer, respectively. Fig_6S1B_errors.csv. This data file contains the test error (measured on 10,000 MNIST images not used for training) across 60 epochs of training, for a one hidden layer network, with synchronized plateau potentials (Regular) and with stochastic plateau potentials. Fig_6S1B_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 60 epochs for each of the two networks described above.

https://doi.org/10.7554/eLife.22901.012

Download elife-22901-fig6-data1-v1.zip

Interestingly, we found that the addition of a second hidden layer further improved learning. The network with two hidden layers learned more rapidly than the network with one hidden layer and achieved an average error rate of 3.2% on the test images by the 60^th epoch, also without hitting a clear asymptote in learning (Figure 6B, red line). However, it should be noted that additional hidden layers beyond two did not significantly improve the error rate (data not shown), which suggests that our particular algorithm could not be used to construct very deep networks as is. Nonetheless, our network was clearly able to take advantage of multi-layer architectures to improve its learning, which is the key feature of deep learning (Bengio and LeCun, 2007; LeCun et al., 2015).

Another key feature of deep learning is the ability to generate representations in the higher layers of a network that capture task-relevant information while discarding sensory details (LeCun et al., 2015; Mnih et al., 2015). To examine whether our network exhibited this type of abstraction, we used the t-Distributed Stochastic Neighbor Embedding algorithm (t-SNE). The t-SNE algorithm reduces the dimensionality of data while preserving local structure and non-linear manifolds that exist in high-dimensional space, thereby allowing accurate visualization of the structure of high-dimensional data (Maaten and Hinton, 2008). We applied t-SNE to the activity patterns at each layer of the two hidden layer network for all of the images in the test set after 60 epochs of training. At the input level, there was already some clustering of images based on their categories. However, the clusters were quite messy, with different categories showing outliers, several clusters, or merged clusters (Figure 6C, bottom). For example, the ‘2’ digits in the input layer exhibited two distinct clusters separated by a cluster of ‘7’s: one cluster contained ‘2’s with a loop and one contained ‘2’s without a loop. Similarly, there were two distinct clusters of ‘4’s and ‘9’s that were very close to each other, with one pair for digits on a pronounced slant and one for straight digits (Figure 6C, bottom, example images). Thus, although there is built-in structure to the categories of the MNIST dataset, there are a number of low-level features that do not respect category boundaries. In contrast, at the first hidden layer, the activity patterns were much cleaner, with far fewer outliers and split/merged clusters (Figure 6C, middle). For example, the two separate ‘2’ digit clusters were much closer to each other and were now only separated by a very small cluster of ‘7’s. Likewise, the ‘9’ and ‘4’ clusters were now distinct and no longer split based on the slant of the digit. Interestingly, when we examined the activity patterns at the second hidden layer, the categories were even better segregated, with only a little bit of splitting or merging of category clusters (Figure 6C, top). Therefore, the network had learned to develop representations in the hidden layers wherein the categories were very distinct and low-level features unrelated to the categories were largely ignored. This abstract representation is likely to be key to the improved error rate in the two hidden layer network. Altogether, our data demonstrates that our network with segregated dendritic compartments can engage in deep learning.

The backpropagation of error algorithm (Rumelhart et al., 1986) is still the primary learning algorithm used for deep supervised learning in artificial neural networks (LeCun et al., 2015). Previous work has shown that learning with random feedback weights can actually match the synaptic weight updates specified by the backpropagation algorithm after a few epochs of training (Lillicrap et al., 2016). This fascinating observation suggests that deep learning with random feedback weights is not completely distinct from backpropagation of error, but rather, networks with random feedback connections learn to approximate credit assignment as it is done in backpropagation (Lillicrap et al., 2016). Hence, we were curious as to whether or not our network was, in fact, learning to approximate the synaptic weight updates prescribed by backpropagation. To test this, we trained our one hidden layer network as before, but now, in addition to calculating the vector of hidden layer synaptic weight updates specified by our local learning rule (${\displaystyle \mathrm{\Delta }{\mathit{W}}^0}$ in Equation (10)), we also calculated the vector of hidden layer synaptic weight updates that would be specified by non-locally backpropagating the error from the output layer, (${\displaystyle \mathrm{\Delta }{\mathit{W}}_{BP}^0}$). We then calculated the angle between these two alternative weight updates. In a very high-dimensional space, any two independent vectors will be roughly orthogonal to each other (i.e. ${\displaystyle \mathrm{\Delta }{\mathit{W}}^0\mathrm{\angle }\mathrm{\Delta }{\mathit{W}}_{BP}^0\approx {90}^{\circ }}$). If the two synaptic weight update vectors are not orthogonal to each other (i.e. ${\displaystyle \mathrm{\Delta }{\mathit{W}}^0\mathrm{\angle }\mathrm{\Delta }{\mathit{W}}_{BP}^0<{90}^{\circ }}$), then it suggests that the two algorithms are specifying similar weight updates.

As in previous work (Lillicrap et al., 2016), we found that the initial weight updates for our network were orthogonal to the updates specified by backpropagation. But, as the network learned the angle dropped to approximately ${65}^{\circ }$, before rising again slightly to roughly ${70}^{\circ }$ (Figure 7A, blue line). This suggests that our network was learning to develop local weight updates in the hidden layer that were in rough agreement with the updates that explicit backpropagation would produce. However, this drop in orthogonality was still much less than that observed in non-spiking artificial neural networks learning with random feedback weights, which show a drop to below ${45}^{\circ }$(Lillicrap et al., 2016). We suspected that the higher angle between the weight updates that we observed may have been because we were using spikes to communicate the feedback from the upper layer, which could introduce both noise and bias in the estimates of the output layer activity. To test this, we also examined the weight updates that our algorithm would produce if we propagated the spike rates of the output layer neurons, ${\displaystyle {\mathit{\varphi }}^1(t)}$, back directly through the random feedback weights, ${\displaystyle \mathit{Y}}$. In this scenario, we observed a much sharper drop in the ${\displaystyle \mathrm{\Delta }{\mathit{W}}^0\mathrm{\angle }\mathrm{\Delta }{\mathit{W}}_{BP}^0}$ angle, which reduced to roughly ${35}^{\circ }$ before rising again to ${40}^{\circ }$ (Figure 7A, red line). These results show that, in principle, our algorithm is learning to approximate the backpropagation algorithm, though with some drop in accuracy introduced by the use of spikes to propagate output layer activities to the hidden layer.

Figure 7

Download asset Open asset

Figure 7—source data 1

Fig_7A.csv.

This data file contains the time-averaged angle (with a sliding window of 100 images) between weight updates prescribed by our local update learning algorithm compared to those prescribed by backpropagation of error, for a one hidden layer network over 10 epochs of training (600,000 training examples). Fig_7C.csv. The first column of this data file contains the maximum Pearson correlation coefficient between each receptive field learned using our algorithm and all 500 receptive fields learned using backpropagation. The second column of this data file contains the maximum Pearson correlation coefficient between a randomly shuffled version of each receptive field learned using our algorithm and all 500 receptive fields learned using backpropagation.

https://doi.org/10.7554/eLife.22901.014

Download elife-22901-fig7-data1-v1.zip

To further examine how our local learning algorithm compared to backpropagation we compared the low-level features that the two algorithms learned. To do this, we trained the one hidden layer network with both our algorithm and backpropagation. We then examined the receptive fields (i.e. the synaptic weights) produced by both algorithms in the hidden layer synapses (${\displaystyle {\mathit{W}}^0}$) after 60 epochs of training. The two algorithms produced qualitatively similar receptive fields (Figure 7B). Both produced receptive fields with clear, high-contrast features for detecting particular strokes or shapes. To quantify the similarity, we conducted pair-wise correlation calculations for the receptive fields produced by the two algorithms and identified the maximum correlation pairs for each. Compared to shuffled versions of the receptive fields, there was a very high level of maximum correlation (Figure 7C), showing that the receptive fields were indeed quite similar. Thus, the data demonstrate that our learning algorithm using random feedback weights into segregated dendrites can in fact come to approximate the backpropagation of error algorithm.

Once we had convinced ourselves that our learning algorithm was, in fact, providing a solution to the credit assignment problem, we wanted to examine some of the constraints on learning. First, we wanted to explore the structure of the feedback weights. In our initial simulations we used non-sparse, random (i.e. normally distributed) feedback weights. We were interested in whether learning could still work with sparse weights, given that neocortical connectivity is sparse. As well, we wondered whether symmetric weights would improve learning, which would be expected given previous findings (Lillicrap et al., 2016; Lee et al., 2015; Liao et al., 2015). To explore these questions, we trained our one hidden layer network using both sparse feedback weights (only 20% non-zero values) and symmetric weights (${\displaystyle \mathit{Y}={{\mathit{W}}^1}^T}$) (Figure 8A,C). We found that learning actually improved slightly with sparse weights (Figure 8B, red line), achieving an average error rate of 3.7% by the 60^th epoch, compared to the average 4.1% error rate achieved with fully random weights. But, this result appeared to depend on the magnitude of the sparse weights. To compensate for the loss of 80% of the weights we initially increased the sparse synaptic weight magnitudes by a factor of 5. However, when we did not re-scale the sparse weights learning was actually worse (Figure 8—figure supplement 1), though this could likely be dealt with by a careful resetting of learning rates. Altogether, our results suggest that sparse feedback provides a signal that is sufficient for credit assignment.

Figure 8 with 1 supplement see all

Download asset Open asset

Figure 8—source data 1

Fig_8B_errors.csv.

This data file contains the test error (measured on 10,000 MNIST images not used for training) across 60 epochs of training, for our standard one hidden layer network (Regular) and a network with sparse feedback weights. Fig_8B_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 60 epochs for each of the two networks described above. Fig_8D_errors.csv. This data file contains the test error (measured on 10,000 MNIST images not used for training) across 60 epochs of training, for our standard one hidden layer network (Regular), a network with symmetric weights, and a network with symmetric weights with added noise. Fig_8D_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 60 epochs for each of the three networks described above. Fig_8S1_errors.csv. This data file contains the test error (measured on 10,000 MNIST images not used for training) across 20 epochs of training, for a one hidden layer network with regular feedback weights, sparse feedback weights that were amplified, and sparse feedback weights that were not amplified. Fig_8S1_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 20 epochs for each of the three networks described above.

https://doi.org/10.7554/eLife.22901.017

Download elife-22901-fig8-data1-v1.zip

Similar to sparse feedback weights, symmetric feedback weights also improved learning, leading to a rapid decrease in the test error and an error rate of 3.6% by the 60^th epoch (Figure 8D, red line). This is interesting, given that backpropagation assumes symmetric feedback weights (Lillicrap et al., 2016; Bengio et al., 2015), though our proof of Theorem 1 does not. However, when we added noise to the symmetric weights any advantage was eliminated and learning was, in fact, slightly impaired (Figure 8D, blue line). At first, this was a very surprising result: given that learning works with random feedback weights, why would it not work with symmetric weights with noise? However, when we considered our previous finding that during the first epoch the feedforward weights, ${\displaystyle {\mathit{W}}^1}$, learn to have the feedforward Jacobian match the inverse of the feedback Jacobian (Figure 5—figure supplement 1) a possible answer emerges. In the case of symmetric feedback weights the synaptic matrix ${\displaystyle \mathit{Y}}$ is changing as ${\displaystyle {\mathit{W}}^1}$ changes. This works fine when ${\displaystyle \mathit{Y}}$ is set to ${\displaystyle {\mathit{W}}^{1^T}}$, since that artificially forces something akin to backpropagation. But, if the feedback weights are set to ${\displaystyle {\mathit{W}}^{1^T}}$ plus noise, then the system can never align the Jacobians appropriately, since ${\displaystyle \mathit{Y}}$ is now a moving target. This would imply that any implementation of feedback learning must either be very effective (to achieve the right feedback) or very slow (to allow the feedforward weights to adapt).

Another constraint that we wished to examine was whether total segregation of the apical inputs was necessary, given that real pyramidal neurons only show an attenuation of distal apical inputs to the soma (Larkum et al., 1999). Total segregation ($g_a=0$) renders the network effectively feed-forward in its dynamics, which made it easier to construct the loss functions to ensure that reducing $L^0$ also reduces $L^1$ (see Figure 5 and Theorem 1). But, we wondered whether some degree of apical conductance to the soma would be sufficiently innocuous so as to not disrupt deep learning. To examine this, we re-ran our two hidden layer network, but now, we allowed the apical dendritic voltage to influence the somatic voltage by setting $g_a=0.05$. This value gave us twelve times more attenuation than the attenuation from the basal compartments, since $g_b=0.6$ (Figure 9A). When we compared the learning in this scenario to the scenario with total apical segregation, we observed very little difference in the error rates on the test set (Figure 9B, gray and red lines). Importantly, though, we found that if we increased the apical conductance to the same level as the basal ($g_a=g_b=0.6$) then the learning was significantly impaired (Figure 9B, blue line). This demonstrates that although total apical attenuation is not necessary, partial segregation of the apical compartment from the soma is necessary. That result makes sense given that our local targets for the hidden layer neurons incorporate a term that is supposed to reflect the response of the output neurons to the feedforward sensory information (${\displaystyle {\mathit{\alpha }}^f}$). Without some sort of separation of feedforward and feedback information, as is assumed in other models of deep learning (Lillicrap et al., 2016; Lee et al., 2015), this feedback signal would get corrupted by recurrent dynamics in the network. Our data show that electrontonically segregated dendrites is one potential way to achieve the separation between feedforward and feedback information that is required for deep learning.

Figure 9

Download asset Open asset

Figure 9—source data 1

Fig_9B_errors.csv.

This data file contains the test error (measured on 10,000 MNIST images not used for training) across 60 epochs of training, for a two hidden layer network, with total apical segregation (Regular), strong apical attenuation and weak apical attenuation. Fig_9B_final_errors.csv. This data file contains the results of repeated weight tests ($n=20$) after 60 epochs for each of the three networks described above.

https://doi.org/10.7554/eLife.22901.019

Download elife-22901-fig9-data1-v1.zip

The network described here consists of an input layer with $\mathrm{\ell }$ neurons, a hidden layer with $m$ neurons, and an output layer with $n$ neurons. Neurons in the input layer are simple Poisson spiking neurons whose rate-of-fire is determined by the intensity of image pixels (ranging from 0 - ${\displaystyle {\varphi }_{\mathrm{m}\mathrm{a}\mathrm{x}}}$). Neurons in the hidden layer are modeled using three functional compartments—basal dendrites with voltages ${\displaystyle {\mathit{V}}^{0b}(t)=[V_1^{0b}(t),V_2^{0b}(t),...,V_m^{0b}(t)]}$, apical dendrites with voltages ${\displaystyle {\mathit{V}}^{0a}(t)=[V_1^{0a}(t),V_2^{0a}(t),...,V_m^{0a}(t)]}$, and somata with voltages ${\displaystyle {\mathit{V}}^0(t)=[V_1^0(t),V_2^0(t),...,V_m^0(t)]}$. Feedforward inputs from the input layer and feedback inputs from the output layer arrive at basal and apical synapses, respectively. At basal synapses, presynaptic spikes from input layer neurons are translated into filtered spike trains ${\displaystyle {\mathit{s}}^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}(t)=[s_1^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}(t),s_2^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}(t),...,s_{\ell }^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}(t)]}$ given by:

where ${\displaystyle t_{jk}^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}}$ is the $k$ th spike time of input neuron ${\displaystyle j}$ is the response kernel given by:

where ${\tau }_s$ and ${\tau }_L$ are short and long time constants, and $\mathrm{\Theta }$ is the Heaviside step function. Since the network is fully-connected, each neuron in the hidden layer will receive the same set of filtered spike trains from input layer neurons. The filtered spike trains at apical synapses, ${\displaystyle {\mathit{s}}^1(t)=[s_1^1(t),s_2^1(t),...,s_n^1(t)]}$, are modeled in the same manner. The basal and apical dendritic potentials for neuron $i$ are then given by weighted sums of the filtered spike trains at either its basal or apical synapses:

where ${\displaystyle {\mathit{b}}^0=[b_1^0,b_2^0,...,b_m^0]}$ are bias terms, ${\displaystyle {\mathit{W}}^0}$ is the $m\times \mathrm{\ell }$ matrix of feedforward weights for neurons in the hidden layer, and ${\displaystyle \mathit{Y}}$ is the $m\times n$ matrix of their feedback weights. The somatic voltage for neuron $i$ evolves with leak as:

where $V^R$ is the resting potential, $g_l$ is the leak conductance, $g_b$ is the conductance from the basal dendrite to the soma, and $g_a$ is the conductance from the apical dendrite to the soma, and $\tau$ is a function of $g_l$ and the membrance capacitance $C_m$:

Note that for simplicity’s sake we are assuming a resting potential of 0 mV and a membrane capacitance of 1 F, but these values are not important for the results. Equations (13) and (14) are identical to the Equation (1) in results.

The instantaneous firing rates of neurons in the hidden layer are given by ${\displaystyle {\mathit{\varphi }}^0(t)=[{\varphi }_1^0(t),{\varphi }_2^0(t),...,{\varphi }_m^0(t)]}$, where ${\varphi }_i^0(t)$ is the result of applying a nonlinearity, $\sigma (\cdot )$, to the somatic potential $V_i^0(t)$. We chose $\sigma (\cdot )$ to be a simple sigmoidal function, such that:

Here, ${\varphi }_{\max }$ is the maximum possible rate-of-fire for the neurons, which we set to 200 Hz. Note that Equation (17) is identical to Equation (3) in results. Spikes are then generated using Poisson processes with these firing rates. We note that although the maximum rate was 200 Hz, the neurons rarely achieved anything close to this rate, and the average rate of fire in the neurons during our simulations was 24 Hz.

Units in the output layer are modeled using only two compartments, dendrites with voltages ${\displaystyle {\mathit{V}}^{1b}(t)=[V_1^{1b}(t),V_2^{1b}(t),...,V_n^{1b}(t)]}$ and somata with voltages ${\displaystyle {\mathit{V}}^1(t)=[V_1^1(t),V_2^1(t),...,V_n^1(t)]}$ is given by:

where ${\displaystyle {\mathit{s}}^0(t)=[s_1^0(t),s_2^0(t),...,s_m^0(t)]}$ are the filtered presynaptic spike trains at synapses that receive feedforward input from the hidden layer, and are calculated in the manner described by Equation (11). $V_i^1(t)$ evolves as:

where $g_l$ is the leak conductance, $g_d$ is the conductance from the dendrite to the soma, and ${\displaystyle \mathit{I}(t)=[I_1(t),I_2(t),...,I_n(t)]}$ are somatic currents that can drive output neurons toward a desired somatic voltage. For neuron $i$, $I_i$ is given by:

where ${\displaystyle {\mathit{g}}_{\mathit{E}}(t)=[g_{E_1}(t),g_{E_2}(t),...,g_{E_n}(t)]}$ and ${\displaystyle {\mathit{g}}_{\mathit{I}}(t)=[g_{I_1}(t),g_{I_2}(t),...,g_{I_n}(t)]}$ are time-varying excitatory and inhibitory nudging conductances, and $E_E$ and $E_I$ are the excitatory and inhibitory reversal potentials. In our simulations, we set $E_E=8$ V and $E_I=-8$ V. During the target phase only, we set ${\displaystyle g_{I_i}=1}$ and $g_{E_i}=0$ for all units $i$ whose output should be minimal, and $g_{E_i}=1$ and $g_{I_i}=0$ for the unit whose output should be maximal. In this way, all units other than the ‘target’ unit are silenced, while the ‘target’ unit receives a strong excitatory drive. In the forward phase, ${\displaystyle \mathit{I}(t)}$ is set to 0. The Poisson spike rates ${\displaystyle {\mathit{\varphi }}^1(t)=[{\varphi }_1^1(t),{\varphi }_2^1(t),...,{\varphi }_n^1(t)]}$ are calculated as in Equation (17).

At the end of the forward and target phases, we calculate plateau potentials ${\displaystyle {\mathit{\alpha }}^f=[{\alpha }_1^f,{\alpha }_2^f,...,{\alpha }_m^f]}$ and ${\displaystyle {\mathit{\alpha }}^t=[{\alpha }_1^t,{\alpha }_2^t,...,{\alpha }_m^t]}$ for apical dendrites of hidden layer neurons, where ${\alpha }_i^f$ and ${\alpha }_i^t$ are given by:

where $t_1$ and $t_2$ are the end times of the forward and target phases, respectively, $\mathrm{\Delta }{t}_s=30$ ms is the settling time for the voltages, and ${\displaystyle \mathrm{\Delta }{t}_1}$ and ${\displaystyle \mathrm{\Delta }{t}_2}$ are given by:

Note that Equation (21) is identical to Equation (5) in results. These plateau potentials are used by hidden layer neurons to update their basal weights.

All feedforward synaptic weights are updated at the end of each target phase. Output layer units update their synaptic weights ${\displaystyle {\mathit{W}}^1}$ in order to minimize the loss function

where ${\displaystyle {\mathit{\varphi }}^{1\ast }={\overline{{\mathit{\varphi }}^1}}^t}$ as in Equation (6). Note that, as long as neuronal units calculate averages after the network has reached a steady state, and the firing-rates of the neurons are in the linear region of the sigmoid function, then for layer $x$,

Thus,

as in Equation (7).

All average voltages are calculated after a delay $\mathrm{\Delta }{t}_s$ from the start of a phase, which allows for the network to reach a steady state before averaging begins. In practice this means that the average somatic voltage for output layer neuron $i$ in the forward phase, ${\displaystyle {\overline{V_i^1}}^f}$, has the property

where $k_d$ is given by:

Thus,

Note that these partial derivatives assume that the activity during the target phase is fixed. We do this because the goal of learning is to have the network behave as it does during the target phase, even when the teaching signal is present. Thus, we do not update synapses in order to alter the target phase activity. As a result, there are no terms in the equation related to the partial derivatives of the voltages or firing-rates during the target phase.

The dendrites in the output layer use this approximation of the gradient in order to update their weights using gradient descent:

where ${\eta }^1$ is a learning rate constant, and $P^1$ is a scaling factor used to normalize the scale of the rate-of-fire function.

In the hidden layer, basal dendrites update their synaptic weights ${\displaystyle {\mathit{W}}^0}$ by minimizing the loss function

We define the target rates-of-fire ${\displaystyle {\mathit{\varphi }}^{0\ast }=[{\varphi }_1^{0\ast },{\varphi }_2^{0\ast },...,{\varphi }_m^{0\ast }]}$ such that

where ${\displaystyle {\mathit{\alpha }}^f=[{\alpha }_1^f,{\alpha }_2^f,...,{\alpha }_m^f]}$ and ${\displaystyle {\mathit{\alpha }}^t=[{\alpha }_1^t,{\alpha }_2^t,...,{\alpha }_m^t]}$ are forward and target phase plateau potentials given in Equation (21). Note that Equation (31) is identical to Equation (8) in results. These hidden layer target firing rates are similar to the targets used in difference target propagation (Lee et al., 2015).

Using Equation (24), we can show that

as in Equation (9). Hence:

where $k_b$ is given by:

Note that although ${\displaystyle {\mathit{\varphi }}^{0\ast }}$ is a function of ${\displaystyle {\mathit{W}}^0}$ and ${\displaystyle {\mathit{b}}^0}$, we do not differentiate this term with respect to the weights and biases. Instead, we treat ${\displaystyle {\mathit{\varphi }}^{0\ast }}$ as a fixed state for the hidden layer neurons to learn to reproduce. Basal weights are updated in order to descend this approximation of the gradient:

Again, we assume that the activity during the target phase is fixed, so no derivatives are taken with respect to voltages or firing-rates during the target phase.

Importantly, this update rule is spatially local for the hidden layer neurons. It consists essentially of three terms, (1) the difference in the plateau potentials for the target and forward phases (${\displaystyle {\mathit{\alpha }}^t-{\mathit{\alpha }}^f}$), (2) the derivative of the spike rate function (${\displaystyle {\varphi }_{\mathrm{m}\mathrm{a}\mathrm{x}}{\sigma }^{\mathrm{\prime }}({\overline{{\mathit{V}}^0}}^f)}$), and (3) the filtered presynaptic spike trains (${\displaystyle {\overline{{\mathit{s}}^{\mathrm{i}\mathrm{n}\mathrm{p}\mathrm{u}\mathrm{t}}}}^f}$). All three of these terms are values that a real neuron could theoretically calculate using some combination of molecular synaptic tags, calcium currents, and back-propagating action potentials.

One aspect of this update rule that is biologically questionable, though, is the use of the term (${\displaystyle {\mathit{\alpha }}^t-{\mathit{\alpha }}^f}$). This requires a difference between plateau potentials that are separated by tens of milliseconds. How could such a signal be used by basal dendrite synapses to guide their updates? Plateau potentials can drive bursts of spikes (Larkum et al., 1999), which can propagate to basal dendrites (Kampa and Stuart, 2006). Since the plateau potentials are similar to rate variables (i.e. a sigmoid applied to the voltage), the number of spikes during the bursts, ${\displaystyle {\mathit{N}}^f=[N_1^f,N_2^f,...,N_m^f]}$ and ${\displaystyle {\mathit{N}}^t=[N_1^t,N_2^t,...,N_m^t]}$, for the forward and target plateaus, respectively, could be sampled from a Poisson distribution with rate parameter equal to the plateau potential level:

If the distinct phases (forward and target) were marked by some global signal, $\varphi (t)$, that was communicated to all of the neurons, for example a neuromodulatory signal, the phase of a global oscillation, or some blanket inhibition signal, then we can imagine an internal cellular memory mechanism in the basal dendrites of the $i^{th}$ neuron, $M_i$ (e.g. a molecular signal like the activity of an enzyme, the phosphorylation level of some protein, or the amount of calcium released from intracellular stores), which could be differentially sensitive to the inter-spike interval of bursts, depending on $\varphi$. So, for example, if we define:

where x indicates the forward or target phase. Then, the change in $M_i$ from before the bursts occur to afterwards would be, on average, proportional to the difference (${\displaystyle {\mathit{\alpha }}^t-{\mathit{\alpha }}^f}$), and could be used to calculate the weight updates.

However, this is highly speculative, and it is not clear that such a mechanism would be present in real neurons. We have outlined the mathematics here to make the reality of implementing the current model explicit, but we would predict that the brain would have some alternative method for calculating differences between top-down inputs at different times, for example by using somatostatin positive interneurons that are preferentially sensitive to bursts and which target the apical dendrite (Silberberg and Markram, 2007). We are ultimately agnostic as to this mechanism, and so, it was not included in the current model.

In order to extend our algorithm to deeper networks with multiple hidden layers, our model incorporates direct synaptic connections from the output layer to each hidden layer. Thus, each hidden layer receives feedback from the output layer through its own separate set of fixed, random weights. For example, in a network with two hidden layers, both layers receive the feedback from the output layer at their apical dendrites through backward weights ${\displaystyle {\mathit{Y}}^0}$ and ${\displaystyle {\mathit{Y}}^1}$. The local targets at each layer are then given by:

where the superscripts ${}^0$ and ${}^1$ denote the first and second hidden layers, respectively, and the superscript ${}^2$ denotes the output layer.

The local loss functions at each layer are:

where $L^2$ is the loss at the output layer. The learning rules used by the hidden layers in this scenario are the same as in the case with one hidden layer.

For each of the three network sizes that we present in this paper, a grid search was performed in order to find good learning rates. We set the learning rate for each layer by stepping through the range $[0.1,0.3]$ with a step size of $0.02$. For each combination of learning rates, a neural network was trained for one epoch on the 60, 000 training examples, after which the network was tested on 10,000 test images. The learning rates that gave the best performance on the test set after an epoch of training were used as a basis for a second grid search around these learning rates that used a smaller step size of $0.01$. From this, the learning rates that gave the best test performance after 20 epochs were chosen as our learning rates for that network size.

In all of our simulations, we used a learning rate of 0.19 for a network with no hidden layers, learning rates of 0.21 (output and hidden) for a network with one hidden layer, and learning rates of 0.23 (hidden layers) and 0.12 (output layer) for a network with two hidden layers. All networks with one hidden layer had 500 hidden layer neurons, and all networks with two hidden layers had 500 neurons in the first hidden layer and 100 neurons in the second hidden layer.

For all simulations described in this paper, the neural networks were trained on classifying handwritten digits using the MNIST database of 28 pixel $\times$ 28 pixel images. Initial feedforward and feedback weights were chosen randomly from a uniform distribution over a range that was calculated to produce voltages in the dendrites between $-6$ - $12$ V.

Prior to training, we tested a network’s initial performance on a set of 10,000 test examples. This set of images was shuffled at the beginning of testing, and each example was shown to the network in sequence. Each input image was encoded into Poisson spiking activity of the 784 input neurons representing each pixel of the image. The firing rate of an input neuron was proportional to the brightness of the pixel that it represents (with spike rates between $[0$ - ${\varphi }_{\max }]$. The spiking activity of each of the 784 input neurons was received by the neurons in the first hidden layer. For each test image, the network underwent only a forward phase. At the end of this phase, the network’s classification of the input image was given by the neuron in the output layer with the greatest somatic potential (and therefore the greatest spike rate). The network’s classification was compared to the target classification. After classifying all 10,000 testing examples, the network’s classification error was given by the percentage of examples that it did not classify correctly.

Following the initial test, training of the neural network was done in an on-line fashion. All 60,000 training images were randomly shuffled at the start of each training epoch. The network was then shown each training image in sequence, undergoing a forward phase ending with a plateau potential, and a target phase ending with another plateau potential. All feedforward weights were then updated at the end of the target phase. At the end of the epoch (after all 60,000 images were shown to the network), the network was again tested on the 10,000 test examples. The network was trained for up to 60 epochs.

For each training example, a minimum length of 50 ms was used for each of the forward and target phases. The lengths of the forward and target training phases were determined by adding their minimum length to an extra length term, which was chosen randomly from a Wald distribution with a mean of 2 ms and scale factor of 1. During testing, a fixed length of 500 ms was used for the forward transmit phase. Average forward and target phase voltages were calculated after a settle duration of $\mathrm{\Delta }{t}_s=30$ ms from the start of the phase.

For simulations with randomly sampled plateau potential times (Figure 5—figure supplement 1), the time at which each neuron’s plateau potential occurred was randomly sampled from a folded normal distribution ($\mu =0,{\sigma }^2=3$) that was truncated ($\text{max}=5$) such that plateau potentials occurred between 0 ms and 5 ms before the start of the next phase. In this scenario, the average apical voltage in the last 30 ms was averaged in the calculation of the plateau potential for a particular neuron.

The time-step used for simulations was $dt=1$ ms. At each time-step, the network’s state was updated bottom-to-top beginning with the first hidden layer and ending with the output layer. For each layer, dendritic potentials were updated, followed by somatic potentials, and finally their spiking activity. Table 1 lists the simulation parameters and the values that were used in the figures presented.

All code was written using the Python programming language version 2.7 (RRID: SCR_008394) with the NumPy (RRID: SCR_008633) and SciPy (RRID: SCR_008058) libraries. The code is open source and is freely available at https://github.com/jordan-g/Segregated-Dendrite-Deep-Learning (Guerguiev, 2017). The data used to train the network was from the Mixed National Institute of Standards and Technology (MNIST) database, which is a modification of the original database from the National Institute of Standards and Technology (RRID: SCR_006440) (Lecun et al., 1998). The MNIST database can be found at http://yann.lecun.com/exdb/mnist/. Some of the simulations were run on the SciNet High-Performance Computing platform (Loken et al., 2010).