Year: 2021

Abstract:
A dynamical system of spiking neurons with only feedforward connections can classify spatiotemporal patterns without recurrent connections. However, the theoretical construct of a feedforward Spiking Neural Network (SNN) for approximating a temporal sequence remains unclear, making it challenging to optimize SNN architectures for learning complex spatiotemporal patterns. In this work, we establish a theoretical framework to understand and improve sequence approximation using a feedforward SNN. Our framework shows that a feedforward SNN with one neuron per layer and skip-layer connections can approximate the mapping function between any arbitrary pairs of input and output spike train on a compact domain. Moreover, we prove that heterogeneous neurons with varying dynamics and skip-layer connections improve sequence approximation using feedforward SNN. Consequently, we propose SNN architectures incorporating the preceding constructs that are trained using supervised backpropagation-throughtime (BPTT) and unsupervised spiking-timing-dependent plasticity (STDP) algorithms for classification of spatiotemporal data. A dual search-space Bayseian optimization method is developed to optimize architecture and parameters of the proposed SNN with heterogeneous neuron dynamics and skip-layer connections.

heading: 1 INTRODUCTION
text: Spiking neural network (SNN) (Ponulak & Kasinski, 2011) uses biologically inspired neurons and synaptic connections trainable with either biological learning rules such as spike-timingdependent plasticity (STDP) (Gerstner & Kistler, 2002) or statistical training algorithms such as backpropagation-through-time (BPTT) (Werbos, 1990). The SNNs with simple leaky integrate-andfire (LIF) neurons and supervised training have shown classification performance similar to deep neural networks (DNN) while being energy efficient (Kim et al., 2020b; Wu et al., 2019; Srinivasan & Roy, 2019). One of SNN’s main difference from DNN is that the neurons are dynamical systems with internal states evolving over time, making it possible for SNN to learn temporal patterns without recurrent connections. Empirical results on feedforward-only SNN models show good performance for spatiotemporal data classification, using either supervised training (Lee et al., 2016; Kaiser et al., 2020; Khoei et al., 2020), or unsupervised learning (She et al., 2021). However, while empirical results are promising, a lack of theoretical understanding of sequence approximation using SNN makes it challenging to optimize performance on complex spatiotemporal datasets.
In this work, we develop a theoretical framework for analyzing and improving sequence approximation using feedforward SNN. We view a feedforward connections of spiking neurons as a spike propagation path, hereafter referred to as a memory pathways (She et al., 2021), that maps an input spike train with an arbitrary frequency to an output spike train with a target frequency. Consequently, we argue that an SNN with many memory pathways can approximate a temporal sequence of spike trains with time-varying unknown frequencies using a series of pre-defined output spike trains with known frequencies. Our theoretical framework aims to first establish SNN’s ability to map frequencies of input/output spike trains within arbitrarily small error; and next, derive the basic principles for adapting neuron dynamics and SNN architecture to improve sequence approximation. The theoretical derivations are then investigated with experimental studies on feedforward SNN for spatiotemporal classifications. We adopt the basic design principles for improving sequence approx-
imation to optimize SNN architectures and study whether these networks can be trained to improve performance for spatiotemporal classification tasks. The key contributions of this work are:
• We prove that any spike-sequence-to-spike-sequence mapping functions on a compact domain can be approximated by feedforward SNN with one neuron per layer using skip-layer connections, which cannot be achieved if no skip-layer connection is used.
• We prove that using heterogeneous neurons having different dynamics and skip-layer connection increases the number of memory pathways a feedforward SNN can achieve and hence, improves SNN’s capability to represent arbitrary sequences.
• We develop complex SNN architectures using the preceding theoretical observations and experimentally demonstrate that they can trained with supervised BPTT and unsupervised STDP for classification on spatiotemporal data.
• We design a dual-search-space option for Bayesian optimization process to sequentially optimize network architectures and neuron dynamics of a feedforward SNN considering heterogeneity and skip-layer connection to improve learning and classification of spatiotemporal patterns.
We experimentally demonstrate that our network design principles coupled with the dual-searchspace Bayesian optimization improve classification performance on DVS Gesture (Amir et al., 2017), N-caltech (Orchard et al., 2015), and sequential MNIST. Results show that the design principles derived using our theoretical framework for sequence approximation can improve spatiotemporal classification performance of SNN.

heading: 2 RELATED WORK
text: Prior theoretical approaches to analyze SNN often focus on the storage and retrieval of precise spike patterns (Amit & Huang, 2010; Brea et al., 2013). There are also works that consider SNN for solving optimization problems (Chou et al., 2018; Binas et al., 2016) and works that analyze the dynamics of SNN (Zhang et al., 2019; Barrett et al., 2013). Those are different topics from the approximation of spike-sequence-to-spike-sequence mappings functions. SNN that incorporates excitatory and inhibitory signal is shown for its ability to emulate sigmoidal networks (Maass, 1997) and is theoretically capable of universal function approximation. Feedforward SNN with specially designed spiking neuron models (Iannella & Back, 2001; Torikai et al., 2008) have been demonstrated for function approximation, while for networks using LIF neurons, function approximation has been shown with only empirical results (Farsa et al., 2015). On the other hand, the existing works that has developed efficient training process for SNN and demonstrated classification performance comparable to deep learning models, have mostly used simpler and generic LIF neuron models (Lee et al., 2016; Kaiser et al., 2020; Kim et al., 2020b; Wu et al., 2019; Sengupta et al., 2019; Safa et al., 2021; Han et al., 2020). Therefore, this paper develops the theoretical basis for function approximation using feedforward SNN with LIF neurons, and studies applications of the developed theoretical constructs in improving SNN-based spatiotemporal pattern classification.
The effectiveness of heterogeneous neurons (She et al., 2021) and skip-layer connections (Srinivasan & Roy, 2019; Sengupta et al., 2019) in SNN has been empirically studied in the past. However, no theoretical approach has been presented to understand why such methods improve learning of spike sequences, and how to optimize a SNN’s architecture and parameters to effectively exploit these design constructs. It is possible to search for the optimal SNN configurations through optimization algorithms, but the large number of hyper-parameters for spiking neurons and network structure creates a high-dimensional search space that is long and difficult to solve. Bayesian optimization (Snoek et al., 2012) uses collected data points to make decisions on the next test point that could provide improvement, thus accelerates the optimization process. Prior works (Parsa et al., 2019; Kim et al., 2020a) have shown that SNN performance can also be effectively improved with Bayesian optimization. While those works consider a single or a few neuron parameters, the dualsearch-space Bayesian optimization proposed in this work optimizes both network architecture and neuron parameters efficiently by separating the discrete search spaces from the continuous search spaces.
t

heading: 3 APPROXIMATION THEORY OF FEEDFORWARD SNN

heading: 3.1 DEFINITIONS AND NOTATIONS
text: Definition 1 Neuron Response Rate γ For a spiking neuron n with membrane potential at vreset and input spike sequence with period tin, γ is the number of input spike n needs to reach vth.
Definition 2 Memory Pathways For a feedforward SNN with m layers, a memory pathway is defined as a spike propagation path from input to the output layer. Two memory pathways are considered distinct if the set of neurons contained in them is different.
Definition 3 Minimal Multi-neuron-dynamic (mMND) Network A densely connected network in which each layer has an arbitrary number of neurons that have different neuron parameters. All synapses from one pre-synaptic neuron have the same synaptic conductance.
Notations Neuron Delay tnd is the time for a spike from pre-synaptic neuron to arrive at its postsynaptic neurons, as shown in Figure 1(a). For a feedforward SNN with m layers, a skip-layer connection can be defined with source layer and target layer pair (ls, lt). The output feature map from source layer is concatenated to the original input feature map of the target layer. For the analysis of spike sequence in temporal space, the notation of Tmax and Tmin are defined as positive real numbers such that Tmax > Tmin. ϵ > 0 is the error of approximation.
Figure 1(a) shows two memory pathways receiving an input spike sequence with time-varying periods. As the neurons have different dynamics, the two memory pathways have different response to the input spike sequence. An example of mMND network with m layers and n neuron dynamics is shown in Figure 1(b). SNN with multilayer perceptron (MLP) structure can be considered a scaled-up mMND network with multiple neurons for each dynamic. A network with convolutional structure can be considered a scaled-up mMND network with duplicated connections in each layer. We analyze the correlation of network capacity and structure based on mMND networks. The design of neuron heterogeneity can also be implemented in MLP-SNN and Conv-SNN as described in Section 4. The analysis for network capacity can be extended to those networks according to their specific layer dimensions.

heading: 3.2 MODELING OF SPIKING NEURON
text: SNN consists of spiking neurons connected with synapses. The spiking neuron model studied in this work is leaky integrate-and-fire (LIF) as defined by the following equations:
τm dv
dt = a+RmI − v; v = vreset, if v > vthreshold (1)
Rm is membrane resistance, τm = RmCm is time constant and Cm is membrane capacitance. a is the resting potential. I is the sum of current from all input synapses that connect to the neuron. A spike is generated when membrane potential v cross threshold and the neuron enters refractory period r, during which the neuron maintains its membrane potential at vreset. The time it take for a pre-synaptic neuron to send a spike to its post-synaptic neurons is tnd. Neuron response rate γ is a property of a spiking neuron’s response to certain input spike sequence. We show how the value of γ can be evaluated below.
Remark For any input spike sequence, individual spike can be described with Dirac delta function δ(t− ti) where ti is the time of the i-th input spike. For the membrane potential of a spiking neuron receiving the input before reaching spiking threshold, with initial state at t = 0 with v = vreset, solving the differential equation (1) leads to:
v(t) = vresete − tτm + a(1− e− t τm ) + Rm τm e− t τm ∑ i G ∫ t 0 δ(t− ti)e t τm dt (2)
Here, G is the conductance of input synapse connected to the neuron, which is trainable. From (2), there exists a value of u such that vm(tu−1) < vthreshold and vm(tu) >= vthreshold. By evaluating (2) for u given neuron parameters and input spike sequence, the neuron response rate γ can be found.

heading: 3.3 APPROXIMATION THEOREM OF FEEDFORWARD SNN
text: To develop the approximation theorem for feedforward SNN, we first aim to understand the range of neuron response rate that can be achieved. We show with Lemma 1 that for any input spike sequence with periods in a closed interval, it is possible to set the neuron response rate γ to any positive integer. Based on this property, we show with Theorem 1 that by connecting a list of spiking neurons with certain γ sequentially and inserting skip-layer connections, approximation of spikesequence mapping functions can be achieved. To understand whether this capability of feedforward SNN relies on skip-layer connections, we develop Lemma 2 to prove that skip-layer connections are indeed necessary. In subsection 3.4 we investigate the correlation between approximation capability and network structures by analyzing the cutoff property of spiking neurons, which can change the network’s connectivity. In our analysis, we focus on two particular designs: heterogeneous network (Lemma 4) and skip-layer connection (Lemma 5), and show their impact on the number of distinct memory pathways in a network. All lemmas are formally proved in the appendix.
Lemma 1 For any input spike sequence with period tin in range [Tmin, Tmax], there exist a spiking neuron n with fixed parameters vth, vreset, a, Rm and τm, such that by changing synaptic conductance G, it is possible to set the neuron response rate γn to be any positive integer.
Proof Sketch. (Formal proof in Appendix A) Given an input spike sequence, it is possible to derive the highest possible membrane potential decay ∆v for any input tin ∈ [Tmin, Tmax] as a function of neuron parameters. We show that it is possible to make ∆v tends to zero by configuring the neuron parameters. Since the decay of v can be negligible, γn can be set to any positive integer by changing G.
Theorem 1 For any input and target output spike sequence pair with periods (tin, tout) ∈ [Tmin, Tmax]× [Tmin, Tmax], there exist a minimal-layer-size network with skip-layer connections that has memory pathway with output spike period function P (t) such that |P (tin)− tout| < ϵ.
Proof Sketch. (Formal proof in Appendix B) With skip-layer connections, there can be multiple memory pathways in a minimal-layer-size network as neurons can be either included or skipped through. Hence it is possible to create memory pathways with different delay times for each input spike in a network and by connecting the output of those memory pathways to a common neuron n′, spike sequence of any arbitrary period tint such that tint <= tin can be generated within ϵ. By setting γn′ > 1, the output from n′ receiving input spike sequence with tint is tn ′
out = γ · tint. Hence it is possible to achieve a network with output spike period P (t) such that |P (tin)− tout| < ϵ.
Lemma 2 With no skip-layer connection, there does not exist a minimal-layer-size network that has output spike period function P (t) such that for any input and target output spike sequence pair with periods (tin, tout) ∈ [Tmin, Tmax]× [Tmin, Tmax], |P (tin)− tout| < ϵ.
Proof Sketch. (Formal proof in Appendix C) A minimal-layer-size network without skip-layer connection has only one memory pathway. For a particular input spike sequence with period tin, different output period P (tin) can be achieved by changing γ of neurons in the memory pathway. We show that there exists two output spike periods P (tin) and P ′(tin), such that P (tin)− P ′(tin) is a constant value independent of network or neuron configurations, and there can be no P ′′(tin) such that P (tin) < P ′′(tin) < P ′(tin). Therefore, for any minimal-layer-size network, there exists tout within the range of (P (tin), P ′(tin)) such that |P (tin)− tout| < ϵ does not hold true.

heading: 3.4 NETWORK STRUCTURE AND MEMORY PATHWAYS
text: Based on Theorem 1, it is possible to approximate an input/output spike sequence mapping function using a minimal-layer-size network with specific configuration, which can be considered as a memory pathway. Since any bounded continuous function on a compact domain can be approximated to arbitrary accuracy using a piece-wise constant function, and it is possible to use a memory pathway to approximate each of the piece-wise constant function, with increasing number of distinct memory pathways, a feedforward SNN can achieve approximation of continuous functions with less error. In this subsection, we show that two SNN structural designs: heterogeneous network i.e. a network having neurons with different dynamics and adding skip-layer connections, a feedforward SNN has the capability to achieve more distinct memory pathways.
Cutoff Frequency of a Memory Path We first show the correlation of cutoff period and spiking neuron parameters with Lemma 3, which is proved in Appendix D.
Lemma 3 A spiking neuron has cutoff period ωc = τm ln( vreset−avreset−a+Rmτm G ) above which input spike sequence cannot cause the spiking neuron to spike.
Remark From Lemma 3, it can be observed that the cutoff period ωc of a neuron can be configured to any positive real number by changing the neuron parameters and synaptic conductance G. Further, with fixed G, ωc can be configured to any positive real number by changing the neuron parameters. Neurons that are in cutoff change the spike propagation path in a network as they send no output spikes. This creates different memory pathways without changing the connections in a network.
Heterogeneous Network If an mMND network has the same parameters for all neurons in each layer, the majority of the neurons are included in the same memory pathway, leading to the upper bound of number of distinct memory pathways to be limited. With Lemma 4, we show the relationship between the upper bound of the number of distinct memory pathways and the number of different neuron dynamics in an mMND network.
Lemma 4 For an mMND network with m layers and {λ1, λ2, ...λm} number of different neuron dynamics in each layer, the least upper bound of the number of distinct memory pathways is ∏m i=1 λi.
Proof Sketch. (Formal proof in Appendix D) For an mMND network, it is possible to have neurons with different ωc in each layer, which creates λi number of different neuron activation states for layer i. Across all network layers the highest possible number of different neuron activation states is therefore the product of λ of each layer. Since neurons in cutoff do not propagate spikes, they can be removed from a memory pathway. This leads to ∏m i=1 λi as the least upper bound of the number of memory pathways.
Compared to a network with homogeneous neuron parameters, in which the upper bound of number of distinct memory pathways is λm, Lemma 4 indicates that heterogeneous network increases the maximum achievable number of distinct memory pathways in a feedforward SNN.
Skip-layer Connection We show that adding skip-layer connection increases the upper bound of the number of memory pathways in a network with Lemma 5.
Lemma 5 For a mMND network with m layers and {λ1, λ2, ...λm} different neuron dynamics in each layer and a skip-layer connection made between layer la and lb, s.t. a, b ∈ {1, 2, ...m} and (b − a) > 1, the least upper bound of the number of memory pathways is ∏m i=1 λi + ( ∏a i=1 λi ·∏m
i=b λi)
Proof Sketch. (Formal proof in Appendix D) Compared to the network considered in Lemma 4, by adding a skip-layer connection, there are additional possible neuron activation states in the network that result from the cutoff of neurons in layers between la and lb. Without layers between la and lb in the spike propagation path, the least upper bound of the number of memory pathways is increased by the maximum number of distinct memory pathways in a network that has layers {l1, l2, ..., la, lb, lb+1, ..., lm} connected sequentially.

heading: 4 SNN ARCHITECTURE DESIGN USING APPROXIMATION THEORY
text: In this section, we discuss design of SNN architectures as inspired by the developed approximation theory for feedforward SNN, with more details in Appendix G.
Network Template for BPTT Training For BPTT training, the network template is shown in Figure 2(a). The heterogeneity in neuron dynamics is implemented by using the multi-neuron-dynamic layers, which can either be convolutional or fully connected. The multi-neuron-dynamic layers use different neuron parameters for spiking neurons in each neuron dynamic module, which contains a certain number of feature maps for convolutional layers or a certain number of neurons for fullyconnected layers. There are two types of synapses between layers: transferred synapses marked as black dashed arrows and learned synapses marked as red solid arrows. The conductance of learned synapses is optimized by the BPTT algorithm during training, and the transferred synapses have the same conductance as the learned synapses from the same pre-synaptic neuron. During forward pass, neurons in each layer receive the same input features and respond differently based on their neuron dynamics to generate different output features. During back-propagation, only conductance of the learned synapses are updated. Skip-layer connection is implemented with the output spike matrix from source layer concatenated to the original input spike matrix of the target layer. The skip-layer connection has the same implementation as the regular connection between consecutive layers, with both learned and transferred synapses (Figure 2(a)). The last layer of the network is a full-connected layer with homogeneous dynamic to generate prediction labels.
Network Template for STDP Learning For networks trained with STDP, the template is shown in Figure 2(b). Each layer contains a learner module and a memory module. Learner modules use homogeneous neuron dynamic that is suitable for STDP learning, and memory modules consist of neurons with different dynamics. There are also two types of synapses: transferred synapses and learned synapses. Between two layers, memory modules are connected with transferred synapses and memory modules are connected to learner modules with learned synapses. Leaner modules be-
tween layers are not directly connected. STDP training proceeds as a layer-by-layer process. During training of the first layer, conductance of synapses connecting neurons in memory module to neurons in learner module are learned with STDP using all training data without supervision. Then, the learned conductance is transferred to the transferred synapses in the same layer. During training of the second layer, layer 1 memory module perceives the same input features and generate different output features with the heterogeneous neurons. This lay-by-layer process is repeated until the layer before the final layer finishes learning. The final linear layer is then fine-tuned using stochastic gradient descent (SGD) based on spike frequency array from the last multi-neuron-dynamic layer generated based on the labeled data. Skip-layer connection is implemented by connecting the memory module of the source layer to the target layer. The connections are made with the two types of synapses and follow the same training process as the consecutive layers.
Dual-search-space Bayesian Optimization Bayesian optimization uses Gaussian process to model the distribution of an objective function, and an acquisition function to decide points to evaluate. For data points in a target dataset x ∈ X and the corresponding label y ∈ Y , an SNN with network structure V and neuron parameters W acts as a function fV,W(x) that maps input data x to predicted label ỹ. The optimization problem in this work is defined as
minV,WP where P = ∑
x∈X,y∈Y L(y, fV,W(x)) (3)
V contains the number of layers Nlayers, the number of memory dynamics Ndynmaic and skiplayer connection configuration variables Nskip, Lstart and Lend, each controlling the number of skip-layer connections, the first layer and last layer to implement skip-layer connections. All of the values are discrete. W contains the values for a, τm and Rm in (1), which are continuous. We separate the discrete and continuous search spaces by implementing a dual-search-space optimization process, where V is first optimized with fixed, manually tuned neuron parameters. After an optimal structure is found, W are optimized for the selected V . Details on the configurations of the optimization process are listed in the appendix. To achieve Bayesian optimization with constraints, we implement a modified expected improvement (EI) acquisition function similar to the one shown by Gardner (Gardner et al., 2014), which uses a Gaussian process to model the feasibility indicator due to its high evaluation cost. In this work, since the constraint function can be explicitly defined, we use feasibility indicator that is directly evaluated. The modified EI function is defined as:
Ic(W) = ∆(W) · max{0,P(W)− P(W+)} (4)
where W is the network configuration containing W and V . W+ is the test point that provided the best result. ∆(W) is the explicitly defined indicator function that takes the value of 1 when all constraints are satisfied and 0 otherwise.

heading: 5 EXPERIMENTS

heading: 5.1 EXPERIMENT SETTINGS
text: Datasets tested in the experiment include the DVS Gesture (Amir et al., 2017), which is a human gesture dataset captured by DVS cameras, and the N-Caltech101 (Orchard et al., 2015), which is an event-based version of the Caltech101 dataset. The proposed method is also tested for MLP-style SNN on the sequential MNIST dataset presented row-by-row. We also vary the amount of labeled data used during training ranging from using 100% labeled data for training down to 10% labeled data (30% for N-Caltech101) during training. Note, during STDP training networks always uses the entire but un-labeled training dataset; however, only the fraction of the labeled data is used for supervised fine-tuning of the last layer. Comparison is made for DVS Gesture and N-Caltech101 with prior works including ConvLSNN, which is a combination of convolutional SNN and recurrent SNN with long and short-term neurons trained with BPTT (Salaj et al., 2020), DECOLLE (Kaiser et al., 2020), which uses surrogate gradient to train a convolutional feedforward SNN, HATS (Sironi et al., 2018), which implements time surfaces and SVM for classification and H-SNN (She et al., 2021) which uses STDP to train a convolutional SNN with two neuron dynamics. Additional function approximation experiments are discussed in Appendix E.

heading: 5.2 EFFECT OF DUAL-SEARCH-SPACE BAYESIAN OPTIMIZATION
text: We compare the proposed dual-search-space Bayesian optimization with regular Bayesian optimization using a single search space for network validation error over 5 runs. The result from the NCaltech101 dataset is shown in Figure 3. It can be observed that the two optimization approaches achieve similar minimum validation error after convergence. By separating the search spaces, the proposed optimization process reaches convergence faster than regular single-search-space optimization. It is also worth noting that, between the two stages in the optimization process for BPTT training, the first stage accounts for more reduction in validation error than the second stage. This indicates that optimizing network structure causes more impact to BPTT training than optimizing neuron parameters, which is potentially due to the reason that network structures more heavily affects the number of memory pathways in the network than neuron parameters. On the other hand, for STDP training where learning behavior is sensitive to the dynamic of spiking neurons, the reduction of validation error is more equally shared between the two optimization stages. Over the 5 runs, among all network configurations achieved after the dual-search-space optimization converges, we compare the configuration with the lowest number of trainable parameters against baseline models. The specific configurations for the optimized networks are listed in Table 1. It can be observed that for BPTT algorithm, the optimized networks have more layers than the STDP trained networks, and the optimal values found for neuron parameters are highly distinct for the two training methods.

heading: 5.3 ABLATION STUDIES
text: To investigate the effect of using multiple neuron dynamics, we apply the same dual-search-space Bayesian optimization process for networks that have homogeneous neuron dynamic for the same
number of evaluations as the proposed design. Similarly, to study the contribution to performance gain from skip-layer connections, the Bayesian optimization process is used for network templates without skip-layer connections. The optimization process runs for the same number of evaluations as the proposed design. From the results shown in Table 2, it can be observed that compared to baselines, the proposed networks achieve the best accuracy for all datasets. Specifically, when homogeneous network is used, the performance of STDP trained network is noticeably lower than the proposed method for DVS Gesture and N-Caltech101. For BPTT training, using heterogeneous network and skip-layer connection shows different level of benefit for each dataset. For sequential MNIST which has less complexity, the improvement from using heterogeneous neurons and skip-layer connections is not as significant.

heading: 5.4 COMPARISON WITH PRIOR WORKS
text: DVS Gesture With 100% labels available during training, the proposed network trained with BPTT demonstrates higher accuracy than all tested networks with the less trainable parameters than all baselines in Table 3 (top). The proposed network trained with STDP has slightly lower accuracy than ConvLSNN and DECOLLE when 100% labels are used; for reduced-label training it outperforms all network (including H-SNN).
N-Caltech101 As shown in Table 3 (bottom), the proposed network trained with BPTT outperforms all baselines with both 70% and 100% training labels and also has less number of trainable parameters. The un-supervised learning models i.e., H-SNN and the proposed network with STDP, show considerably lower performance (more than what was observed for DVS Gesture) than supervised ones when 100% labels are available; However, the proposed network with STDP shows better performance than H-SNN, and outperforms all network when available labels are below 50%.

heading: 6 CONCLUSION
text: We develop a theoretical basis to understand and optimize the ability of a feedforward SNN to approximate a temporal sequence. We analytically show how heterogeneity and skip-layer connections can improve sequence approximation with SNN, and empirically demonstrate their impact on spatiotemporal learning. It is well-known in neuroscience that, heterogeneity (De Kloet & Reul, 1987) and irregular connectivity (Eickhoff et al., 2018) are intrinsic properties of human brains. Our analysis shows that incorporating such concepts within artificial SNN is beneficial for designing high-performance SNN for classification of spatiotemporal data.

heading: A SNN DYNAMICS
text: Remark For a sequentially connected neuron list with m neurons all with γ = 1 and neuron delay tnd, an input spike at time t leads the neuron list to generate an output spike at time t+mtnd
Remark For any input sequence with period tin to a spiking neuron with response rate γ such that γ > 1, if refractory period is set to r < tin, the neuron can exit refractory period before the next spike arrives.
Lemma 1 For any input spike sequence with period tin in range [Tmin, Tmax], there exist a spiking neuron n with fixed parameters vth, vreset, a, Rm and τm, such that by changing synaptic conductance G, it is possible to set the neuron response rate γn to be any positive integer.
Proof. For a given input spike sequence period tin, consider the maximum possible membrane potential decay that can be reached within a period of tin. From (1), when I = 0, dvdt < 0 and | dv dt | increases with higher v. Hence the maximum decay of v is reached when initial membrane potential v(t = 0) → v−th and the neuron decays for period tin = Tmax. The decayed membrane potential v(t = Tmax) can be derived by solving the differential equation (1) for v(t) = vth at t = 0:
v(t = Tmax) = vthe −Tmaxτm − ae− Tmax τm + a (5)
It is possible to have a spiking neuron with Rm, a and τm such that ∆v, defined as
∆v = v(t = Tmax)− v(t = 0) = vthe− Tmax τm − ae− Tmax τm + a− vth (6)
,
tends to zero. With this configuration, since the the highest possible decay of membrane potential is negligible, for any target γ, it is possible to set G such that
G = vth − vreset
γ (7)
The proof is complete.
Note, input sequence that has no spike has tin tending to infinity, thus violating the bounded constraint mentioned in Lemma 1.

heading: B PROOF OF THEOREM 1
text: Theorem 1 For any input and target output spike sequence pair with periods (tin, tout) ∈ [Tmin, Tmax]× [Tmin, Tmax], there exist a minimal-layer-size network with skip-layer connections that has memory pathway with output spike period function P (tin) such that |P (tin)− tout| < ϵ. Proof. For any given (tin, tout), first consider the condition where tin > tout. It is possible to construct a minimal-layer-size network N connecting m spiking neurons with neuron response rate γ = 1 sequentially, denoted as a m-tuple of neurons {n1, n2, ..., nm}. Since any configuration of skip-layer connection with source layer and target layer pair (ls, lt), such that ls ∈ [1,m − 2], and lt ∈ [ls + 2,m], can be added, it is possible to add a (m− 2)-tuple of skip-layer connections
Ssl = {(i,m) ∀ i ∈ {1, 2, 3, ...,m− 2}} (8)
Denote the synaptic conductance for all the skip-layer connections as a (m− 2)-tuple
SGsl = {Gsl1 , Gsl2 , Gsl3 , ..., Gslm−2} (9)
For any tout < tin, it is possible to find a k-tuple of synaptic conductance
S′Gsl = {G sl i , G sl 2i, G sl 3i, ..., G sl ki} s.t. i = ⌊ tout tnd ⌋ and k = ⌊m− 2 i ⌋ (10)
Set synaptic conductance in SGsl \S′Gsl to 0. Then set the conductance of synapse connecting nm−1 and nm to 0. In such way, The output spikes from network N has period
P (tin) = ⌊ tout tnd ⌋ · tnd (11)
For given ϵ, it is possible to choose tnd such that tnd < 2ϵ, therefore satisfying |P (tin)− tout| < ϵ. m can be chosen as
m = Tmax − Tmin
tnd (12)
or equivalently:
m = Tmax − Tmin
2ϵ (13)
Since Tmax−Tmin2ϵ is finite, m is finite.
For tin < tout, using N as described above, it is possible to achieve output spike with period within ϵ of any period in (0, tin]. For a given tout, assume the configuration in neuron list N has output spike interval t′int such that kt ′ int = tout, where k is a positive integer. From Lemma 1, it is possible to set G for a neuron nm+1 such that its neuron response delay satisfies γnm+1 = k for input spike period t′int. A new network, denoted as N
′, can be formed by connecting nm+1 to the output of N . N ′ has output spike with period P (tin) = kt′int = tout. Hence, to reach the given ϵ, it requires neuron list N to have output spike interval tint such that
|tint − t′int| < ϵ
k (14)
Since k is finite, (14) can be achieved.
For tin >= tout, it is possible to configure network N ′ such that tint satisfies |tint − tout| < ϵ, and the value of γnm+1 set to 1, hence |P (tin)− tout| < ϵ can be achieved. The proof is complete.

heading: C PROOF OF LEMMA 2
text: Lemma 2 With no skip-layer connection, there does not exist a minimal-layer-size network that has output spike period function P (tin) such that for any input and target output spike sequence pair with periods (tin, tout) ∈ [Tmin, Tmax]× [Tmin, Tmax], |P (tin)− tout| < ϵ. Proof. A minimal-layer-size network N with m layers can be denoted as a m-tuple of neurons {n1, n2, ..., nm} connected sequentially. Since no skip-layer connection exists, there is only one distinct memory pathway that contains all neurons {n1, n2, ..., nm}. Denote the set of neuron response rate corresponding to each neuron in N as
Γ = {γ1, γ2, ..., γm} (15)
For a given input sequence with tin, denote the timing of the first spike as tin1 , consider the output spike sequence for network with
γi = 1 ∀ γi ∈ Γ (16)
The first output spike from N has timing t̃out1 = t in 1 +mtnd, and the second output spike has timing t̃out2 = t in 1 + tin +mtnd. It can be easily derived that the period of the output spike sequence is
P (tin) = tin (17)
Also consider the output spike sequence for network with
γj = 2 for any j ∈ {1, 2, 3, 4, ...,m} and γi = 1 ∀ i ∈ ({1, 2, 3, 4, ...,m} \ {j}) (18)
Following the same process, the period of the output spike sequence is
P (tin) = 2tin (19)
Since the smallest increase to any γi is by 1, there is no set of values for γ such that the network output spike sequence has period P (tin) satisfying tin < P (tin) < 2tin. Since within the range (tin, 2tin), there exists values of tout such that |P (tin) − tout| < ϵ does not hold. The proof is complete.

heading: D MEMORY PATHWAYS IN SNN
text: In this section we analyze the increase to the upper bound of the number of memory pathways in a network by using heterogeneous networks and skip-layer connections.
Lemma 3 A spiking neuron has cutoff period ωc = τm ln( vreset−avreset−a+Rmτm G ) above which input spike sequence cannot cause the spiking neuron to spike.
Proof. Consider (2), since the membrane potential increases at time of ti and decays otherwise, solving for t = ti and the equation can be expanded:
vm(t i) = vresete
ti τm + a(1− e ti τm ) + Rm τm Ge ti−t1 τm + Rm τm Ge ti−t2 τm + ...+ Rm τm G (20)
For input with frequency f , ti+1 − ti = ∆t = 1f , subtracting membrane potential values at two consecutive ti provides:
∆vm = vm(t i+1)− vm(ti) = vreset(e
ti+1 τm − e ti τm )− a(e ti+1 τm − e ti τm ) + Rm τm Ge ti+1−t1 τm (21)
setting time of first input spike t1 to zero leads to:
∆vm = e ti τm ((e ∆t τm − 1)(vreset − a) + Rm τm Ge ∆t τm ) (22)
As e ti τm > 0, and the term ((e ∆t τm )−1)(vreset−a)+Rmτm Ge ∆t τm ) does not depend on ti, the polarity of ∆vm does not change with time. vm is either strictly increasing, staying the same or decreasing with higher ti. This indicates that, when ∆vm ≤ 0 the post-synaptic neuron can never spike regardless of how many pre-synaptic spike it receives. ∆vm ≤ 0 when input spike period tin satisfies
tin ≥ τm ln( vreset − a
vreset − a+ Rmτm G ) (23)
Therefore, the cutoff period of the neuron is
ωc = τm ln( vreset − a
vreset − a+ Rmτm G ) (24)
The proof is complete.
In the following proof, we consider cutoff frequency, fc = 1ωc of spiking neurons. From (24) it can be easily derived that fc can be configured to any positive real number by changing the neuron parameters.
Lemma 4 For an mMND network with m layers and {λ1, λ2, ...λm} number of different neuron dynamics in each layer, the least upper bound of the number of memory pathways is ∏m i=1 λi.
Proof.
Denote the set of neurons in layer l with distinct neuron dynamics as
Sln = {nl1, nl2, nl3, ..., nlλl} (25)
Since the network is mMND, |Sln| = λl. Denote the set of cutoff frequency corresponding to each neuron in Sln as
F lc = {f nl1 c , f nl2 c , f nl3 c , ..., f nlλl c } (26)
Since neurons in Sln can have different neuron parameters, from Lemma 3, it is possible to set the parameters such that all entries of F lc are distinct. Hence, there exists a permutation π such that
f nlπ(1) c < f nlπ(2) c < f nlπ(3) c ... < f nlπ(λl) c (27)
Denote the input spike frequency to layer l as f lin, neuron n l i is a part of a valid memory pathway in the network if
f nli c ≤ f lin (28)
If the input spike frequency is f lin ≥ f nlπ(λl) c , all neurons in Sln can be a part of a valid memory pathway of the network. For input spike frequency such that f nlπ(i) c ≤ f lin < f nlπ(i+1) c , i neurons can be a part of a valid memory pathway of the network.
For any input to the network fin ∈ [Fmin, Fmax], denote the number of ways different neurons in Sln can be part of valid memory pathways as kl. The total number of distinct memory pathways K in the network is
K = m∏ l=1 kl (29)
Since 0 ≤ kl ≤ λl,
K ≤ Kmax = m∏ l=1 λl (30)
For any layer l ∈ {1, 2, 3, ...,m}, the input f lin is bounded by [F lmin, F lmax]. kl = λl can be achieved by setting f nlπ(1) c and f nlπ(λl) c such that:
f nlπ(1) c = F l min and f nlπ(λl) c = F l max (31)
Hence the bound is tight for (30). The proof is complete.
Lemma 5 For a mMND network with m layers and {λ1, λ2, ...λm} different neuron dynamics in each layer and a skip-layer connection made between layer la and lb, s.t. a, b ∈ {1, 2, ...m} and (b− a) > 1, the least upper bound of the number of memory pathways is ∏m i=1 λi+( ∏a i=1 λi · ∏m i=b λi)
Proof. Denote the mMND network with skip-layer connection between layer la and layer lb as P . Denote the set of all neurons in P as
SP = {n11, n12, n13, ..., n21, n22, n23, ...} (32)
where n1i is a neuron in layer l1, and n 2 i is a neuron in layer l2, etc. The activation state of a neuron onji can be denoted with binary values 0 and 1 with onji = 1 representing that n j i receives input frequency that is above its cutoff frequency fn j i c .
The set of all possible neuron activation states O in P that generates non-zero network output feature vector can be partitioned into two subsets denoted as OA and OB .
The set OA contains all states where the input frequency fkin to any layer lk such that a < k < b satisfies
fkin < f nki c ∀ i ∈ {1, 2, 3, ..., λk} (33)
The set OB contains all the remaining neuron activation states in O, where all layers receive input frequency higher than cutoff frequency of at least one neuron in each layer.
For all the states in OA, no spike signal is sent from layer b − 1 to layer b, since at least one layer between la and lb generates no output. Hence, output from P is not affected if connections between layer li and li+1, such that i ∈ {a, a+1, ..., b− 1}, are removed. Network P is therefore equivalent to network P ′ that has layers {l1, l2, ..., la, lb, lb+1, ..., lm} connected sequentially. According to Lemma 4, it can be derived that the least upper bound of the number of distinct memory pathways in P ′ is ∏a i=1 λi · ∏m i=b λi.
Hence, for all states in set OA, the least upper bound of the number of distinct memory pathways in P is also ∏a i=1 λi · ∏m i=b λi. For all states in set OB , since the activation of neurons in the source layer of the skip-layer connection is already accounted for when considering layer la, the least upper bound of the number of distinct memory pathways is the same as network P that has no skip-layer connection, which is ∏m i=1 λi according to Lemma 4.
For the set of memory pathways from states in OA, denoted as MA, and the set of memory pathways from states in OB , denoted as MB , the number of memory pathways of network P is
K = |MA ∪MB | (34)
Since
|MA ∪MB | ≤ m∏ i=1 λi + ( a∏ i=1 λi · m∏ i=b λi) (35)
,
K ≤ m∏ i=1 λi + ( a∏ i=1 λi · m∏ i=b λi) (36)
. From Lemma 5, the bound is tight for |MA| ≤ ∏a i=1 λi · ∏m i=b λi and for |MB | ≤ ∏m
i=1 λi. It also satisfies that MA ∩MB = ∅, since all elements in MA have (m− (b− a− 1)) layers, and all elements in MB have m layers. Hence the bound is tight for (36). The proof is complete.

heading: E TIME-VARYING FUNCTION APPROXIMATION
text: A time-varying function f(t) can be approximated using piece-wise constant function, such as illustrated in Figure 4. It is therefore possible to approximate the time-varying function with a feedforward SNN by approximating each of the constant function with a memory pathway. In this section,
we test the approximation capability of feedforward SNN for time-varying functions using this principle. The target functions to approximate has the form of:
f(x) = xn
x−m (37)
Here, x is variable; n and m are function parameters. For discrete-time simulation of the network, we approximate the target function with
x = tin and f(x) = tout (38)
where tin is the input spike period and tout is the output spike period. We test approximation performance of a small-scale feedforward SNN with 6 fully-connected layers, skip-layer connections {(2, 5), (3, 5)} and 4 neuron dynamics. The network is trained with BPTT to minimize mean squared error (MSE) loss between the spike period of network output t′out and target spike period tout.
To evaluate the network’s approximation performance, we construct 6 networks with the same structure as discussed above, and change each to have different trainable parameter numbers by scaling the layer size. For baseline comparison, networks with 6 layers, no skip-layer connection and using homogeneous neurons, configured to have the same trainable parameter numbers as the proposed
networks, are tested. All networks are trained with BPTT. The loss function measures the difference between the network output spike trains and the target spike trains with MSE. Two sets of function parameters: (m = 1, n = 3.3) and (m = 2, n = 2.1) are tested for the target functions f(x) on domain [3, 10]. The resulting MSE loss for different network scales are shown at the bottom of Figure 5. It can be observed that the proposed networks can approximate target functions with less error than the baseline networks at all network scales. The smallest tested networks have relatively high losses while performance increase quickly with more trainable parameters. The rate of performance improvement decreases when trainable parameter numbers is above 4000.
To understand the impact of the target function parameters to approximation performance of SNN, we test the baseline and proposed network with 4167 trainable parameters for different pairs of function parameters m and n. The resulting MSE loss is shown with a heat map in Figure 6. It can be observed that for all m and n value pairs, the proposed network can achieve lower loss than the baseline network. Another observation is that, there is no clear correlation between the value of m and approximation error. On the other hand, for higher values of n, the approximation error generally increases for both baseline and the proposed networks.

heading: F NETWORK SPIKING ACTIVITY
text: To investigate neuron spiking activity in the proposed network, for the function approximation task as described in Appendix E, with parameters (m = 1, n = 2.3) and a randomly selected approximation point f(x = 19), we plot the timing of neuron spike at the last layer over training epochs in Figure 7 (a). It can be observed that the network initially generates spikes at widely distributed timings. After the first 50 training epochs spike timing starts to converge and remains stable at around 200 epochs. The final output spike period is t′out = 48, which matches the target output period.
Next, we investigate spiking activity of the BPTT trained network as described in Section 5.4 for N-Caltech101 classification task. Three different test data points are presented to the network, and the spikes from neurons in the first depth of layer 8 and neurons in the final layer, are recorded and shown in Figure 7 (b). Here, (i), (iii) and (v) are from layer 8; (ii), (iv) and (vi) are from the final layer. (i) and (ii) are from the observation of a data point with class label “5”; (iii) and (iv) are from the observation of another data point also with class label “5”; (iii) and (iv) are from the observation of a data point with class label “1”.
It can be observed that neurons in layer 8 exhibit similar activity in (i) and (iii) as the network is presented with data points from the same class. (ii) and (iv) show that most spiking activity is from the neuron with index 5, indicating correct prediction for the two data points. When a data point from a different class is presented, spiking activity in (v) shows different patterns than those in (i) and (iii). In the final layer, neuron with index 1 generates the most spike, leading to correct prediction. However, for this data point, neuron with index 58 is also generating a considerable number of
spikes, indicating that the network is more likely to mis-classify this data point, compared to the other two tested data points.
G IMPLEMENTATION OF HETEROGENEOUS CONV-SNN
Multi-neuron-dynamic (MND) networks with convolutional layers can be considered a scaled-up version of the mMND network, where each neuron dynamic now contains a matrix of neurons that receive input features from different spatial locations. To implement heterogeneous Conv-SNN, first consider a regular Conv-SNN layer with no heterogeneity. The spiking neuron matrix has dimension {W,H,D}, where D is the depth of the layer. The convolution filter has dimension {C,w, h,D} where C is the number of input channels and D is the number of output channels.
Based on this, a heterogeneous layer l can be constructed, by concatenating heterogeneous neuron matrices with the same W and H along layer depth. The resulting spiking neuron matrix has dimension {W,H, λD}, where λ is the number of neuron dynamics. Hence, layer depth i ∈ {k + 1, k + 2, k + 3, ..., k +D} have the same neuron parameters, where 0 ≤ k ≤ λ− 1 is an integer representing the index of neuron dynamics. The convolution filter has the same dimension as the regular Conv-SNN layer: {C,w, h,D}, which is shared by layer depth with different k values. During forward pass of this layer, a convolution operation is applied to generate an input signal matrix with dimension {W,H,D}. For each value of k, neurons in depth {k+1, k+2, k+3, ..., k+d} are simulated based on the neuron parameters for such neuron dynamic index k using the signal matrix as input. For the next convolutional layer l′ receiving input from layer l, the filter dimension is {λD,w′, h′, D′} with λD input channels.

heading: H DETAILS ON THE BAYESIAN OPTIMIZATION PROCESS
text: During the first stage of the dual-search-space optimization process, the parameters to optimize include: Nlayer, Lstart, Lend, Nskip, Ndynamic, all of which are positive integers. Specifically, Nlayer is the number of convolutional layers. For skip-layer connection, there are three configuration parameters to optimize: starting layer Lstart, which is the source layer of the first skip-layer connection; ending layer Lend, which is the target layer of the last skip-layer connection; skip-layer connection number Nskip, which defines how many connections to implement. The source layer of the Nskip skip-layer connections are placed evenly between Lstart and Lend, each with skip length equal to
⌊Lend − Lstart Nskip
⌋ .
In the case of
Lend − Lstart Nskip ̸= ⌊Lend − Lstart Nskip
⌋ ,
the value of Lend is reduced to the maximum value that satisfies
Lend − Lstart Nskip = ⌊Lend − Lstart Nskip
⌋ .
For heterogeneity, the number of different dynamic Ndynamic in all layers are optimized jointly. The constraints for the parameters are defined as:
Nlayer ∈ [4, 15] 2 ≤ Lstart < (Nlayer − 1)
(Lstart + 1) < Lend ≤ Nlayer
0 ≤ Nskip ≤ (Lend − Lstart)/2
and
Ndynamic ∈ [1, 10] .
The manually configured neuron parameters are, τm = 100 and Rm = 300 for all neuron dynamics, and a ∈ [−30,−5] is distributed evenly for each neuron dynamics. Due to the exponential increase of search space with the number of neuron dynamics in each layer, it is highly inefficient to search for every neuron parameters in each dynamic. In the second stage of the optimization process, we choose to apply Bayesian optimization for the parameter a of each neuron dynamic separately, while τm and Rm are optimized jointly with the same values shared by all neuron dynamics. a values are taken to the precision of 100, and τm and Rm values are taken to the precision of 101. The constraints are a ∈ [−30,−5], τm ∈ [50, 200] and Rm ∈ [200, 400]. The value of tnd for all networks is set to 1. The parameters of each optimized networks are shown in Table 1. Note, the skip-layer connections are listed as source and target layer pairs.

Review summaries:
In this publication, the authors present an interesting framework for approximation the mapping between sequences by feedforward spiking neural networks (SNN), providing insights on the computational capabilities of feedforward SNNs for the approximation of the mapping of input to output spike trains. And how these are influenced by hyperparameters such as the network architecture, heterogenic properties of neurons and skipping connections. The authors also provide an ansatz how these hyperparameters can be optimised in a two-steps Bayesian optimisation fashion. The performance of the approach is demonstrated by several spatiotemporal classification tasks, on datasets like DVS Gesture, N-Caltech 101 and sequential MNIST. 

This paper presents some theoretical understanding of sequence approximation using feedforward SNNs. The main conclusions are two folds: (1) a feedforward SNN with one neuron per layer and skip-layer connections can approximate the mapping rate-coding function; (2) a feedforward SNN constructed by heterogeneous neurons with varying dynamics and skip-connections can improve sequence approximation. Besides, the authors provide the DSBO to optimize the architecture and parameters of their proposed SNNs.  

The authors introduce a theory on how to learn arbitrary input spike train - output spike train mappings, using a feedforward SNN with heterogeneous neurons and skip connections. Then this theory is used to train a deep SNN using either BPTT or STDP. The SNN is then evaluated on IBM DVS Gesture, N-Caltech101. The accuracy that they get with BPTT is beyond the SOTA.