1 Keywords:
Spiking Neural Networks, Leaky Integrate and Fire, Generalization, Hausdorff Dimension, logSTDP, addSTDP, multSTDP, Bayesian Optimization
2 Introduction
A Spiking Neural Network (SNN) (maass1997networks; gerstner2002spiking; pfeiffer2018deep) is a neuroinspired machine learning (ML) model that mimics the spikebased operation of the human brain (bi1998synaptic). The Spike Timing Dependent Plasticity (STDP) is a policy for unsupervised learning in SNNs (bell1997synaptic; magee1997synaptically; gerstner2002mathematical). The STDP relates the expected change in synaptic weights to the timing difference between postsynaptic spikes and presynaptic spikes (feldman2012spike). Recent works using STDP trained SNNs have demonstrated promising results as an unsupervised learning paradigm for various tasks such as object classification and recognition (she2021heterogeneous; diehl2015unsupervised; kheradpisheh2018stdp; masquelier2009competitive; lee2018deep; mozafari2019bio).
The generalizability is a measure of how well an ML model performs on test data that lies outside of the distribution of the training samples (kawaguchi2017generalization; neyshabur2017exploring)
. The generalization properties of Stochastic Gradient Descent (SGD) based training for deep neural network (DNN) have received significant attention in recent years
(poggio2019theoretical; allen2018learning; allen2019can). The dynamics of SGD have been studied via models of stochastic gradient Langevin dynamics with an assumption that gradient noise is Gaussian (simsekli2020fractional). Here SGD is considered to be driven by a Brownian motion. Chen et al. showed that SGD dynamics commonly exhibit exhibits rich, complex dynamics when navigating through the loss landscape (chen2020anomalous). Recently Gürbüzbalaban et al. (gurbuzbalaban2020heavy) and Hodgkinson et al. (hodgkinson2021multiplicative)have simultaneously shown that the law of the SGD iterates can converge to a heavytailed stationary distribution with infinite variance when the stepsize
is large and/or the batchsize is small. These results form a theoretical basis for the origins of the observed heavytailed behavior of SGD in practice. The authors proved generalization bounds for SGD under the assumption that its trajectories can be wellapproximated by the Feller Process (capocelli1976transformation), a Markovbased stochastic process. Modeling the trajectories of SGD using a stochastic differential equation (SDE) under heavytailed gradient noise has shed light on several interesting characteristics of SGD.In contrast, the generalizability analysis of STDP trained SNNs, although important, has received much less attention. Few studies have shown that, in general, the biological learning process in the human brain has significantly good generalization properties (zador2019critique; sinz2019engineering)
. However, none of them have characterized the generalization of an STDPtrained SNN using a mathematical model. There is little understanding of how hyperparameters of the STDP process impact the generalizability of the trained SNN model. Moreover, the generalization of STDP cannot be characterized by directly adopting similar studies for SGD. For example, SGD has been modeled as a Feller process for studying generalizability.
Rossum et al. showed that random variations arise due to the variability in the amount of potentiation (depression) between the preand postsynaptic events at fixed relative timing (van2000stable)
. At the neuron level, fluctuations in relative timing between preand postsynaptic events also contribute to random variations
(roberts2000computational). For many STDP learning rules reported in the literature, the dynamics instantiate a Markov process ((bell1997synaptic); (markram1997regulation); (bi1998synaptic); (han2000reversible); (van2000stable)); changes in the synaptic weight depend through the learning rule only on the current weight and a set of random variables that determine the transition probabilities. However, recent literature have shown that weight update using STDP is better modeled as an OrnsteinUhlenbeck process
(cateau2003stochastic), (aceituno2020spiking), (legenstein2008learning).As described by Camuto et al. (camuto2021fractal), fractals are complex patterns, and the level of this complexity is typically measured by the Hausdorff dimension (HD) of the fractal, which is a notion of dimension. Recently, assuming that SGD trajectories can be wellapproximated by a Feller Process, it is shown that the generalization error can be controlled by the Hausdorff dimension of the trajectories of the SDE simsekli2020hausdorff. That is, the ambient dimension that appears in classical learning theory bounds is replaced with the Hausdorff dimension. The fractal geometric approach presented by Simsekli et al. can capture the low dimensional structure of fractal sets and provides an alternative perspective to the compressionbased approaches that aim to understand why over parametrized networks do not overfit.
This paper presents a model to characterize the generalizability of the STDP process and develops a methodology to optimize hyperparameters to improve the generalizability of an STDPtrained SNN. We use the fact that the sample paths of a Markov process exhibit a fractallike structure (xiao2003random). The generalization error over the sample paths is related to the roughness of the random fractal generated by the driving Markov process which is measured by the Hausdorff dimension (simsekli2020hausdorff) which is in turn dependent on the tail behavior of the driving process.
Normally, validation loss of a model on a testing set is used to characterize the accuracy of that model. However, the validation loss is dependent on the choice of the test set, and does not necessarily give a good measure of the generalization of the learning process. Therefore, generalization error a model is generally measured normally measured by comparing the difference between training and testing accuracy  a more generizable model has less difference between training and testing accuracy (goodfellow2016deep). However, such a measure of generalization error requires computing the validation loss (i.e. testing accuracy) for a given test dataset. To optimize the generalizability of the model, we need an objective function that measures the generalizability of the learning process. If the validation loss is used as a measure, then for each iteration of the optimization, we need to compute this loss by running the model over the entire dataset, which will significantly increase the computation time. We use Hausdorff Dimension as a measure of generalizability of the STDP process to address the preceding challenges. First, the Hausdorff measure characterizes the fractal nature of the sample paths of the learning process itself and does not depend on the testing dataset. Hence, HD provides a better (i.e. testset independent) measure of the generalization of the learning process. Second, HD can be computed during the training process itself and does not require running inference on the test data set. This reduces computation time per iteration of the optimization process.
We model the STDP learning as an OrnsteinUhlenbeck process which is a Markov process and show that the generalization error is dependent on the Hausdorff dimension of the trajectories of the STDP process. We use the SDE representation of synaptic plasticity and model STDP learning as a stochastic process that solves the SDE.
Using the Hausdorff dimension we study the generalization properties of an STDP trained SNN for image classification. We compare three different STDP processes, namely, logSTDP, addSTDP, and multSTDP, and show that the logSTDP improves generalizability. We show that modulating hyperparameters of the STDP learning rule and learning rate changes the generalizability of the trained model. Moreover, using logSTDP as an example, we show the hyperparameter choices that reduce generalization error increases the convergence time, and training loss, showing a tradeoff between generalizability and the learning ability of a model. Motivated by the preceding observations, we develop a Bayesian optimization technique for determining the optimal set of hyperparameters which gives an STDP model with the least generalization error. We consider an SNN model with 6400 learning neurons trained using the logSTDP process. Optimizing the hyperparameters of the learning process using Bayesian Optimization shows a testing accuracy of and a generalization error of
on the MNIST dataset. This shows a mean increase of almost 40% in generalization performance for a mean drop of about 1% in testing accuracy in comparison to randomly initialized training hyperparameters. In order to further evaluate the learning methodologies, we also evaluated them on the more complex Fashion MNIST dataset and observed a similar trend.
3 Materials and Methods
3.1 Background
3.1.1 Spiking Neural Networks
We chose the leaky integrateandfire model of a neuron where the membrane voltage is described by
where is the resting membrane potential; and
are the equilibrium potentials of excitatory and inhibitory synapses, respectively; and
and are the conductances of excitatory and inhibitory synapses, respectively. The time constant is longer for excitatory neurons than for inhibitory neurons. When the neuron’s membrane potential crosses its membrane threshold, the neuron fires, and its membrane potential is reset. Hence, the neuron enters its refractory period and cannot spike again for the duration of the refractory period.Synapses are modeled by conductance changes, i.e., synapses increase their conductance instantaneously by the synaptic weight when a presynaptic spike arrives at the synapse, otherwise, the conductance decays exponentially. Thus, the dynamics of the conductance can be written as:
(1) 
If the presynaptic neuron is excitatory, the dynamics of the conductance is with the time constant of the excitatory postsynaptic potential being . On the other hand, if the presynaptic neuron is inhibitory, it’s synaptic conductance is given as and the time constant of the inhibitory postsynaptic potential as .
3.1.2 STDP based Learning Methods
Spiketimingdependent plasticity is a biologically plausible learning model representing the time evolution of the synaptic weights as a function of the past spiking activity of adjacent neurons.
In a STDP model, the change in synaptic weight induced by the preand postsynaptic spikes at times are defined by:
(2) 
where the learning rate
determines the speed of learning. The Gaussian white noise
with zero mean and variance describes the variability observed in physiology. The function describes the long term potentiation (LTP) and depression (LTD) based on the the relative timing of the spike pair within a learning window , and is defined by:(3) 
The shape of the weight distribution produced by STDP can be adjusted via the scaling functions in 3 that determine the weight dependence. We study three different types of STDP processes, namely, addSTDP, multSTDP, and logSTDP. All STDP models follow the equations 2 and 3, however, they have different scaling functions () in 3 as discussed below. The weight distributions of these three different STDP processes at the end of the last training iteration are shown in Fig. 1.At the beginning of the training iterations, the distribution is uniform for all three reflecting on the weight initialization conditions. Additive STDP gives rise to strong competition among synapses, but due to the absence of weight dependence, it requires hard boundaries to secure the stability of weight dynamics. To reach a stability point for the addSTDP, we followed the analysis done by Gilson et al. (gilson2011stability) and Burkitt et al. (burkitt2004spike) and chose the fixed point
. Fig. 1 approximates the probability density function based on the weight distributions of the different STDP models. We are using a Gaussian KDE to get this pdf from the empirical weight distribution obtained. Given
independent realizations from an unknown continuous probability density function (p.d.f.) on, the Gaussian kernel density estimator is defined as
(4) 
where
(5) 
is a Gaussian p.d.f. (kernel) with location and scale . The scale is usually referred to as the bandwidth. Much research has been focused on the optimal choice of , because the performance of
as an estimator of
depends crucially on its value (jones1992progress), (sheather1991reliable).Logarithmic STDP (logSTDP) (gilson2011stability). The scaling functions of logSTDP is defined by:
(6)  
(7) 
In this study, is chosen such that LTP and LTD in logSTDP balance each other for uncorrelated inputs, namely Therefore, will also be referred to as the ’fixed point’ of the dynamics. Since depression increases sublinearly (blue solid curve for in Fig.1 ), noise in logSTDP is weaker than that for multSTDP for which depression increases linearly (orange dashed curve for in Fig.1).
The weight dependence for LTD in logSTDP is similar to multSTDP for , i.e., it increases linearly with However, the LTD curve becomes sublinear for and determines the degree of the loglike saturation. For larger , LTD has a more pronounced saturating loglike profile and the tail of the synaptic weight distribution extends further.
We choose the function for LTP to be roughly constant around , such that the exponential decay controlled by only shows for . Thus, the scaling functions , ,and are the hyperparameters that can be tuned to model an efficient logSTDP learning model.
Additive STDP (addSTDP) (song2000competitive). It is weight independent and is defined by the following scaling functions:
(8) 
with such that LTD overpowers LTP. The drift due to random spiking activity thus causes the weights to be depressed toward zero, which provides some stability for the output firing rate. (gilson2011stability). For the simulations, we are using a fast learning rate that is synonymous to a high level of noise, and more stability. This requires a stronger depression. Thus, we use and .
Multiplicative STDP (multSTDP) (van2000stable). The multiplicative STDP has a linear weight dependence for LTD and constant LTP:
(9)  
(10) 
The equilibrium mean weight is then given by For the simulations we use and . This calibration corresponds to a similar neuronal output firing rate to that for logSTDP in the case of uncorrelated inputs.
3.1.3 Generalization  Hausdorff Dimension and Tail Index Analysis
Recent works have discussed the generalizability of SGD based on the trajectories of the learning algorithm. Simsekli et al. (simsekli2020hausdorff) identified the complexity of the fractals generated by a Feller process that approximates SGD. The intrinsic complexity of a fractal is typically characterized by a notion called the Hausdorff dimension (le2019hausdorff; lHorinczi2019multifractal), which extends the usual notion of dimension to fractional orders. Due to their recursive nature, Markov processes often generate random fractals (xiao2003random)
. Significant research has been performed in modern probability theory to study the structure of such fractals
(bishop2017fractals; khoshnevisan2009fractals; khoshnevisan2017macroscopic; yang2018multifractality). Thus, the STDP learning method follows an OrnsteinUhlenbeck (OU) process which is a special type of Lévy process. Again, the Hausdorff Dimension measures the roughness of the fractal patterns of the sample paths generated by the stochastic process which is measured using the tail properties of the Lévy measure of the OU process. Lévy measure is a Borel measure on satisfying . The OrnsteinUhlenbeck process which is a Lévy process can thus be characterized by the triplet where denotes a constant drift, is a positive semidefinite matrix and is the Lévy measure as defined above. Thus, taking Lévy processes as stochastic objects, their sample path behavior can be characterized by the Hausdorff dimension which is in turn measured using the BGindices. Thus, the generalization properties of the STDP process can be modeled using the Hausdorff dimension of the sample paths of the OU process. We formally define the Hausdorff dimension for the OrnsteinUhlenbeck process modeling the STDP learning process in Section 4.2.3.2 STDP as a Stochastic Process
In this paper, we evaluate the generalization properties of an STDP model using the concept of the Hausdorff dimension. In this section, we discuss the learning methodology of STDP and how the plasticity change can be modeled using a stochastic differential equation. The state of a neuron is usually represented by its membrane potential which is a key parameter to describe the cell activity. Due to external input signals, the membrane potential of a neuron may rise until it reaches some threshold after which a spike is emitted and transferred to the synapses of neighboring cells. To take into account the important fluctuations within cells, due to the spiking activity and thermal noise, in particular, a random component in the cell dynamics has to be included in mathematical models describing the membrane potential evolution of both the preand postsynaptic neurons similar to the analysis shown by Robert et al., (robert2020stochastic). Several models take into account this random component using an independent additive diffusion component, like Brownian motion, of the membrane potential . In our model of synaptic plasticity, the stochasticity arises at the level of the generation of spikes. When the value of the membrane potential of the output neuron is at , a spike occurs at rate where
is the activation function
(chichilnisky2001simple). In particular, we consider the instants when the output neuron spikes are represented by an inhomogeneous Poisson process as used by Robert et al. (robert2020stochastic). Thus, in summary, (1) The presynaptic spikes are modeled using a Poisson process and hence, there is a random variable added to membrane potential. (2) The postsynaptic spikes are generated using a stochastic process based on the activation function. Hence, STDP, which depends on the preand postsynaptic spike times can be modeled using a stochastic differential equation (SDE). Hence, we formulate the STDP as a SDE. The SDE of a learning algorithm shares similar convergence behavior of the algorithm and can be analyzed more easily than directly analyzing the algorithm.Mathematical Setup: We consider the STDP as an iterative learning algorithm which is dependent on the dataset and the algorithmic stochasticity . The learning process returns the entire evolution of the parameters of the network in the time frame where being the parameter value returned by the STDP learning algorithm at time . So, for a given training set , the learning algorithm will output a random process indexed by time which is a trajectory of iterates.
In the remainder of the paper, we will consider the case where the STDP learning process is chosen to be the trajectories produced by the OrnsteinUhlenbeck (OU) process , whose symbol depends on the training set . More precisely, given , the output of the training algorithm will be the random mapping , where the symbol of is determined by the drift , diffusion matrix , and the Lévy measure , which all depend on . In this context, the random variable represents the randomness that is incurred by the OU process. Finally, we also define the collection of the parameters given in a trajectory, as the image of , i.e., and the collection of all possible parameters as the union . Note that is still random due to its dependence on .
We consider the dynamics of synaptic plasticity as a function of the membrane potential and the synaptic weight .
Time Evolution of Synaptic Weights and Plasticity Kernels As described by Robert et al., (robert2020stochastic), the time evolution of the weight distribution depends on the past activity of the input and output neurons. It may be represented using the following differential equation:
(11) 
where are two nonnegative processes where the first one is associated with potentiation i.e., increase in W and the latter is related to the depression i.e., decrease in W. The function needs to be chosen such that the synaptic weight stays at alltime in its definition interval . The function can thus be modified depending on the type of implementation of STDP that is needed. Further details regarding the choice of for different types of STDP is discussed by Robert et al. (robert2020stochastic)
When the synaptic weight of a connection between a presynaptic neuron and a postsynaptic neuron is fixed and equal to , the time evolution of the postsynaptic membrane potential is represented by the following stochastic differential equation (SDE) (robert2020stochastic):
(12) 
where is the left limit of at , and is the exponential decay time constant of the membrane potential associated with the leaking mechanism. The sequence of firing instants of the presynaptic neuron is assumed to be a Poisson point process on with the rate . At each presynaptic spike, the membrane potential is increased by the amount If the synapse is said to be excitatory, whereas for the synapse is inhibitory. The sequence of firing instants of the postsynaptic neuron is an inhomogeneous Poisson point process on whose intensity function is . The drop of potential due to a postsynaptic spike is represented by the function When the postsynaptic neuron fires instate , its state just after the spike is .
Uniform Hausdorff Dimension: The Hausdorff dimension for the training algorithm is a notion of complexity based on the trajectories generated by . Recent literature has shown that the synaptic weight update using an STDP rule can be approximated using a type of stochastic process called the OrnsteinUhlenbeck process which is a type of Markov process (cateau2003stochastic), (aceituno2020spiking), (legenstein2008learning). Hence, we can infer that the STDP process will also have a uniform Hausdorff dimension for the trajectories. We use the Hausdorff Dimension of the sample paths of the STDP based learning algorithm which has not been investigated in the literature.
Let be the class of functions which are right continuous, monotone increasing with and such that there exists a finite constant such that
(13) 
A function in is often called a measure function. For , the Hausdorff measure of is defined by
(14) 
where denotes the open ball of radius centered at The sequence of balls satisfying Eq. 14 is called an covering of . We know that is a metric outer measure and every Borel set in is measurable. Thus, the function is called the Hausdorff measure function for if .
It is to be noted here that in Eq. 14, we only use coverings of by balls, hence is usually called a spherical Hausdorff measure in the literature. Under Eq. 13, is equivalent to the Hausdorff measure defined by using coverings by arbitrary sets. The Hausdorff dimension of is defined by
(15) 
The STDP learning process is modeled using the SDEs for the temporal evolution of the synaptic weights and the membrane potential given in Eqs. 11, 12. Considering the empirical observation that STDP exhibits a diffusive behavior around a local minimum (baity2018comparing), we take to be the local minimum found by STDP and assume that the conditions of Proposition hold around that point. This perspective indicates that the generalization error can be controlled by the BG index of the Lévy process defined by ; the subsymbol of the process 11 around . The choice of the SDE 11 imposes some structure on , which lets us express in a simpler form. This helps us in estimating the BG index for a general Lévy process. As shown by Simsekli et al., there is a layerwise variation of the tailindex of the gradient noise in a DNNbased multilayer neural network(simsekli2019tail). Thus, for our STDP model we assume that around the local minimum , the dynamics of STDP will be similar to the Lévy motion with frozen coefficients: . Also assuming the coordinates corresponding to the same layer have the same tailindex around , the BG index can be analytically computed as (meerschaert2005dimension). We note here that and thus, determines the order of the generalization error. Using this simplification, we can easily compute , by first estimating each by using the estimator proposed by Mohammadi et al. (mohammadi2015estimating), that can efficiently estimate by using multiple iterations of the STDP learning process.
As explained by Simsekli et al. (simsekli2020hausdorff), due to the decomposability property of each dataset , the stochastic process for the synaptic weights given by behaves like a Lévy motion around a local point . Because of this locally regular behavior, the Hausdorff dimension can be bounded by the BlumenthalGetoor (BG) index (blumenthal1960some), which in turn depends on the tail behavior of the Lévy process. Thus, in summary, we can use the BGindex as a bound for the Hausdorff dimension of the trajectories from the STDP learning process. Now, as the Hausdorff dimension is a measure of the generalization error and is also controlled by the tail behavior of the process, heaviertails imply less generalization error (simsekli2020hausdorff), (simsekli2020fractional).
To further demonstrate the heavytailed behavior of the OrnsteinUhlenbeck process (a type of stable Lévy process ) that characterizes the STDP learning mechanism, we plot its trajectories and their corresponding pdf distributions. We plot these for varying values of the stability factor of the Lévy process , . We hence also plotted the probability density function of the Lévy processes to show the heavytailed nature of the Lévy trajectories as the tail index decreases. Figure 2 shows a continuoustime random walk performed using the OU random process in 3D space. In the Figure, X, Y, Z are random variables for the stable distributions generated using the OU process. Fig. 3, shows the corresponding probability density function of the OU process for varying values of corresponding to the different trajectories shown in Fig. 2. From Figs. 2,3, that as the OU process becomes heaviertailed (i.e., decreases), and the Hausdorff dimension gets smaller.
3.3 Optimal Hyperparameter Selection
Using the Hausdorff Dimension as a metric for the generalizability of the learning process, we formulate an optimization process that selects the hyperparameters of the STDP process to improve the generalizability of the models. The Hausdorff dimension is a function of the hyperparameters of the STDP learning process. Thus, we formulate an optimization problem to select the optimal hyperparameters of the STDP using the Hausdorff dimension of the STDP learning process as the optimization function. Now, since the BGindex is the upper bound of the Hausdorff dimension, as discussed earlier, we in turn aim to optimize the BGindex of the STDP stochastic process. The optimization problem aims to get the optimal set of hyperparameters of the STDP process that can give a more generalizable model without looking at the test set data. Now, given an STDP process, we aim to tune its hyperparameters to search for a more generalizable model. Let us define where is the set of hyperparameters of the STDP process, . Let denote the domain of the hyperparameter . The hyperparameter space of the algorithm is thus defined as . Now, we aim to design the optimization problem to minimize the Hausdorff Dimension of the learning trajectory for the STDP process. This is calculated over the last training iteration of the model, assuming that it reaches near the local minima. When trained with on the training data , we denote the algorithm’s Hausdorff dimension as . Thus, using Kfold cross validation, the hyperparameter optimization problem for a given dataset is to given as follows:
(16) 
We choose the Sequential Modelbased Bayesian Optimization (SMBO) technique for this problem (feurer2015initializing). SMBO constructs a probabilistic model of based on point evaluations of and any available prior information. It then uses that model to select subsequent configurations to evaluate. Given a set of hyperparameters for an STDP learning process , we define the point functional evaluation as the calculation of the BG index of with the hyperparameters . The BG index gives an upper bound on the Hausdorff dimension of the learning process. In order to select its next hyperparameter configuration using model SMBO uses an acquisition function which uses the predictive distribution of model at arbitrary hyperparameter configurations . This function is then maximized over to select the most useful configuration to evaluate next. There exists a wide range of acquisition functions (mockus1978application), all of whom aim to tradeoff between exploitation and exploration. The acquisition function tries to balance between locally optimizing hyperparameters in regions known to perform well and trying hyperparameters in a relatively unexplored region of the space.
In this paper, for the acquisition function, we use the expected improvement (mockus1978application) over the best previouslyobserved function value attainable at a hyperparameter configuration where expectations are taken over predictions with the current model :
(17) 
4 Results
4.1 Experimental Setup
We empirically study the generalization properties of the STDP process by designing an SNN with 6400 learning neurons for handwritten digit classification using the MNIST dataset. The MNIST dataset contains training examples and test examples of pixel images of the digits . It must be noted here that the images from the ten classes in the MNIST dataset are randomized so that there is a reinforcement of the features learned by the network.
Architecture. We use a twolayer SNN architecture as done by the implementation of Diehl et al., (diehl2015unsupervised) as shown in Figure 4. The first layer is the input layer, containing neurons with one neuron per image pixel. The second layer is the processing layer, with an equal number of excitatory and inhibitory neurons. The excitatory neurons of the second layer are connected in a onetoone fashion to inhibitory neurons such that each spike in an excitatory neuron will trigger a spike in its corresponding inhibitory neuron. Again, each of the inhibitory neurons is connected to all excitatory ones, except for the one from which it receives a connection. This connectivity provides lateral inhibition and leads to competition among excitatory neurons. There is a balance between the ratio of inhibitory and excitatory synaptic conductance to ensure the correct strength of lateral inhibition. For a weak lateral inhibition, the conductance will not have any influence while an extremely strong signal would ensue that one dominant neuron suppresses the other ones.
The task for the network is to learn a representation of the dataset on the synapses connecting the input layer neurons to the excitatory layer neurons. The excitatoryinhibitory connectivity pattern creates competition between the excitatory neurons. This allows individual neurons to learn unique filters where the single most spiked neuron in each iteration updates its synapse weights to match the current input digit using the chosen STDP rule. Increasing the number of neurons allows the network to memorize more examples from the training data and recognize similar patterns during the test phase.
Homeostasis. The inhomogeneity of the input leads to different firing rates of the excitatory neurons, and lateral inhibition further increases this difference. However, all neurons should have approximately equal firing rates to prevent single neurons from dominating the response pattern and to ensure that the receptive fields of the neurons differentiate. To achieve this, we employ an adaptive membrane threshold resembling intrinsic plasticity (zhang2003other). Specifically, each excitatory neuron’s membrane threshold is not only determined by but by the sum , where is increased every time the neuron fires and is exponentially decaying (querlioz2013immunity). Therefore, the more a neuron fires, the higher will be its membrane threshold and in turn, the neuron requires more input to a spike soon. Using this mechanism, the firing rate of the neurons is limited because the conductancebased synapse model limits the maximum membrane potential to the excitatory reversal potential , i.e., once the neuron membrane threshold is close to (or higher) it will fire less often (or even stop firing completely) until decreases sufficiently.
Input Encoding. The input image is converted to a Poisson spike train with firing rates proportional to the intensity of the corresponding pixel. This spike train is then presented to the network in an alltoall fashion for 350 ms as shown in Fig. 4. The maximum pixel intensity of 255 is divided by 4, resulting in input firing rates between 0 and 63.75 Hz. Additionally, if the excitatory neurons in the second layer fire less than five spikes within 350 ms, the maximum input firing rate is increased by 32 Hz and the example is presented again for 350 ms. This process is repeated until at least five spikes have been fired during the entire time the particular example was presented.
Training and STDP Dynamics Analysis. To train the network, we present digits from the MNIST training set to the network. It is to be noted that, before presenting a new image, no input to any variable of any neuron is given for a time interval of 150 ms. This is done to decay to their resting values. All the synaptic weights from input neurons to excitatory neurons are learned using the STDP learning process as described in Sec. 3.1.2. To improve simulation speed, the weight dynamics are computed using synaptic traces such that every time a presynaptic spike arrives at the synapse, the trace is increased by 1 (morrison2007spike). Otherwise, decays exponentially. When a postsynaptic spike arrives at the synapse the weight change is calculated based on the presynaptic trace as described in section 3.1.2. To evaluate the model, we divide the training set into 100 divisions of 600 images each and check the model performance after each such batch is trained using the STDP learning model. In the remainder of the paper, we call this evaluation of the model after 600 images as one iteration.
Inference. After the training process is done, the learning rate is set to zero and each neuron’s spiking threshold is fixed. A class is assigned to each neuron based on its highest response to the ten classes of digits over one presentation of the training set. This is the first time labels are used in the learning algorithm, which makes it an unsupervised learning method. The response of the classassigned neurons is used to predict the digit. It is determined by taking the mean of each neuron response for every class and selecting the class with the maximum average firing rate. These predictions are then used to measure the classification accuracy of the network on the MNIST test set
Hyperparameter  logSTDP  addSTDP  multSTDP 
Synaptic Delay  0.75ms  0.75ms  0.75ms 
Synaptic epsp  1ms  1ms  1ms 
Synaptic epsp  5ms  5ms  5ms 
Number of correlated pools  4  4  4 
Number of neurons per pool  50  50  50 
Spiking Rate of inputs  10Hz  10Hz  10Hz 
Learning rate ()  
STDP Apre (LTP) time constant  17ms  17ms  17ms 
STDP Apre (LTD) time constant  34ms  34ms  34ms 
Increase in (LTP), on prespikes  1.0  1.0  1.0 
Increase in (LTD), on postspikes  0.5  0.55  100 
LTD Curvature Factor ()  5  N/A  N/A 
Exponential LTP Decay factor ()  50  N/A  N/A 
Threshold Fixedpoint weight ()  0.006  N/A  N/A 
Computation of Generalization Error and Hausdorff Dimension. We empirically study the generalization behavior of STDP trained SNNs. We vary the hyperparameters of the STDP learning process which controls the LTP/LTD dynamics of the STDP learning algorithm. Table 1 shows the hyperparameters for various STDP processes. We trained all the models for 100 training iterations. In this paper, we consider the synaptic weight update to follow a Lévy process, i.e., continuous with discrete jumps similar to the formulation of Stein et al. (stein1965theoretical) and Richardson et al. (richardson2010firing). As discussed in Section 3.2, the generalizability can be measured using the Hausdorff dimension which is bounded by BGindex.
Therefore, the BGindex is computed in the last iteration when all the neurons have learned the input representations. We also compute the generalization error as the difference between the training and test accuracy. we study the relations between BGindex, generalization error, testing accuracy, and convergence behavior of the networks.
4.2 Analysis of Generalizability of STDP Processes
Impact of Scaling Functions. Kubota et al. showed that the scaling functions play a vital role in controlling the LTP/LTD dynamic of the STDP learning method (kubota2009modulation). In this subsection, we evaluate the impact of scaling functions (i.e. in the equation 3) on the generalizability properties of the STDP methods. We define the ratio of the postsynaptic scaling function to the presynaptic one (i.e. in add, mult, and log STDP equations), hereafter referred to as the scaling function ratio (SFR), as our variable. Kubota et al. has shown that the learning behavior is best when this SFR lies between the range of 0.9 to 1.5. Hence, we also modulate the SFR within this set interval. Table 2 shows the impact of scaling function on Hausdorff dimension (measured using BGindex), generalization error, and testing accuracy. We observe that a smaller SFR leads to a lower Hausdorff dimension and a lower generalization error, while a higher ratio infers a less generalizable model. However, a higher SFR marginally increases the testing accuracy. The analysis shows confirms that a higher Hausdorff dimension (i.e. a higher BGindex) corresponds to a higher generalization error, as discussed in section 3.2, justifying the use of BGindex as a measure of the generalization error.
logstdp  addstdp  multstdp  











2.1  1.352  6.8  89.92  1.969  9.7  88.17  1.824  8.1  89.26  
1.7  1.294  6.2  89.98  1.911  9.3  88.12  1.797  7.6  89.15  
1.2  1.209  5.9  89.79  1.875  8.9  88.09  1.702  7.0  88.99  
0.9  1.174  5.7  89.26  1.799  8.6  88.10  1.633  6.5  88.87 
Impact of the Learning Rate. One of the major parameters that control the weight dynamics of the STDP processes is the learning rate i.e. the variable in equation 2. In this subsection, we evaluate the effect of the learning rate on the generalizability of the STDP process. We have summarized the results in Table 3. We observe that a larger learning rate converges to a more generalizable solution. This can be attributed to the fact that a higher learning rate inhibits convergence to sharper minimas; rather facilitates convergence to a flatter one resulting in a more generalizable solution. We also observe the monotonic relation between the BGindex and the generalization error.
logstdp  addstdp  multstdp  











0.2  1.312  6.6  89.12  1.844  9.1  87.69  1.769  7.9  88.38  
0.15  1.255  6.1  89.03  1.783  8.9  87.53  1.648  7.4  88.19  
0.1  1.112  5.3  88.35  1.698  8.5  87.11  1.596  6.8  87.95  
0.05  1.068  5.0  88.01  1.632  8.2  87.02  1.512  6.3  87.82 
Impact of the STDP models on Generalizability. In this section, we compare the three different STDP models, namely, addSTDP, multSTDP, and logSTDP to its generalization abilities with changing SFR (scaling function ratio) and learning rate. The results are summarized in Tables 2, 3. In all the above cases we see that the logSTDP process has a significantly lower generalization error compared to the other two STDP methods. The difference between the generalizability of various STDP models comes from the nature of the stochastic distribution of weights generated by different models.
Gilson et al. (gilson2011stability) has discussed that addSTDP (gutig2003learning) can rapidly and efficiently select synaptic pathways by splitting synaptic weights into a bimodal distribution of weak and strong synapses. However, the stability of the weight distribution requires hard bounds due to the resulting unstable weight dynamics. In contrast in multSTDP (rubin2001equilibrium)
, weightdependent update rules can generate stable unimodal distributions. However, multSTDP weakens the competition among synapses leading to only weakly skewed weight distributions. The probability distributions of the three different STDP models are shown in Fig.
1. On the other hand, logSTDP proposed by Gilson et al. (gilson2011stability) bypass these problems by using a weightdependent update rule while does not make the other synapses weak as in multSTDP. The logSTDP results in a lognormal solution of the synaptic weight distribution as discussed by Gilson et al. (gilson2011stability). A lognormal solution has a heavier tail and thus a smaller Hausdorff dimension leading to a lower generalization error. A detailed comparison of the weight distributions of the three types of STDP processes can be found in the paper by Gilson et al. (gilson2011stability). We further evaluated the training loss for iterations for the different STDP models. The results are plotted in Fig. 8. From the figures we see that the logSTDP outperforms the addSTDP and the multSTDP in terms of training loss for either case, thus demonstrating the performance ofImpact on Different Datasets: To demonstrate the generalizability of the STDP models, we also tested its performance on FashionMNIST (xiao2017fashion) which is an MNISTlike fashion product dataset with 10 classes. FashionMNIST shares the same image size and structure of the training and testing splits as MNIST but is considered more realistic as its images are generated from front look thumbnail images of fashion products on Zalando’s website via a series of conversions. Therefore, FashionMNIST poses a more challenging classification task than MNIST. We preprocessed the data by normalizing the sum of a single sample gray value because of the high variance among examples. The results for SFR and learning rate is shown in Table 4. We see that as seen in MNIST datasets, for the FashionMNIST also, the logSTDP method has a lower generalization error corresponding to a lower BG Index.
logSTDP  addSTDP  multSTDP  











c+/c_ =0.9  1.322  10.02  82.43  1.788  13.87  79.16  1.573  11.38  81.92  
1.204  9.93  81.75  1.692  11.73  79.76  1.489  10.11  80.35 
4.3 Generalizability vs Trainability Tradeoff
In this section, we study the relations between the generalizability and trainability of a learning model. For the sake of brevity, we only focus on the logSTDP process as it has shown better generalizability compared to addSTDP and multSTDP.
We plot the training loss as a function of the time evolution of the synaptic weights trained with the STDP learning method.
Since STDP is an online unsupervised learning algorithm, there is no formal concept of training loss. So, to evaluate the performance of the model, we define the training loss as follows: We divide the MNIST dataset into 100 divisions, with each division consisting of 600 images. We evaluate the model after training the model on each subset of the full training dataset and this is considered as one training iteration. We train the SNN model using STDP with this limited training dataset. After the training is done, we set the learning rate to zero and fix each neuron’s spiking threshold. Then, the image of each class is given as input to the network, and the total spike count of all the neurons that have learned that class is calculated. Hence, the spike counts are normalized by dividing the number of spikes by the maximum number of spike counts and the crossentropy loss was calculated across all the classes. This is defined as the training loss. To show the confidence interval of each training iteration, we evaluated each of the partially trained models on the test dataset 5 times, randomizing the order of images for each of the test runs. We see from Figure 6, that initially, as the model was trained with fewer images, the variance in training loss was high demonstrating a low confidence interval. However, as the model is trained on a larger training set, the variance decreases as expected.
Table 2 and Figure 4 show the training loss versus the number of iterations for the logSTDP process for various SFR. We see that the shows a lower generalization error and almost similar testing accuracy, compared to the other SFRs. The results show that increasing the SFR increases the generalization error. If the presynaptic scaling function is stronger than the postsynaptic scaling function (i.e. is lower), it implies that the synaptic weights of the neurons gradually decay. Since we have the images in the MNIST dataset randomized over the ten classes, the more important features which help in the generalization ability of the network over unknown data are reinforced so that the network does not forget these filters as shown by Panda et al. panda2017asp. Thus, the network only forgets the less important features and preserves the more important ones, hence making it more generalizable. Since the neuron forgets some features which would help to fit better into the current dataset, it affects its training/testing accuracy as can be seen in Tables 2,3. Thus, the model learns the more important features and is essentially more generalizable.
Note here that the training loss for the STDP processes all reach their convergence around iteration 60  i.e., images after that added little information for the training of the model. The models here are not optimized and hence optimizing the hyperparameters can also help in reducing the number of images required for extracting enough information from the training dataset. Thus, if SFR is too high, training gets messed up since a neuron starts spiking for multiple patterns, in which case there is no learning. As the SFR value increases from 1, the SNN tends to memorize the training pattern and hence the generalization performance is poor. On the other hand, if when SFR is less than 1 but is close to 1, it is hard to memorize the training patterns as the STDP process tends to forget the patterns which are nonrepeating, leading to better generalization.
On the other hand, if the postsynaptic scaling function is stronger than the presynaptic one (i.e. is higher), then the neurons tend to learn more than one pattern and respond to all of them. Similar results can be verified from Figure 7 where the learning rate was varied instead of the SFR. In this study as well we observe that a higher learning rate, although leads to faster convergence and lower training loss, leads to a less generalizable model. Hence, we empirically observe that hyperparameters of STDP models that lead to better generalizability can also make an SNN harder to train.
4.4 Results of Hyperparameter Optimization
In Section 4.2, we empirically showed that the Hausdorff dimension is a good measure of the generalizability of the model and it can be very efficiently controlled using the hyperparameters of the STDP learning process. In this section, we show the application of our Bayesian optimization strategy to search for the optimal hyperparameters to increase the generalizability of an STDPtrained SNN model. For the sake of brevity, we demonstrate the application of Bayesian optimization on the logSTDP process. Table 5 shows the set of hyperparameters that are optimized and their optimal values obtained by our approach. It should be clarified here that the hyperparameters are not necessarily the absolute global optimum but a likely local optimum found in the optimization algorithm. The optimized logSTDP model results in a training accuracy of 93.75%, testing the accuracy of 90.49%, and a BG Index of 0.718 for the MNIST dataset.
We study the behavior of Bayesian optimization. Each iteration in the Bayesian optimization process corresponds to a different set of hyperparameters for the logSTDP model. Each such iteration is called a functional evaluation. For each functional evaluation, the Bayesian Optimization trains the SNN with the corresponding set of hyperparameters of the logSTDP model and measures the BGindex of the weight dynamics of the trainedSNN. Figure 9(a) shows the change in the BGIndex as a function of a number of the function evaluations of the search process. It is to be noted here that at each functional evaluation, we train the network with the STDP learning rule with the chosen hyperparameters and estimate the Hausdorff dimension from the trained network. We see that for the optimal set of hyperparameters, the BG Index converges to . Figure 9(b) shows the corresponding training accuracy of the model with the change in iteration number. We see that the logSTDP configurations during Bayesian optimization that have a higher BG Index (i.e. a higher generalization error) also have has a higher training accuracy. These results further validate our observations on the generalizability vs trainability tradeoff for a logSTDP trained SNN.
Hyperparameter  Domain  logSTDP  addSTDP  
Before BO  After BO  Before BO  After BO  

[0.05, 0.2]  0.1  0.063  0.1  0.017  

[0.1, 1]  0.5  0.581  0.5  0.632  

3  5  N/A  N/A  

45  57  N/A  N/A  

[0, 1]  0.5  0.244  N/A  N/A  

(0,1](0,1] 






[10,20][20, 40]  15,30  17, 36  15,30  18,22  
Testing Accuracy  91.41  90.65  88.68  89.77  
Generalization Error  5.79  3.17  8.29  4.61 
Comparison with AddSTDP: In order to compare the performance of the logSTDP, we performed a similar analysis using the addSTDP model. The results of the Bayesian Optimization for the addSTDP and the logSTDP are plotted in Figure 10. We see that the logSTDP process outperforms the addSTDP model in terms of both training/testing accuracy and the generalization error thus showing the robustness of the logSTDP process. We see that though the addSTDP has a higher training accuracy, and a comparable test accuracy, its generalization error is higher compared to the logSTDP method.
5 Discussion
In this paper, we presented the generalization properties of the spike timingdependent plasticity (STDP) models. A learning process is said to be more generalizable if it can extract features that can be transferred easily to unknown testing sets thus decreasing the performance gap between the training and testing sets. We provide a theoretical background for the motivation of the work treating the STDP learning process as a stochastic process (an OrnsteinUhlenbeck process) and modeling it using a stochastic differential equation. We control the hyperparameters of the learning method and empirically study their generalizability properties using the Hausdorff dimension as a measure. From Tables 2,3 and corresponding Figures 5, we observed that the Hausdorff Dimension is a good measure for the estimation of the generalization error of an STDPtrained SNN. We compared the generalization error and testing error for the logSTDP, addSTDP, and multSTDP models, as shown in Tables 2, 3. We observed that the lognormal weight distribution obtained from the logSTDP learning process leads to a more generalizable STDPtrained SNN with a minimal decrease in testing accuracy. In this paper, when we refer to a model as more generalizable, we mean there is a smaller difference between the training and testing performance, i.e., the generalization error. The objective of the paper was to get a model which is more generalizable in the sense that the performance of the network on unknown datasets should not differ much from its performance in the training dataset. It is to be noted that in this paper we are using the generalization error as the metric of generalizability of the network. Generalization error is not a measure of absolute accuracy, but rather the difference between training and testing datasets. As such, we see that models which have lower generalization error extract lesser and more important features compared to less generalizable models. However, we see that with this reduced set of features, the model has almost no drop in the testing accuracy, showing the generalizability of the model at comparable accuracy. Thus, we get a model which is more generalizable in the sense that the performance of the network on unknown datasets does not differ much from its performance in the training dataset. As such, these ’more generalizable’ models, extract lesser and more important features compared to less generalizable models. However, we see that with this reduced set of features, the model has almost no drop in the testing accuracy, showing the generalizability of the model as we can see from Tables 2, 3. This phenomenon can be explained using the observations of Panda et al. panda2017asp on how the ability to forgets boosts the performance of spiking neuron models. The authors showed that the final weight value towards the end of the recovery phase is greater for the frequent input. The prominent weights will essentially encode the features that are common across different classes of old and new inputs as the preneurons across those common feature regions in the input image will have frequent firing activity. This eventually helps the network to learn features that are more common with generic representations across different input patterns. This extraction of more generalizable features can be interpreted as a sort of regularization wherein the network tries to generalize over the input rather than overfitting such that the overall accuracy of the network improves. However, due to this regularization, we see that the training performance of the network decreases. However, since the model is more generalizable, the testing performance remains almost constant as seen in Figures 10, 11. We further observe that the logSTDP models which have a lower Hausdorff dimension and hence have lower generalization error, have a worse trainability i.e., takes a long time to converge during training and also converges to a higher training loss. The observations show that an STDP model can have a tradeoff between generalizability and trainability. Finally, we present a Bayesian optimization problem that minimizes the Hausdorff dimension by controlling the hyperparameter of a logSTDP model leading to a more generalizable STDPtrained SNN.
Future work on this topic will consider other models of STDP. In particular, the stochastic STDP rule where the probability of synaptic weight update is proportional to the time difference of the arrival of the pre and postsynaptic spikes has shown improved accuracy over deterministic STDP studied in this paper. The trajectories of such a stochastic STDP model will lead to a Feller process as shown by Kuhn (helson2017new). Hence, in the future, we will perform a similar Hausdorff dimensionbased analysis for generalization for the stochastic STDP model. Moreover, in this work, we have only considered the hyperparameters of the STDP model to improve the generalizability of the SNN. An important extension is to consider the properties of the neuron dynamics, which also controls the generation of the spikes and hence, weight distribution. The choice of the network architecture will also play an important role in the weight distribution of the SNN. Therefore, a more comprehensive optimization process that couples hyperparameters of the STDP dynamics, neuron dynamics, and network architecture like convolutional SNN (kheradpisheh2018stdp) and heterogeneous SNN (she2021heterogeneous) will be interesting future work.
Funding
This material is based on work sponsored by the Army Research Office and was accomplished under Grant Number W911NF1910447. The views and conclusions contained in this document are those of the authors and should not be interpreted as representing the official policies, either expressed or implied, of the Army Research Office or the U.S. Government.
Comments
There are no comments yet.