Accepted Manuscript Synaptic dynamics: Linear model and adaptation algorithm Ali Yousefi, Alireza A. Dibazar, Theodore W. Berger PII: DOI: Reference:
S0893-6080(14)00085-9 http://dx.doi.org/10.1016/j.neunet.2014.04.001 NN 3318
To appear in:
Neural Networks
Received date: 7 January 2014 Revised date: 10 April 2014 Accepted date: 20 April 2014 Please cite this article as: Yousefi, A., Dibazar, A. A., & Berger, T. W. Synaptic dynamics: Linear model and adaptation algorithm. Neural Networks (2014), http://dx.doi.org/10.1016/j.neunet.2014.04.001 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
*Title Page (With all author details listed)
SYNAPTIC DYNAMICS: LINEAR MODEL AND ADAPTATION ALGORITHM Authors: 1. Ali Yousefi PhD Candidate, Department of Electrical Engineering, University of Southern California Address: DRB 140, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9359 Email: [email protected] 2. Alireza A. Dibazar Research Assistant Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 384, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9315 Email: [email protected] 3. Theodore W. Berger Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 140, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9315 Email: [email protected]
Corresponding Author: Alireza A. Dibazar Research Assistant Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 384, 1042 Downey Way, Los Angeles, California, 90089 Phone: 213-740-9315 Email: [email protected]
Address: Alireza A. Dibazar 1042 Downey Way, DRB 140, Los Angeles, California, 90089
*Manuscript Click here to view linked References
SYNAPTIC DYNAMICS: LINEAR MODEL AND ADAPTATION ALGORITHM Ali Yousefi, Alireza A. Dibazar, Theodore W. Berger Abstract—In this research, temporal processing in brain neural circuitries is addressed by a dynamic model of synaptic connections in which the synapse model accounts for both pre- and post-synaptic processes determining its temporal dynamics and strength. Neurons, which are excited by the post-synaptic potentials of hundred of the synapses, build the computational engine capable of processing dynamic neural stimuli. Temporal dynamics in neural models with dynamic synapses will be analyzed, and learning algorithms for synaptic adaptation of neural networks with hundreds of synaptic connections are proposed. The paper starts by introducing a linear approximate model for the temporal dynamics of synaptic transmission. The proposed linear model substantially simplifies the analysis and training of spiking neural networks. Furthermore, it is capable of replicating the synaptic response of the non-linear facilitation-depression model with an accuracy better than 92.5%. In the second part of the paper, a supervised spike-in-spike-out learning rule for synaptic adaptation in dynamic synapse neural networks (DSNN) is proposed. The proposed learning rule is a biologically plausible process, and it is capable of simultaneously adjusting both pre- and post-synaptic components of individual synapses. The last section of the paper starts with presenting the rigorous analysis of the learning algorithm in a system identification task with hundreds of synaptic connections which confirms the learning algorithm's accuracy, repeatability and scalability. The DSNN is utilized to predict spiking activity of cortical neurons and pattern recognition tasks. The DSNN model is demonstrated to be a generative model capable of producing different cortical neuron spiking patterns and CA1 Pyramidal neurons recordings. A single-layer DSNN classifier on a benchmark pattern recognition task outperforms a 2-Layer Neural Network and GMM classifiers while having fewer numbers of free parameters and decides with a shorter observation of data. DSNN performance in the benchmark pattern recognition problem shows 96.7% accuracy in classifying three classes of spiking activity. I.
INTRODUCTION
Synaptic transmission is a dynamic process; it is regulated by a variety of short-term and long-term processes in which interplay of the processes defines the synaptic response (Liaw & Berger, 1996; Liaw & Berger, 1999; Dittman, Kreitzer, & Regehr, 2000; Zucker & Regehr, 2002; Tsodyks, Pawelzik, & Markram, 1998; Abbott, Varela, Sen, & Nelson, 1997; Mongillo, Barak, & Tsodyks, 2008; Maass & Zador, 1999). The synaptic response can change on a short time scale by an order of magnitude depending on the past pre-synaptic activity where its temporal dynamics forms the basis of neural computation (Liaw & Berger, 1996; Maass & Zador, 1999). Synaptic facilitation and depression are the core components of any mathematical model describing neurotransmitter release/recovery kinetics for the chemical synapses. In fact, synaptic models are mainly implemented through static weights and axonal delays failing to address fine spatio-temporal dynamics observed in spiking neural networks (SNNs). In this research, the facilitation depression (FD) model of the synaptic transmission is modified, and an approximate linear state space representation of the synaptic transmission is introduced. The synaptic response is modeled by a dynamic process addressing both STP and LTP characteristics where the mutual regulation of STP and LTP factors shape synapse temporal dynamics. There is biological evidence that in almost any neural system which requires a high processing-speed, the timing of a single spike is remarkably precise and reliable (Bohte, 2003; Mainen & Sejnowski, 1995; Berry, Warland, & Meister, 1997). Using simulation analysis, it was shown that neural spike trains have a fine temporal structure with up to one millisecond resolution in which the reliability of spike train is greater than the Renewal process (Hatsopoulos, Geman, Amarasingham, & Bienenstock, 2003; Tiesinga, Fellous, & Sejnowski, 2002). From the collected evidence, there is a well-established fact that precise and reliable spike time coding is a key component for neural processing and information conveyance in the brain. All of the aforementioned facts are in accordance with the
1
synaptic temporal dynamics; synapses are provoked by spiking activity and builds a spike-domain processing unit. In this research, a new family of artificial neural models with hundreds of synaptic connections is built which has the necessary components of a fast spike-domain computational engine. The synaptic plasticity has been accepted for over 60 years. Synapses with repeated use demonstrate different forms of plasticity. At some synapses the plasticity is in the form of facilitation dominance, and in others, the result is a decrement in synaptic strength and depression. Multiple short-term and long-term processes are present in synaptic plasticity, and its efficacy emerges by a combination of facilitation and depression processes which is highly dependent upon the timing of incoming action potentials (AP). In other words, synaptic plasticity is not merely the synaptic strength regulation; it is a dynamic process built upon multiple processes occurring in both the pre- and post-synaptic terminals of a synapse. Biologically realistic learning algorithms have emerged in parallel with experimental evidence of synaptic adaptation to address the learning mechanism in synapses and eventually artificial models of neural circuitries. STDP and its variants introduce the spike-time formulation of the Hebbian learning rule which are biologically plausible frameworks of synaptic adaptation (Song, Miller, & Abbott, 2000; Strain, Daid, Maguire, & Ginnity, 2006; Pfister, Toyoizumi, Barber & Gerstner, 2006; Bohte & Mozer, 2004; Masquelier, Guyonneau, & Thorpe, 2008; Senn, Markram, & Tsodyks, 2000; Morrison, Diesmann, & Gerstner, 2008; Izhikevich &. Desai, 2003; Bush, Philippides, Husbands, & Oshea, 2010; Burkitt, Meffin, & Grayden, 2004). Additionally, machine learning techniques including supervised (Bohte, 2003; Bohte, Poutr, & Kok, 2002; Boij, 2004; Davis, Erdogmus, Rao, & Principe, 2003; Gutig & Sompolinsky, 2006; Ghosh-Dastidar & Adeli, 2009; Ponulak & Kasinski, 2010; Fang, Wang, & He, 2010; Rostro, Cessac, Vasquez, & Vieville, 2009) and statistical (Fiete, 2003; Xie & Seung, 2004;Fiete & Seung, 2006; Williams, 1992; Legenstein, Pecevski, & Maass, 2008; Lee & Kwon, 2008; Deco & Schurmann,1999; Hinton, 2006) approaches are widely utilized for synaptic adaptation. Combinatory techniques plus heuristic search have also been extensively applied in learning of complex dynamic models (Belatreche, Maguire, & McGinnity, 2007; Jin, Wen, & Sendhoff, 2007). Reinforcement learning is the other biologically plausible synaptic adaptation framework which is categorized as a statistical learning approach. Except for one supervised learning algorithm which performs synaptic facilitation adaptation (Storck, Jakel, & Deco, 2001), other supervised learning algorithms, proposed thus far, fail to adjust the temporal dynamics of synapses. Though the learning rule proposed by Storck et al. demonstrates facilitation adaptation in DSNNs, its capability is limited to generating a single spike in the output neuron. Gupta & Long (Gupta & Long, 2009) argue that the STDP algorithm is the optimal learning approach for SNNs, even though it has worked well on a number of applications. One issue with STDP is noncausality of the learning process; when a post-synaptic neuron fires in a time-marching code, it is unknown whether one of the pre-synaptic neurons will fire in the near future. It is also unclear how a STDP rule operating in a time scale of milliseconds can justify most of the behavioral/learning events occurring at much higher time scale. Evidence for triplet and quadruplet STDP rules implies that paired-STDP is not a generalized learning rule for synaptic adaptation of neural circuitries. Statistical learning algorithms and specifically reinforcement learning suffer from slow convergence in large size networks. In reinforcement learning and other sampling methods, synaptic parameters of the neural system are updated with a noisy estimate of the expected reward's gradient. Variance of this estimate diverges with the increase in the network size. In fact, the correlation of reward with the fluctuations in spiking activity of a single neuron becomes weaker as the network size grows. Statistical learning maximizes the probability of a specific activity pattern rather than its precise timing. This is in contrast with the biological finding of precise and reliable spike timing. Heuristic search algorithms fail to demonstrate a learning mechanism in neural circuitries. Moreover, the processing time in finding optimum synaptic parameters increases exponentially as the network size and processing time increases.
2
In this research, a family of learning algorithms is introduced satisfying neuro-biological constraints including simultaneous pre- and post-synaptic adaptation and precise spike timing. A supervised learning algorithm capable of simultaneously adjusting STP and LTP factors of individual synapse of a DSNN is introduced to perform a spike-in-spike-out functional mapping (Yousefi, Dibazar, & Berger, 2011). The assessment of neural circuitries has generated considerable research interest. One of the frontiers of neuroscience research focuses on replacing a damaged human brain region with a neural prosthesis in which the prosthesis duplicates the functionality of the damaged region. Synaptic plasticity is the core component of any cortical neural circuitries. Reproducing the neural functionality requires a precise artificial model capable of performing the spike-domain processing task. Though functional approximation and detailed compartmental methods have been applied in reproducing neural region functionality, the learning scalability and computational complexity of both approaches limit their practicality in the prosthesis application. The Volterra series and its simplified kernel counterpart encounter an exponential growth of the model parameters in accurate modeling of the neural functionality. Though detailed conductance models demonstrate kinetics of different compartments of a neural system, these models are exhaustive and require large computational resources. Similar to the functional approximation technique, free parameters of the detailed model grow rapidly making their tuning process a complicated task. In this research, an intermediary technique for modeling the dynamics of cortical neural circuitries is proposed; the phenomenological models show a balance of biological feasibility and computation complexity. Moreover, structure of the model has the flexibility of expanding in both synaptic connections and neuron layers making the model a superior technique in reproducing the neural system functionality. It will be shown that the phenomenological model of the manuscript along with the proposed learning algorithm are able to predict the spiking activity of Hippocampal neurons significantly better than non-parametric approaches while being able to explain the possible learning mechanisms that occur in the corresponding biological system. The efficiency of the proposed model and training algorithm will be demonstrated in a dynamic pattern recognition task and its performance will be compared with a 2-layer sigmoid neural network. II.
LINEAR MODEL OF SYNAPSE DYNAMICS
The facilitation depression (FD) model of synaptic response (Tsodyks, Pawelzik, & Markram, 1998; Dittman, Kreitzer, & Regehr, 2000; Mongillo, Barak, & Tsodyks, 2008) is described by: ⁄ ⁄ 1 ∗ ∆ ∗
(1.a)
⁄ 1 ⁄ ∗ ∗
(1.b)
⁄ ⁄ ∗ ∗ ∗
(1.c)
⁄ ⁄ ∗ ∗ ∗
(1.d)
(1.e)
Variable describes facilitation dynamics in response to incoming APs, and variable represents the portion of release ready vesicles. Auxiliary variables and model an -synapse function representing the PSP - - elicited due to synaptic neurotransmitter release (Gerstner, 2006). The PSP of a single synapse in equation set (1) is a non-linear function of synapse state variables - ∗ -; and it is the interplay of facilitation and vesicle recovery/release that determines the temporal dynamics of synapse. The ∗ term defines the portion of vesicle sites that release their neurotransmitter due to an impinging AP time while parameter defines the population of neurotransmitter receptor sites. Parameter ∆ controls calcium influx in the presynaptic terminal of a synapse, and it is the main biological factor in modulating synaptic dynamics. Parameter is a post-synaptic factor representing the number of neurotransmitter receptors or simply synaptic strength. Facilitation increase and vesicle release dynamics are both AP-driven process where the AP time is defined by . The FD model proposed in equation (1) defines the fundamental biological components shaping
3
synapse temporal dynamics. Appendix A defines a more intricate synapse model including multiple pre- and post-synaptic facilitation processes in addition to the calcium dependent vesicle recovery time constant proposed by Regehr et al. The discrete counterpart state space model of synapse dynamics is built by utilizing the Euler discretization technique. The discretization resolution ∆ is defined in milliseconds with a nominal value of one millisecond. The millisecond resolution is required to accurately trace the fast kinetics of synaptic facilitation and vesicle release process. Facilitation increase and vesicle release in response to incoming AP occurs in about the one millisecond time period; thus, in the discrete model both facilitation increase and vesicle release evolve at the same time of impinging APs. Equation set (2) defines the discrete state space model of the synaptic response: !"# $1 ∆ ∗ 1 ' ( ∗ ! ∗ ∆ ∗ 1 '
% % ∆
∆
!"# $1
&
∆ %)
&
* ∗ ∆+ ∗ 1 ' ( ∗ ! *1 ∆ %&
∆ %&
(2.a)
∆+ ∗ 1 ' ∗ ! ∗ !
!"# *1 + ∗ ! ∗ * ∗ ∆+ ∗ 1 ' ∗ ! ∗ *1 ∆
%,
∆ %&
!"# -1 / ∗ ! ∗ * ∗ ∆+ ∗ 1 ' ∗ ! ∗ *1
! ! !
∆
%.
∆ %&
∆ %&
∆ %&
∆ %)
∆+ ∗ 1 ' ∗ ! ∗ !
∆+ ∗ 1 ' ∗ ! ∗ !
(2.b) (2.c) (2.d)
Variable 1 ' denotes AP-time, and it is equal to one at spike time. The state variables of equation set (2) define discrete samples of the continuous FD model state variables. Table 1 defines each parameter of the model and its nominal value. Figure (1) shows state variables of the FD model in response to a sample spike train. (2.e)
Table 1 Range of parameters in synapse FD model Parameters ∆ ∆
Facilitation time constant Resting facilitation Facilitation increment factor Vesicle recovery time constant Maximum number of release sites Rise time of PSP Rise time of PSP Processing time resolution
Value 100…500 milliseconds 0…1 0…1 100…1500 milliseconds 1...20 3…15 milliseconds 3…15 milliseconds 1 millisecond
The synapse model defined in equations (1) and (2) are a non-linear system, and it is the AP occurrence that injects non-linearity to the synaptic response (Morrison, Diesmann, & Gerstner, 2008; Yousefi, Dibazar, & Berger, 2011). Non-linearity in the synapse model appears in the ! ∗ ! term, and it becomes linear if the ! ∗ ! term is replaced by a surrogate variable. The new variable 0! ! ∗ ! is defined for the linearization process. The variable denotes the minimum vesicle release for a possible AP at time ' 1. To build the synapse linear state space model, 0!"# needs to be defined by a linear combination of the ! , ! and 0! state variables. Variable 0!"# !"# ∗ !"# can be rewritten as a linear function of ! , ! , 0! plus a new non-linear term - ! 1 ∗ ! . 0!"# !"# ∗ !"# 2# ∗ ! 21 ∗ 3# ∗ ! 31 ∗ ! ∗ ! 34
2# ∗ 3# ∗ ! ∗ ! 2# ∗ 31 ∗ ! ∗ ! ∗ ! 2# ∗ 34 ∗ ! 21 ∗ 3# ∗ ! 21 ∗ 31 ∗ ! ∗ ! 21 ∗ 34
2# ∗ 3# 21 ∗ 31 ∗ 0! 2# ∗ 31 ∗ ! ∗ ! ∗ ! 2# ∗ 34 ∗ ! 21 ∗ 3# ∗ ! 21 ∗ 34
4
(3)
0.5
100
200
300
400
0
500
1
0.5
0.5
200
300
400
0
500 G
2 0
K
100
100
200
300
400
2 0
100
200
300
400
1 0.5 0
100
200
300
400
1.5 1 0.5 0
AP
1 0.5 0
100
200 300 Time (msec)
400
500
(a)
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200 300 Time (msec)
400
500
1 0.5
500
200
2 0
500
100
2 0
500 K
G
0
N
1
Y /Y m a x
N
F
0.5 0
Y /Y m a x
∆F = 0.8, Nmax = 1.0 1
AP
F
∆F = 0.1, Nmax = 1.0 1
0
(b)
Figure 1 Synapse FD model response to a 500 millisecond spike train with 50 milliseconds ISI. a) Synapse with ∆F=0.1, b) Synapse with ∆F=0.8
Parameter 2# , 21 , 3# , 31 and 34 are defined using equation set (2). The linearization process can generate higher order terms of ! 6 ∗ ! and ! 6 , 8 1 whose values shrink to zero as increases. This is because both ! and ! vraibles are normalized positive terms with a maximum value of 1.0. The linearization process is completed by truncating higher-order non-linear terms of the synapse FD model, and then defining a set of new surrogate variables - along with the recursive equations for updating each new variable. For the first order linear model, the only non-linear term is ! ∗ ! and other terms of ! 6 , 8 1 is set to zero. In the second order linear model, non-linear terms include ! ∗ ! , ! 1 ∗ ! and ! 1 and other extra terms containing ! 6 , 8 2 is set to zero. The same procedure is applied for a higher order linear model in which any non-linear term consisting of ! 6 , 8 : is set to zero; : defines the order of linear model. Figure (2) shows the normalized error of vesicle release for different linearization orders - : ∈ - within the possible range of interspike interval (ISI) and ∆. Independent of linearization order and facilitation factor, the maximum vesicle release error increases with smaller ISIs. The vesicle release error decreases as the linearization order increases. Meanwhile the maximum of vesicle error is delayed to later APs of the impinging spike train. Similar to the ISI scenario, the normalized error gets larger for lower facilitation factors. Synapses with a lower facilitation factor correspond to potentiating synapses. Therefore linearization error accumulates through multiple impinging APs. At lower ISIs, the synapse vesicle pool is depleted completely after ten to fifteen
5
impinging spikes due to the lack of recovery reco time between incoming APs. Thus, hus, the maximum error is more reasonable until til a specific number of consecutive APs. The simulation results show that the maximum vesicle release error for the linear model with a third-order order approximation is below 15% following a prolonged spiking activity. For spiking activity with a maximum of ten consecutive impinging spikes, the maximum vesicle release error drops below 7.5%. The simulation resultss show that for spiking activity with ISIs of more than 100 milliseconds, the first-order first order linear approximate is an accurate model for synapse dynamics. The results result confirm that the third-order order linear approximate state space model is an accurate model of synapse synapse dynamics independent of ISI variation and the facilitation factor. In the third-order order linear model, other surrogate variables of the model - besides the 0! variable- are: ?! ! 1 ∗ ! , @! ! 1 , P! ! 4 ∗ ! and H! ! 4 . Higher igher order term terms of ! 6 ∗ ! and ! 6 , 8 3 are set to zero. Equation (4)) presents synapse linear state space model: ∗ D A!"# BC 1 ' ∗ B ∗ A! DC 1 '
! E0 0 0 0 0 0 0 1 1G ∗ A! Q ∗ A!
(4.a)
A! E! @! ! 0! ?! H! P! ! ! is the model state vector, and ! corresponds to the PSP. Elements of BC , B , DC and D are constant functions of synapse internal parameters, and 1 ' is one at AP time. Equation (4) ( introduces a generalized linear model for synaptic response response. The model also provides a compact tool for temporal and steady state analysis of the synaptic response. The following subsections analyze synapse temporal dynamics as a function of synapse internal parameters and AP variability. GR
(a)
(4.b)
(b)
Figure 2 Maximum vesicle release error in the first, second and third order linear approximate model of synapse FD model model. The error is equal to normalized vesicle release between the linear and non-linear non models defined by OOOOOOOOOOOOO IJJ ∗ |LM"I ∗ NM L M"I ∗ NM |⁄LM"I ∗ NM ; in which LM"I ∗ NM determines vesicle release at impinging spike time, and OOOOOOOOOOOOO L M"I ∗ NM is the linear estimate of vesicle release using the approximate linear model. a) Maximum vesicle error for different ISIs. For each specific ISI, the parameter space - facilitation factor and AP index - is searched for the maximum error, b) Maximum vesicle error for different facilitation tation factors. factor Number of APs is limited to 10 consecutive spikes. Synapse time constants - ST , SU - are set at IVJ, WVJ milliseconds.
6
A. Parameter Sensitivity Analysis
Each element of the state matrix - BC , B - and input vector - DC , D - is a differentiable function of the synapse internall parameters. Equation (5) (5 denotes the derivative of a state vector relative to any of the synapse parameters: XYZ[ X\
X`Z X\
-
X]^
Q∗
X\
1 ' ∗
XYZ X\
X]_ X\
/ ∗ A! BC 1 ' ∗ B ∗
XYZ X\
-
Xf^ X\
1 ' ∗
Xf_ X\
/
(5.a) (5.b)
(b)
(a)
Figure 3 Synapse temporal dynamics. a) Synaptic PSP response for ∆L J. JV, J. h and J. i plus PSP sensitivity in response to ∆L changes. Input APs has ISI of 50 milliseconds illiseconds and last for 400 milliseconds. Smaller ∆L leads to a potentiating synapse. Higher values of ∆L turns to a depressing synapse, b) Synaptic temporal selectivity in response to ISI variation for different ∆L. Smaller ∆L has a higher sensitivity for lower ISIs, ISI and the peak of synaptic efficacy move to lower ISI for higher values of ∆L. Synapse time constants - ST , SU - are set at IVJ, WVJ milliseconds.
Variable a can be any synapse parameter; in figure (3.a) .a) the PSP derivative with respect to facilitation factor - ∆ - at different values is shown. The PSP derivative with respect to ∆ is not a monotonic function, and its sign and amplitude changes depending on AP time and internal synapse parameters. The PSP derivative with respect to is equal to a normalized PSP. Thus, the PSP increases with the growth of . Equation (5)) will be used in the learning algorithm; the equation defines a recursive function to calculate the derivative of the neuron membrane potential relative to the synapse parameters. B. AP Pattern Sensitivity Analysis
Synaptic efficacy is defined by the ratio rat of vesicle release in the jk and the first impinging AP time. This ratio for different ISIs builds a quantitative measure of the synaptic temporal dynamic dynamics in response to the variability of the impinging spike train. Equation (6) denotes the efficacy ratio for a spike train with b-millisecond ISIs. Figure (3.b) shows the ratio for different ISIs and different facilitation factors. c A∗ "# 4 ⁄A "# 4 1 1 j ∗ b 1
(6)
AO!"# BC e ∗ B ∗ AO! DC e ∗ D
(7.a)
Synaptic efficacy is one possible measure to quantify synaptic temporal dynamics. s. It is also possible to identify synapse temporal dynamics using its maximum release lease in response to impinging spikes spikes. Similarly, synaptic transmission in response to the mean firing rate of incoming incoming spike trains can be defined using the linear model.. For a homogenous firing rate - e -, equation set (7) represents the synaptic response and the state variables' derivatives with respect to change in firing rate at its steady state.
7
XYOZ[ Xl
BC e ∗ B ∗
XYOZ Xl
B ∗ AO! D
(7.b)
AOmnopqr s BC e ∗ B t# ∗ DC e ∗ D
(7.c)
uOmnopqr Q ∗ AOmnopqr
(7.d)
Temporal dynamics of the synaptic response is addressed in equations (4-7). The next section extends the linear state representation to a neuron model consisting of multiple synaptic connections. C. Neuron Membrane Dynamics Model The neuron model integrates the PSP response of multiple synaptic connections plus the refractory process emerging from its previous firing activity. The neuron model is capable of receiving external stimuli in which both the external stimuli and refractory process are defined by a linear state space representation. Equation set (8) defines neuron membrane potential dynamics and its firing activity criterion: ' ∑ ' Ref z ' {o '
(8.a)
{| ' 1 1 ∆⁄ ∗ {| ' ∗ ∆ ∗ s| '
(8.b)
Ref z ' 1 1 ∆⁄ | ∗ Ref z ' Refp}~ ∗ 1 '
If ' 8 k ∩ ' ' 1 ⇒
1 '
(8.c)
1
(8.d)
Variable ' is neuron membrane potential at time ', and variable { | ' is the impact of external stimuli - s | ' - on the neuron membrane potential. Variable ' is the PSP of the jk synapse in response to incoming APs, and variable Ref z ' corresponds to the refractory effect in response to neuron spiking activity - 1 ' . Equation (8.d) defines the firing criteria of the neuron; the neuron fires an AP whenever its membrane potential passes the firing threshold. Parameter and are defined by neuron membrane characteristics, and | and Refp}~ - a negative term - model refractory dynamics. The neuron membrane dynamics follows the same generalized form proposed for synaptic response by equation (4) in which 1 ' is replaced by a matrix replicating all parent’s plus neuron's spiking activity at time '. Figure (4) shows the motif of neuron models with multiple synaptic connections. For a Type I neuron model, the network consists of feed-forward synapses projecting input spiking activities on the neuron membrane while the Type II neuron has an extra synapse projecting its firing activity to itself. The feedback synapse has different dynamics compared to the refractory process. The feedback synapse can be either excitatory or inhibitory with a dynamic response. Similarly, feed-forward synaptic connections can be excitatory or inhibitory. The dynamic synapse neural model and its morphology Type I and Type II - are defined; the next section will focus on the learning algorithm for the DSNN. Ref
Ref
(a)
(b)
Figure 4 Type I and Type II DSNN. a) Type I DSNN with multiple feed-forward synaptic connections, b) Type II DSNN with multiple feed-forward synaptic connections plus a feed-back synaptic connection. Feed-back synaptic connection has different dynamics comparing to the refractory process, and it can be either excitatory or inhibitory.
III.
SYNAPTIC ADAPTATION ALGORITHM
A biologically plausible learning rule capable of simultaneously adjusting long term plasticity - LTP and short term plasticity - STP - parameters of individual synaptic connections of a Type I/II DSNN is introduced. Facilitation factor - ∆ - and synaptic strength factor - - are adapted simultaneously
8
to shape desired neuron spiking activity. The learning rule is a supervised learning rule in which the similarity between a neuron’s spiking activity and a desired spike train is increased. Moreover, a new STDP learning rule capable of training Type I/II DSNNs with hundreds of synaptic connections is introduced. In the following sections, a similarity measure between spike trains and supervised learning rules are introduced. The new STDP and its formulation will be described in detail. A. Spike Train Similarity Measure The spike train similarity measure quantifies the coincidence between APs in the test and the desired spike trains (Kreuz, Haas, Morelli, Abarbanel, & Politi, 2007; Yousefi, Dibazar, & Berger, 2011). The proposed similarity measure quantifies the similarity between two spike trains, the test and desired spike, in two different time scales. The fine temporal scale identifies each AP of the test spike train as a: 1) similar, 2) missing, or 3) extra AP. Meanwhile, the global scale measure returns a normalized value determining the overall similarity between two spike trains. An AP is considered similar if it happens in the time vicinity of one of the APs in the desired spike train. This vicinity time interval is called the similarity range and determines acceptable jitter between spike pairs in the test and desired spike trains. It is defined to be less than half of refractory period and is generally limited to a couple of milliseconds. The overall similarity of two spike trains is defined by: C ⁄j 1, ,
C is the number of similar spikes; and denote the number of spikes in the desired and test spike trains respectively. The similarity measure - - returns a positive normalized value with a maximum of one corresponding to completely similar spike trains. The measure can be defined at different time scales, and as the time scale increases temporal resolution of the measure drops. Figure (5) shows similar, missing and extra spikes for a sample pair of test and desired spike trains in which gets value of 0.6. Similarity between two spike trains is increased by reducing the number of missing and extra spikes. To increase the similarity, the membrane potential at extra spikes needs to be dropped below firing threshold, and at the same time, the membrane potential of missing spikes needs to be increased above the firing threshold. The learning algorithm adjusts membrane potential for increasing the overall similarity between two spike trains. The following section defines the learning algorithm in DSNN. (9)
Target B'
Test D'
s
m
s
e
m
s
'
j
AP Refractory Period Similarity Range Similar Extra Missing
Figure 5 Similarity measure. Similarity measure classifies APs of a test spike train to three classes of : 1) Similar (s), 2) missing (m), and 3) extra (e) spikes. ' is the relaxed spike train defined by the comparison between B' and D' spike trains. The S value is equal to 0.6.
B. Off-line Learning Algorithm The learning algorithm modifies synaptic parameters to maximize the similarity measure between the training and desired spike trains. AP firing is a discontinuous function of synapse membrane potential, and it is approximated with a continuous function to build a tractable learning algorithm. The learning objective function approximates the similarity measure; it is defined by: E ∑ 1qp~ n ∗ Ep~ n 1 1qp~ n
∗ Eztp~ n
Ep~ n log -1 exp -α ∗ Vn Vn // log *1 exp -β ∗ Vn 1 Vn /+
Eztp~ n log1 expα ∗ Vn Vn ∗ log *1 exp -β ∗ Vn Vn 1 /+
9
(10.a) (10.b) (10.c)
Ep~ n and Eztp~ n determine the network membrane potential deviation from its desired behavior, and 1qp~ n is either zero or one depending on the desired spiking activity, or the similarity measure output, at processing time '. Both Ep~ n and Eztp~ n return a small positive, almost zero, value whenever the training network neuron membrane potential satisfies firing/non-firing criteria. The objective function E returns a small positive value at similar spike times, and its returned value grows as the number of dissimilar spike grows. It is theoretically shown - in Appendix B - that the expected value of c is proportional to the number of missing and extra spikes. Thus, reducing the objective function value yields a higher similarity measure for the training network spiking activity. In the objective function definition, free parameters α and β are set to positive large values that penalize neuron membrane potential deviation from its desired value. The optimum direction for synaptic parameter adjustment is defined by the learning objective function gradient. Synaptic parameters are adjusted in the gradient descent direction to increase the similarity between the test and the desired spike train. The learning algorithm is defined by: a a ∗ ~ c
(11.a)
P maxLower Bound, minUpper Bound, P η ∗ ~ E
(11.b)
Parameter ∆ is a positive number with a maximum value of one, and similarly is constrained to be a positive number. Excitatory or inhibitory synapses are defined independent of the sign of . Thus parameter adjustment - equation (11.b) - is defined by considering parameters' range and sign. The learning objective function returns a value close to zero for processing bins with similar spiking activities. Thus the objective function gradient is defined by gradient of either Ep~ n at missing spikes or Eztp~ n at extra spike time. In fact, the gradient of the learning objective function is defined by the sum of sparse samples of the neuron membrane potential gradient at missing and extra spike time. It is theoretically shown - in Appendix B - that the membrane potential of Type I/II DSNNs is a unique function of its synaptic parameters. The uniqueness is extended to spiking activity in Type I/II DSNNs in which the free synaptic parameters of two DSNN networks with the same number of synaptic connections and similar spiking activity converge to each other. The gradient descent learning algorithm reduces the returned value of the learning objective function, and this reduction corresponds to increases of the similarity measure. Thus, the learning algorithm converges to optimum parameter adjustment. Algorithm I describes processing steps for the learning algorithm in Type I/II DSNNs. Algorithm I Supervised Off-line Learning Algorithm in Type I/II DSNN
1. 2. 3. 4. 5.
6. 7.
Initialize synaptic parameters of the training network with random values Set the desired similarity range - Set the learning iterations – L Set j to zero – there are M set of input and output spike pairs. ¤¥' - and desired spiking activity - B' . ' ∈ For jk pair of input - 1 ¤¥, 0 a. Set § 0, D' 0 b. If § 8 @ ∪ 8 → Jump to step 6 c. ª c 0, § § 1 ¤¥' input stimuli (Equation 8.b) d. Run the network with 1 e. Calculate and build ' - ' is built by matching B' - network spiking activity- with B'
f. Calculate Ref z using ' (Equation 8.c) - Refractory compensation g. Calculate output neuron membrane potential dynamics (Equation 8.a) h. Calculate optimum parameter update (Equation 11.a) i. Update network parameters (Equation 11.b) j. Jump to step b m m 1, if j ¬ H jump to step 5 Stop learning
10
The similarity signal - ' - generates a relaxed desired spike train emerging from the spike train similarity process. Figure (5) shows the ' signal for a sample pair of training and desired spike trains. The relaxed spike train - ' - will be utilized in refractory process compensation. For a Type II DSNN, this will also determine the impinging spike train fed into the feedback synaptic connection. The refractory process of a training network initiates at each spike of the output neuron. Any mismatch in the neuron spiking activity changes the dynamics of the refractory process. Using the missing and extra spikes of similarity signal, the refractory process of a training network is compensated to replicate the desired refractory process. The feedback synaptic connection of the Type II DSNN can be broken by replacing its impinging spiking activity with the desired spiking activity. Replacing '
with the desired spiking activity - B' - makes the learning algorithm insensitive to spiking jitter, and the learning convergence is boosted as both of the desired spike train and synapse parameters are adjusted at each learning step. Practically in Type II DSNNs, the learning algorithm is repeated multiple times for each training sample. Although Algorithm I is an off-line learning algorithm, the optimum parameter adjustment is defined by observing only a sparse set of spiking activity or a partial observation of the neuron membrane potential. Regularity in the parameter update rule at missing and extra spike time suggests the possibility of building a more generalized learning rule which is online and is independent of observing the desired spiking activity. The next section will introduce a gradient-gain learning rule which is an on-line learning algorithm applicable in the synaptic adaptation of Type I/II DSNNs. C. Gradient-Gain Learning Rule The gradient of the post-synaptic membrane potential with respect to synapse parameters is positively gained at missing spike times and negatively gained at extra spike times. The gain is either a constant value, or it reflects the neuron membrane potential deviation from the firing threshold. The gain can be generalized to any function of membrane potential and its firing threshold - including a constant with alternating sign -, and this gain function - defined by ' - adapts the spiking activity of the training network. In the Algorithm I, the gain value during similar spiking activity and any silent period is set to zero, and it is equal to a constant at missing and extra spike time. It is assumed that all synaptic connections are excitatory, and desired spiking activity is generated by the same neural topology. Refractory compensation as defined in the algorithm I will be dropped if the gain function predicts spiking activity - in the next milliseconds where defines the similarity range - of the training network rather than its instantaneous spiking activity. For example, the training algorithm expects to see a spike in next milliseconds and the training network activity is impotent to fire an AP shortly. Thus, a positive gain is injected to the network boosting its membrane potential to fire an AP. If the network responds to the gain injection by firing an AP, the refractory process is already adjusted by the network response. The same scenario happens whenever the network receives a negative gain, and it reduces its membrane activity below the firing threshold. The gradient gain with a predictive behavior resembles the temporal difference - TD - learning algorithm in which the system behavior is adapted to reduce upcoming error. Algorithm II Gradient-gain Learning Algorithm in Type I/ II DSNN 1. 2. 3. 4. 5.
Initialize synaptic parameters of training network with random values Set ' 0 and ª~® 0 ( is the PSP of mth synapse, and p is the synapse free parameter) Build network state space model (Equation set 8) Run the network with ¤¥ 1' input stimuli (Equation set 8) Update network parameters with the received gain - ¯ is the forgetting factor and ' is the gain function a. ª~® ∶ ¯ ∗ ª~® ' ∗ ª~® '
p} ∶ maxLower Bound, minUpper Bound, p} η ∗ ª~® ' ' 1 Jump to step 3
b. 6. 7.
11
The generalized learning rule for a Type I/II DSNN is defined by Algorithm II. The learning rule is defined independent of refractory compensation and uses the gain function to adjust the DSNN synaptic parameters for a desired spike-domain functional mapping. The gain is any evaluative function reflecting the similarity between the training network spiking activity and a desired one. Learning algorithm II is similar to a family of approximate dynamic programming algorithms called Gradient Policy (Williams, 1992). The gain function corresponds to the reward function of the gradient policy algorithm, and the rewarding process measures the training network’s ability to perform a defined dynamic task. Similarly, the gradient-gain can be generalized to any evaluative function of spiking activity of the training network. An evaluative gain function can be defined by the instantaneous spiking activity of the training network, its average spiking activity through time, or any other combination of its temporal response. Thus, the gain function is shaping neural system response for performing a spike-domain functional task. Gradient policy learning algorithms - including REINFORCE and GOMPD (Baxter & Bartlett, 2001; Williams, 1992)- build a stochastic gradient ascent approach to adjust training network policy parameters with a guarantee of convergence. In contrast, the gradient-gain learning proposed here is a deterministic approach in which optimum parameter adjustment is determined by the variability in spiking activity of the training network. Optimum parameter adjustment of the gradient-gain learning rule is the same as in off-line learning; thus, the off-line learning convergence can be extended to this gradient-gain learning. The learning objective function is a convex function with a minimum value of zero. The learning objective function converges to zero as the network converges to its optimum parameters. Optimum parameter update is defined by two terms: 1) the local gradient of each individual synapses, and 2) the gain function. The gradient term is generally non-zero; thus, receiving a non-zero gain implies an incomplete learning process. The gain function definition requires convergence to zero as the training network learns the desired spiking activity in response to incoming stimuli. The similarity signal generated by the spike similarity process returns zero for a similar spike train; thus, the similarity signal is a viable gain function. Either defining the gain function with a stopping mechanism or dropping its time average from the instantaneous gain are possible choices to build a learning stopping criteria. The gradient-gain learning is a function of the pre- and post-synaptic spiking activity of the network, and it can be extended to a STDP learning rule. The next subsection discusses a new STDP rule for synaptic parameter adjustment in DSNNs. D. STDP Learning Rule In the gradient-gain learning, synapse parameters are updated at the time indexes - '± - at which the network receives the gain - '± . ª c'± is a function of pre-synaptic AP time and synapse internal parameters which involve the history of impinging APs for synaptic adaptation. The gain factor is a function evaluating the network post-synaptic spiking activity. Thus, the gradient-gain is a function of pre- and post-synaptic spiking activity of each individual synapse leading to a new class of STDP learning rule. In contrast to pair-wise STDP learning rules (Morrison, Diesmann, & Gerstner, 2008) in which synaptic weight is adjusted by the relative timing of a pair of pre- and post-synaptic AP time, ª c'± at gain time is determined by a history of AP times and synapse internal parameters. For the pair-wise STDP rule - plus their extension of triplet and higher order STDP rule -, synaptic weight are strengthened whenever a pre-synaptic AP precedes post-synaptic AP and weakened when the postsynaptic AP leads a pre-synaptic AP. The pair-wise STDP learning rule presents a uni-modal learning behavior. However, the STDP learning rule emerging from the gradient-gain alternates the synaptic weight strengthening/weakening mechanism - along with STP factor - depending on the gain function. The pair-wise STDP requires a matching process to find corresponding pre- and post-synaptic APs, whereas the new STDP rule is independent of an AP matching process as the gain function is defined independent of pre- and post-synaptic AP activity sequencing. In contrast to the standard STDP learning rule (Morrison, Diesmann, & Gerstner, 2008), the gradient gain STDP rule is a causal learning
12
process. The gain value can be defined efined without a matching process; therefore, the STDP learning rule emerging from gradient gain learning can be an on-line on learning rule. Figure (6)) describes the STDP learning rule defined by gradient-gain gain learning for different sets set of synapse parameters. Appendix C describes the components of the gain-STDP gain learning rule. Through this section, off-line line and on on-line learning rules for synaptic adjustment in DSNNs DSNN are introduced. The STDP plasticity learning rule addresses both STP- and LTP-based synaptic adjustment in DSNNs. The following sections will focus on the computational capability and verification of the proposed learning algorithms for DSNNs.. 1st
2nd
3rd
4th
Incoming Spike Train Gain Function
5th ∆t
Pre-synaptic neuron Post-synaptic synaptic neuron
(a)
(b) ∆ 0.05, 1
(c) ∆ 0.3, 1
(d) ∆ 0.9, 1
Figure 6 Gradient-gain gain learning rule and STDP. A positive gain corresponds to an expected AP at gain injection time, meanwhile the negative gain reflects an unexpected AP. STDP graphs shown here define STP and LTP adjustment for different facilitation factors with a fixed positive gain injected inject with different delays. ISI of impinging AP is 50 millisecondss and it includes 5 consecutive APs. a) Schematic of gradient-gain gradient learning rule, b) STP and LTP learning curve for a potentiating synapse, c) STP and LTP learning curve for a synapse with medial facilitation factor, d) STP and LTP learning curve for a depressing ressing synapse. The learning curve is a non-linear non linear function of facilitation factor and impinging AP's ISI. Although the learning curve shown here is only limited to 5 APs, practically larger number of APs are involved in synaptic adaptation. Synapse time constants - ST , SU - are set to IVJ, WVJ milliseconds.
IV.
SIMULATION RESULT
In this section, different aspects of DSNN networks and the proposed learning learning rule will be addressed through a set of benchmark problems plus replicating in vivo neural responses. First, a system identification analysis is performed on Type I/II DSNNs. Type I/II DSNNs with different number numbers of synaptic connections are trained by observing the input stimuli and the output spike train of a desired DSNN network having the same me topology. Algorithm I’s convergence and accuracy in parameter estimation in response to network size growth and similarity range variation will be studied in detail. The second set of results examines the FD model’s model capacity to predict the temporal dynamics mics of the
13
Schaffer collateral synaptic pathway which projects CA3 neuron spiking activity onto the CA1 pyramidal neurons. Free parameters of the synapse model are tuned to replicate the CA1 post-synaptic membrane potential. Furthermore, Algorithm II’s ability to estimate the optimum synaptic parameters will be demonstrated. Next, a generative model - a Type II DSNN with one feed-forward synaptic connection - with dynamic synapses is introduced to produce the different spiking and burst patterns of various cortical neurons. The gain function for the on-line learning process is defined, and capability of gradient-gain learning algorithm in generating different spiking patterns will be analyzed. For the last portion, a spatio-temporal pattern recognizer is built utilizing the DSNN model. The classifier consists of three output neurons recognizing three different classes of spiking stimuli in which each output neuron matches with one of the input classes. The network is trained by the gradient-gain learning rule to fire APs with a shorter onset delay - or a higher number of APs - for the neuron that corresponds to the class to which the incoming spike train belongs. In the system identification problem, the gain is an estimate of the membrane error between the test and target neurons. In contrast, synaptic adaptation for the generative model utilizes the spike time difference for the gain function. The gain function for the spatio-temporal pattern classifier and the CA1 neuron model switches between constant positive and negative values at the desired and undesired spiking times to regulate the firing time/rate of the output neurons. This set of problems demonstrate versatility of the gain function and how the spatio-temporal processing capability of synaptic responses can be merged to build a powerful biological computational unit. In the following subsections, the different simulation problems are discussed in detail. A. System Identification in DSNN
Algorithm I is applied for synaptic parameter estimation - Δ and - for a DSNN network by observing its input and output spiking activity. Four different properties of the proposed learning algorithm will be analyzed: 1) convergence, 2) repeatability, 3) scalability, and 4) similarity range variation. The learning properties are quantified by their convergence to unknown synaptic parameters, similarity measure improvement, and the square error in membrane potential prediction. The components and structure of simulation test are as follows: a. Both the training and target network consists of Nm excitatory synaptic connections. b. The facilitation factor for each synapse is generated independently by a uniform distribution within the range of E0.05 0.95G. c. The synaptic strength for each synapse is defined by uniform distribution in the range of E0.05,0.05 K ∗ UE0,1 G. K is dependent on the number of synaptic connections. Its minimum value is 0.2, and its nominal value is defined as 20⁄Nm . UE0,1 is a uniform distribution ranging from 0.0 to 1.0. The synaptic strength is chosen by the assumption that six synapses with the same facilitation factor of 0.5 will satisfy firing criteria at the time of their first impinging spikes. d. The firing threshold is set to be 1.0. e. The refractory process parameters τ·o¸ , Refp}~ are 2.5 msec, 5 . f. The input stimuli to the network consists of an array of Nm spike trains. Each spike train is generated independently by a Renewal process with a random rate uniformly chosen from 5 to 20 Hz. g. There are at most 1500 training spike trains and 100 test spike trains. Training progress is examined on the test data set. h. The simulation time is 400 milliseconds, and the time resolution ∆t is set to 0.2 msec. i. The nominal similarity range is 2 milliseconds unless otherwise mentioned. j. The learning rate is set to be 0.01. k. In Type II DSNNs, the learning process is repeated 5 times per each training sample. For Type I DSNNs, the learning process is applied once per training sample. l. , are set to 150,250 milliseconds.
14
1. Learning Convergence Table 2 and 3 show the learning progress for a Type I and II DSNN with 10 synaptic connections. Results in both tables show a significant improvement in the similarity measure and membrane error reduction. The similarity measure between the desired and training network grows from 0.21 - 4 similar spikes - to 0.97 - completely similar spikes - in the Type I DSNN while the membrane potential error decreases by 34.5 dB. Similarly, the same improvements are observed in the Type II DSNN. The similarity grows from 0.07 - 1 similar spike - to 0.97 - completely similar spikes - while there is more than a 31 dB improvement in the membrane potential error. The result implies an accurate estimation of the network parameters. Figure (7.a) and (7.b) show the initial guesses of the unknown parameters ∆, - and the estimated parameters after training with 1500 samples using a Type I DSNN. Figure (7.c) and (7.d) show the membrane potential of the trained and target networks in response to a sample input using the initial and trained parameters. The numerical result for the sample task verifies convergence to an optimum parameter. The next essential questions are repeatability and scalability of the learning algorithm in synaptic parameter adjustment. The following simulations will focus on Type I DSNN as the learning convergence to optimum parameter in Type II DSNN has been already observed by the convergence simulation result. Table 2 Learning progress for a Type I DSNN network with 10 synaptic connections - in average there are 19 spikes per input sample in 400 milliseconds. Training Sample 1 100 500 1000 1500
Similarity (Mean)
Similarity (Variance)
0.22 0.80 0.93 0.96 0.97
0.007 0.008 0.003 0.002 0.002
Membrane (Mean) 27.06 dB 14.35 dB 2.12 dB -2.47 dB -7.70 dB
Error
Membrane Error (Variance) 0.69 2.98 4.55 2.46 2.30
Table 3 Learning progress for a Type II DSNN network with 11 synaptic connections - in average there are 14 spikes per input sample in 400 milliseconds. Training Sample 1 100 500 1000 1500
Similarity (Mean)
Similarity (Variance)
0.07 0.85 0.95 0.95 0.98
0.005 0.007 0.003 0.003 0.001
Membrane Error (Mean) 29.86 dB 13.32 dB 3.29 dB 5.80 dB -2.36 dB
Membrane Error (Variance) 1.10 4.96 11.12 6.56 12.28
2. Repeatability Analysis For repeatability analysis, the training process is repeated 5 times each time with a new set of parameters. Table 4 shows the minimum and maximum of similarity measure plus the membrane error improvement for the trained network. The results confirm a consistent convergence to an optimum parameter set. Greater accuracy is achievable by increasing number of training samples. The average number of spikes for a trial is about 14 spikes in 400 milliseconds. The results demonstrate that in average 13 out of 14 possible spikes coincide within a 2 milliseconds accuracy to the spike train of the desired DSNN. Similarly, the membrane error measure shows a significant improvement in neuron membrane potential prediction. For the next simulation, scalability of the learning algorithm will be studied. Table 4 Learning progress for a Type I DSNN network with 10 synaptic connections, 1500 training samples for 5 different training/test set Similarity (Min) 0.907
Similarity (Max) 0.971
Membrane Error Drop (Min) 19.01 dB
15
Membrane Error (Max) 34.7 dB
(b) 3
2
2
PSP
PSP
(a) 3
1 0
0
50
100
150
200
250
300
350
1 0
400
2
0
50
100
150
200
250
300
350
400
2
AP
AP
Target Test 1
0
0
50
100
150 200 250 Time(msec)
300
350
1
0
400
0
50
(c)
100
150 200 250 Time(msec)
300
350
400
(d)
Figure 7 Parameter convergence and membrane error in system identification task. task a) Value of ∆, in the target and training networks before training, b) Value of ∆, in the target and training networks after training with 1500 samples, c) Membrane embrane and AP activity of the target and training networks network before training, d) Membrane and AP activity ty of the target and training networks network after training with 1500 samples.
3. Scalability Analysis The learning progress is analyzed for larger numbers of synaptic connections. The previous network structure was replicated but with 5, 10, 20, 40, 80 and 160 1 synaptic connections. Table 5 shows the similarity measure and membrane error for each of these networks network after training with 1500 samples. Table 5 Similarity measure and membrane error improvement for different number of synaptic connections co in a Type I DS DSNN Synaptic connections
Similarity (Mean)
Spike Count (Mean)
Membrane Error (Mean)
5
0.27→0.93
7.97
28.90→1.54 dB
10
0.37→0.91
13.13
25.62→6.60 dB
20
0.19→0.93
9.63
22.12→-3.70 dB
40
0.48→0.93
13.7
17.85→-5.3 dB
80
0.56→0.91
21.7
17.23→-0.81 dB
160
0.72→0.84
54.3
18.28→7.94 dB
Improvements to similarity measure and membrane error are observed independent of the number of synaptic connections. In spite of the fact that a larger network requires more training samples, the improvement in membrane prediction ction and similarity measure is significant for different network size sizes using the same number of training samples. Figure (8) shows the initial guess of synaptic parameters of
16
the Type I network with 80 synaptic connections connection and the trained parameter set after ter training with 1500 samples. The simulation results identify that the convergence to an optimum parameter set is repeatable and scalable. The simulation results also suggest that the neuron membrane potential is a unique function of incoming spike train and synapse parameters. Thus, the learning progress in predicting synapse parameters corresponds to the simultaneous improvement of similarity measure and membrane potential prediction. These results validate the learning capability for accurate predictions of synaptic parameters for both Type I and Type II DSNNs.
(a)
(b) Figure 8 Training progress in a Type I DSNN with 80 synaptic connections. connections a) Desired synapse parameters plus initial guess of parameters, b) Desired and trained synapse parameters trained by 1500 training samples.
The next question is the analysis of similarity range variation in the learning process. Similarity range is the core parameter in quantifying similarity between spike trains, and the following simulation simulations examine its impact on synaptic parameter estimation. 4. Similarity Range Analysis The similarity range defines the acceptable jitter in comparing the generated spiking train to the desired spike train.. Algorithm I utilizes the similarity range for matching the desired and tr training spikes for membrane prediction and eventually synaptic parameter adaptation. Larger values of the similarity range correspond to larger er errors in membrane prediction. Thus, it it is expected that the learning process fails to converge to the optimum parameter set as the similarity range increases increases.
17
30
0.95
25
0.9
0.85
0.8
1
P a ra m e te r E rr.
1
M e m b ra n e E rr (d B )
S im . M e a s u re
Spiking activity in a biological system is precise with an accuracy of about one millisecond. Thus, defining a much larger similarity range is not a viable tool for spike trains similarity assessment. In the following simulations, the similarity range for the learning algorithm changes from 2 to 20 milliseconds with one millisecond steps for a Type I DSNN with 10 synapses. The network is trained with only 500 input stimuli - Table 2 result indicates that 500 training samples is enough for an accurate estimation of the network parameters - , and its learning progress is analyzed on 100 test samples. Figure (9) shows the average similarity measure for the trained network with different similarity ranges. In addition to similarity measure, average membrane error and parameter estimation error are shown. The reference similarity range for performance analysis in the trained network is 2 milliseconds. Similarity measure shows a significant drop as the learning similarity range increases. This drop is predictable as the learning algorithm faces a larger error in membrane prediction. Membrane error grows along with similarity range growth; in fact, the membrane error grows to a level of an untrained network - refer to first row of Table 2 and 3. Similarly, the sum of square error in parameter estimation grows significantly. The simulation results imply that an accurate synaptic parameter adjustment requires a fine similarity range. It was already discussed that the maximum value for similarity range should be smaller than half of the refractory period. The similarity range only within a couple milliseconds - a range from 0 to 7 milliseconds - is a viable choice for the learning algorithm. For the similarity measure, it is possible to maintain the same level of similarity when the similarity range for both the training and the test process is kept equal. Despite maintaining a high similarity measure, the error in parameter estimation and membrane prediction diverges significantly - Figure (9.c). The result suggests that two networks with different synaptic parameters - here desired and training networks - may show a reasonable similarity measure while they have completely different temporal dynamics. Thus, similarity measure and the learning algorithm converge to optimum synaptic parameters for only a fine similarity range- a value smaller than half of refractory period.
20
15
0.8
0.6
0.4
10 0.75
0.2
5 0.7 2
4
6
8 10 12 14 Training Sim. Range
(a)
16
18
20
2
4
6
8 10 12 14 Training Sim. Range
(b)
16
18
20
0 2
4
6
8 10 12 14 Training Sim. Range
16
18
20
(c)
Figure 9 Learning progress as a function of the similarity range. a) Similarity measure versus similarity range. Similarity measure of the trained network is measured with a fixed value of 2 milliseconds, b) Membrane error versus similarity range, c) Error in parameter estimation for different similarity range. A lower similarity range shows a lower membrane and parameter estimation error, and it correlates with higher similarity measure. A longer similarity range is neither a viable choice to evaluate similarity between two spike trains nor an optimum parameter estimation rule for the learning process. Despite the fact that large similarity range maintains a reasonable level of similarity measure, it fails in accurate membrane potential and synaptic parameter prediction.
A complete set of simulations has been applied in examining off-line learning in Type I/II DSNNs. Different aspects of the learning algorithm including convergence, repeatability and scalability were examined. The simulation results verify an acceptable consistency in the system identification task independent of the network size and type. The simulation result along with the theoretical finding Appendix B - verify that the learning algorithm and its extension, gradient-gain learning, are a viable tool for simultaneous STP and LTP adaptation in DSNNs. Through the following simulations other
18
aspects of the DSNN model and the on-line gradient-gain learning will be analyzed. The first simulation examines the FD model’s capability in the prediction of CA1 pyramidal neuron spiking activity. B. FD Model and CA1 Neuron Spiking Activity Prediction A single pyramidal cell of the CA1 hippocampus area was excited by injecting APs on its presynaptic pathway with a 2-Hz Poisson spike timing. The patch-clamp technique was utilized for recording the somatic membrane potential of the CA1 neuron in response to the spiking train impinging on its Schaffer collateral synaptic pathway – the data was collected in the USC Center of Neural Engineering, Theodore W. Berger Lab -. Figure (10.a) shows 10 seconds of the CA1 pyramidal neuron membrane potential and its AP activity in response to the impinging spike train. Using an exhaustive search with a weighted least square criterion, free parameters of a single synapse model are adjusted to replicate the PSP dynamics and its spiking activity. In order to have a better estimate of the neuron's PSP response and spiking activity, the synapse model PSP response was replaced with the neuron's PSP extracted from the recording. Figure (10.b) presents the PSP prediction using the FD model plus its projected spiking activity. The average normalized error in PSP amplitude prediction is 27%. The model predicts CA1 neuron spiking activity with 6.6 milliseconds accuracy while generating only one extra spike. The parameter adaptation was applied for the whole recording time - 215 seconds - where the recorded neuron’s PSP dynamics showed other nonlinear dynamic processes including facilitation regulation and vesicle depression. For a short recording period, adjusting only ∆, parameters suffice to build a precise predictive model of the recorded neuron’s PSP response whereas a longer period was required to capture other dynamics such as facilitation regulation and vesicle depression. Appendix A describes a more intricate FD model for synapse dynamics. The on-line learning algorithm was applied in parameter prediction with a 7 milliseconds similarity range. Facilitation and recovery time constants were fixed by the outcome of the exhaustive search in which both ∆ and parameters converged to the values predicted by the exhaustive search algorithm. The membrane potential dynamics of the CA1 neuron was modeled by a single synapse with fixed time constants. However, multiple synaptic connections with different time constants - and probably synaptic delays - generate a more precise prediction of the neuron PSP dynamics. Table 6 shows the release prediction error using single and double synaptic connections. Increasing synaptic connections correlates with a better prediction of CA1 membrane potential dynamics and, correspondingly, its spiking activity prediction. It was already shown in the system identification section that the learning algorithm is a scalable process for multiple synaptic connections. Thus, the proposed learning algorithm will address synaptic parameter adjustment for multiple synaptic connections as well. Table 6 Performance of the FD model with multiple synaptic connections in predicting a CA1 neuron PSP Synapse parameters ∆T, N¼½¾ , SU , ST 0.115,10.132,500,200
0.115,10.132,500,200
0.115,10.132,500,200
Membrane error
Release error (%)
2.66
27.3%
Spike similarity (milliseconds) 6.6 + 1 extra spike
2.76
22.4%
6.3 + 1 extra spike
Through modeling the PSP response of the FD model was replaced by the realistic PSP shape of the CA1 neuron to provide a more accurate AP timing; this change was independent of the learning process or synapse temporal dynamics. Thus, it is replaced with its phenomenological approximate in the next simulations. Though it is shown that the multiple FD model with different time constants is a better approximate of the pyramidal neuron membrane potential dynamics, through the following simulations each synapse is modeled with only one FD model. However, it is possible to increase the number of synaptic connections to create a more intricate model where a set of synapses receives the same input spiking activity. The next section introduces a generative model to generate different spiking patterns of cortical neurons.
19
1
AP
0.8 0.6 0.4 0.2 0 0
1
2
3
4
5 Time(sec)
6
7
8
9
10
(a) 120
CA1 PSP/AP
100 80 60 40 20 0 -20
0
1
2
3
4
5
6
7
8
9
10
(b) 50 FD CA1 Threshold
CA1/FD PSP
40 30 20 10 0 -10
0
1
2
3
4
5
6
7
8
9
10
5
6
7
8
9
10
(c)
CA1/FD AP
2 CA1 FD
1.5 1 0.5 0 0
1
2
3
4
(d) Figure 10 CA1 PSP dynamics and synapse FD model response. a) Input spike train to a CA1 pyramidal neuron synaptic pathway, b) A single CA1 neuron PSP response and its AP activity, c) Predictive model of CA1 neuron PSP dynamics plus a synapse FD model PSP response, d) Spiking activity of the CA1 neuron plus the synapse FD model. The FD model with only one synaptic connection predicts a CA1 neuron spiking activity with an average of 6.6 milliseconds jitter accuracy plus firing an extra spike.
C. Cortical Neuron Spike Generative Model The Hodgkin-Huxley model is the most well-known biological model capable of generating the different spiking patterns of cortical neurons. Izhikevich proposed a mathematical model capable of generating different AP patterns. A Type II DSNN model can imitate the spiking patterns claimed in Izhikevich model. In the following simulations, synaptic parameters in a Type II DSNN network with only one feed-forward synaptic connection are adjusted to generate different spiking patterns of cortical neurons. The proposed DSNN model builds a spike-domain computation unit capable of mapping incoming spike trains to different spiking and bursting patterns. The network parameters - ∆ , , ∆ ¿ , ¿ - are adjusted to generate different spiking patterns in response to the input spike train. The learning process is on-line, and the model parameters are updated at the gain injection time. Input stimuli to the network have an ISI of 50 milliseconds and lasts for 300 milliseconds. The learning process is repeated 500 times, and in compliance with the on-line learning algorithm, the input stimuli are repeated with a 300 milliseconds silent gap - with this time gap, the network synapses retrieve from their previous release activities. The following gain function adjusts the synaptic parameters to generate the desired spiking pattern:
20
' b ∗ , b 0,1,2, … . 1 missing spike ∈ E' , ' , ' À1 extra spike ∈ E' , ' ,
' b ∗ , b 0,1,2, … .Å 0 Other time
(12 12)
definition is set to 3 milliseconds; thus, the network is trained to The similarity range - - in the '
generate the desired spike train within with 3 milliseconds accuracy. The training raining process for different spiking patterns start from the same initial value for the synaptic parameters. The he facilitation factor ∆ and ∆ ¿ - are set to 0.2 and 0.5, and synaptic strength - and ¿ - are set to 5.0 and 2.0. Training stops whenever the sum of absolute gain g is less than or equal to 1; the stopping ing criterion implies that the trained network is capable of generating at least all but one of the desired spikes within 3 millisecond accuracy. Figure (11) shows the learning progress and gain pattern for chattering - CH spiking activity. Chattering neurons urons can fire high-frequency bursts of 3 to 5 spikes with a relatively short inter-burst burst period in response to incoming stimuli. Figure (12) shows the neural model and the trained parameter set for different spiking activities of cortical neuron. In contrast to the system identification problem, the gain function is determined independent of membrane potential error. The gain function reflects the coincidence of the desired and training APs in a specific time interval; therefore, synaptic parameters are adjusted justed by an evaluative function of matched AP time rather than the network membrane error. Neural models with static synaptic connections fail to emulate ate the diverse spiking activities of a cortical neuron whereas dynamic synapses allow the neural model too generate different spiking and bursting patterns. The model introduced here builds ann adjustable spike-domain computational unit. The feed-forward forward and feed-back synapses of the spike generative model are able to introduce an adjustable delay in their response to the impinging APs.. By utilizing the synaptic delay, the model is capable of generating a diverse spiking pattern. The model odel response to incoming spike activity shows a diverse set of neural processing capabilities; the he processing consists of variable spike time delay, temporal selectivity of incoming activity, and neural synchrony. The next simulation utilizes a Type I network for a temporal pattern recognition task. 1 0.8 AP
0.6 0.4 0.2 0 0
50
100
150 Time (msec)
200
250
300
200
250
300
200
250
300
(c) 3
P S P /AP
2
(a)
1 0 -1 -2 -3
50
100
150 Time (msec)
(d) 1
G a in
0.5 0 -0.5 -1 0
50
100
150 Time (msec)
(e)
(b)
Figure 11 Parameter adaptation in the neuron model for the chattering - CH - spiking activity. a) Feed-forward forward synapse parameter adaptation, b) Feed-back back synapse parameter adaptation, c) Input spike train to cortical neuron, d) PSP and AP in neuron output, e) Gain pattern in response to incoming APs.
21
R S / G a in
Ref
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
T C 1 / G a in
50
100
150
200
250
300
T C 2 / G a in
50
100
150
200
250
300
R Z 1 / G a in
∆ ¿ , ¿
50
100
150
200
250
300
R Z 2 / G a in
2 1 0 -1
50
100
150
200
250
300
R Z 3 / G a in
∆ ,
50
100
150 Time (msec)
200
250
300
I B / G a in
2 1 0 -1
C H / G a in
(a)
F S / G a in
2 1 0 -1
L S / G a in
2 1 0 -1 2 1 0 -1
(b)
2 1 0 -1 2 1 0 -1 2 1 0 -1 2 1 0 -1
(c) 1
AP
0.8 0.6
2 1 0 -1
0.4 0.2 0
0
50
100
150 Time (msec)
200
250
300
(e)
(d)
Figure 12 Cortical neuron model plus different spiking patterns. patterns a) Cortical neuron model. The model is a Type II DSNN with one feed-forward synaptic connection, b) Parameters for the feed-forward forward synapse. There are 10 different parameter set, each one corresponds to one of the desired spiking activity. The RS spiking activity refers to a regular spiking, in which cortical neuron model generates an AP in response to each impinging AP. The The IB spiking activity refers to intrinsically bursting spiking activity; the cortical neuron fires multiple APs for the first incoming AP and it responds with regular spikes afterward. The CH spiking activity corresponds to a chattering neuron; the cortical cortical neuron responds to each impinging AP with multiple APs. Neuron generates more APs in responds to the first AP, and number of APs drops in the next spikes. The FS spiking activity corresponds to fast response neurons, in which neuron fires multiple spikes spike per impinging APs. The LTS refers to a low threshold spiking activity, which is similar to FS expect the neuron fires lower number of APs per incoming AP. The TC1 and TC2 spiking activity refers to Thalamo-Cortical Cortical neuron response which has a temporal selectivity. lectivity. The TC1 neuron responds to incoming AP with a delay, and the TC2 is responding to the first couple of APs. The RZs neuron define other possible temporal selectivity of cortical neurons,, each neuron is responsive to specific speci APs of incoming spike train, c) Parameters arameters for feed feed-back synapse, d) Input stimuli to cortical neuron model, e) Cortical ortical neuron model response for different spiking activity plus gain pattern. Cortical neuron training stops whenever there is an absolute gain sum less than or equal to 1.
22
B. Temporal Pattern Classifier A DSNN based neural classifier is built for recognizing different classes of spiking activity - 3 -. A Type I DSNN consisting of neurons with H synaptic connections per neuron - H 60 in this problem - is built in which synapses project incoming spiking activity to their corresponding neurons. Each neuron is responsible for encoding a specific class and it is trained to fire a higher number of APs for its corresponding class. In fact, the synaptic connections of the network bind APs in multiple spike trains to build the desired post-synaptic firing activity. The synaptic parameters of each neuron are adjusted to be sensitive to temporal features of the corresponding activity class while being insensitive to the other activity classes. The gain function of the learning algorithm is responsible for observing neurons' spiking activity to adjust network synaptic parameters by injecting the proper gain. Equation (13) defines the gain function for each output neuron at different processing time: 1 ∑ÉÊ APÆ n t k ∈ Winning Class Å Æ ' Ç ∑ÉÊ APÆ n t
k ∈ Other Classes
(13)
Parameter defines the ISI for the spiking activity in the winning class. The learning process utilizes competitive learning among the neurons. Each neuron tries to fire more APs in response to its matched incoming activity while at the same time forcing the other neurons to be silent. It is also possible to define new gain functions built upon other temporal characteristics of incoming spiking activity. The gain function can train individual neurons to fire APs sooner in response to their matched incoming spiking activity while postponing firing activity in other neurons to a later time. This second gain function is defined by: Æ ' Í
min 1,2 ∑ÉÊ APÆ n t k ≡ Winning Class, n ¬ TÐ ∑ÉÊ APÆ n t
k ∈ Other Classes, n ¬ TÐ Å 0 n TÐ
TÐ Time of Ñirst AP in winning class D,
D deÑines time gap margin
(14)
Using two different gain functions, the synaptic parameters of the network are adjusted to capture different features of the incoming spiking trains. The first gain function tends to adjust synaptic parameters of the network for a prolonged spiking activity in the winning neurons while the second gain function adjusts synaptic parameters for an abrupt AP firing. To maintain a prolonged spiking activity, the synaptic parameters of a neuron generally move toward potentiating synapses with strong synaptic connections. In contrast, synaptic parameters of a fast firing neuron move toward a depressing synapse making the synapses sensitive to isolated spiking activity. The optimum synaptic parameters are dependent on the temporal characteristics of the incoming spiking trains as well. The discriminative information of a class are encoded in the spike time of individual spike trains and the mutual spiking activity of different spiking trains. The rule of the gain function is to correlate incoming spike train features with the synaptic parameters of the neurons to generate a desired spike pattern in the network output. Thus, optimum synaptic parameters of a network become a function of the incoming spike train, the synaptic parameters, and the gain function. Network training stops whenever the gradient gain converges to zero. The gain function converges to zero when the neurons generate the desired spiking activity. Figure (13) shows the network model plus a sample spiking train for each class. The time period of the input stimuli is 400 milliseconds. There are 500 different samples of spiking activities per each activity class. The network is trained with the first 400 samples, and its performance is examined using the remaining 100 samples. The free parameters include 180 pairs of ∆, while the parameters of the inhibitory synaptic connections are fixed at 1.0,1 . The function of the inhibitory synapse is to reduce the activity of all the neurons except the winning neuron whenever there is an uncertainty about the incoming spiking activity. Whenever there is a distinct class of incoming activity, the winning neuron fires with a higher firing rate in comparison to the other neurons without the help of inhibitory
23
synaptic connection. Facilitation factor in the inhibitory synaptic connections are manually set to 1.0 with synaptic strength in order of 1.0. Therefore, inhibitory synapse has depressing dynamics; synapses with depressing dynamics are sensitive to sparse neural activity which matches the dynamics of the network spiking activity in response to indistinct input stimuli. Ref
30
25
Synapse
20 Ref
15
10 Ref
5
0 0
150 200 250 Time (msec)
300
350
400
300
350
400
30
25
25
20 Synapse
20 Synapse
100
(c)
(a) 30
15
15
10
10
5
5
0 0
50
0 50
100
150 200 250 Time (msec)
300
350
400
0
50
100
150 200 250 Time (msec)
(d)
(b)
Figure 13 Temporal pattern classifier and sample spiking activities. a) Classifier architecture. Synaptic connections to other neurons have the same structure of the first neuron, b) Sample spiking activity of class A - there are 30 input spike trains, each spike train feeds two synaptic connections, c) Sample spiking activity of class B, d) Sample spiking activity of class C.
For training purposes, input stimuli are concatenated with a 400 millisecond silent time gap between different training samples. The initial values for the synaptic parameters are fixed for all three classes, and learning rate is set to 0.001. The learning process repeats 10 times on the training data before the classification performance is examined on the test samples. Table 7 shows the network performance using both higher spike count and time of the first AP gain functions. Network performance using spike count on trained data reaches 95.6% at the end of the training process, and its performance on the test sample is 96.5%. Network performance using spike time reaches 76.6% with an average delay time of 38 milliseconds - less than 16 samples - on correctly classified samples, and the network performance in the test samples is equal to 80.8% with an average delay time of 36 milliseconds - less than 15 samples. The delay time determines the time period until the first AP fires in the wining class. The delay time defines the time window that the network requires to reach a decision about the class of activity. The performance with the spike time classifier implies that the trained network only requires 38 milliseconds - less than 16 samples as input data has been sampled with 400 samples per second of data to decide the activity class with an 80% accuracy. The spiking activity provided in this section belongs to a chain link fence intrusion activity database the fence activity database was collected by USC Neural Dynamics Laboratory in 6 different fence panels which include 300 recordings of 25 seconds - with three different classes of sequential data. The database includes 3-axis vibration signals recorded with 12 bit samples with a sampling rate of 400 Hz. In generating spike trains, the input signal is filtered by a filter bank - consisting of 5 different frequency bands -, and it is then passed to the absolute value function. Each pre-processed output is magnified and then passed to an integrate-fire neuron with a fixed firing threshold. The neuron fires an AP whenever its membrane potential passes the firing threshold. There are two gain levels for each
24
pre-processed output; thus, there are in total 10 neurons mapping input stimuli to spike trains for each axis. In all, there are 30 different spike patterns for each 400 milliseconds of a fence activity. Each spike train feeds two synaptic connections where the initial values of the synaptic parameters in these two synapses are different. To compare the network performance with well-established artificial neural networks, a 2-layer neural network with 15 input, 50 hidden, and 3 output neurons are trained to classify the same activity classes. The feature set of the artificial neural network are defined by energy of the 5 filter banks for each vibration axis. The time period of each activity is partitioned in two parts, and the energy of each partition is added to the feature set. Table 7 shows the performance of the neural network and a Gaussian mixture model - GMM - in classifying the same activity classes. The DSNN performance is comparable with GMM and 2-layer neural network classifiers; however, the number of free parameters for the DSNN are much lower than the other classifiers. The DSNN performance is almost independent of feature selection whereas the performance of other classifiers is dependent on the choice of features. Table 7 Performance of different classifiers on the test dataset of the fence intrusion database Classifier
Performance (%)
Parameters
GMM
97.4
GMM
90.7
Neural Network
97.8
Neural Network
93.9
DSNN (AP Count)
96.5
Gaussian mixture model with 2 components Energy in the first and second half of input are considered as the feature set 30 features per input sample 1862 free parameters (900+900+30+30+1+1) Gaussian mixture model with 2 components 15 features per data 482 free parameters (225+225+15+15+1+1) 2 hidden layers, each layer has 40 neurons Energy in the first and second half of input are considered as the feature set 30 features per data More than 4000 free parameters 2 hidden layers, each layer has 40 neurons 15 features per data More than 2000 free parameters 360 free parameters (60×2×3)
DSNN (First AP)
80.8
360 free parameters (60×2×3)
The DSNN network is also capable of deciding about the incoming activity class in a shorter time period than the other classifiers. Feature extraction in the other classifiers requires a pre-defined time window; otherwise, variability in the feature set significantly drops the classifiers' performance. Neural network and GMM classifiers are incapable of performing on original data; in fact, their learning algorithms are incapable of maintaining a reasonable performance on large feature sets with extensive numbers of free parameters. The DSNN classifier performance along with its superior features: 1) lower number of free parameters, 2) optimum decision in minimal time, and 3) a feature-free classifier make the model the superior choice for temporal pattern classification. V.
DISCUSSION
This research introduces an artificial model incorporating two progressive transitions to develop a more realistic neural model: 1) synaptic response is not a static process; it is a dynamic process changing in the time-course of impinging spiking activity, and 2) synaptic adaptation is not merely a synaptic weight regulation; it is a dynamic process regulating both long-term and short-term characteristics of each synapses. The novel neural model introduces a family of neuro-computational processing engines applicable in modeling a cortical neural circuitry's spike-domain functionality and a variety of spatio-temporal processing tasks. A linear state space model of synaptic transmission was introduced providing an accurate and compact mathematical representation of the synaptic dynamics. A comprehensive analytical study of the
25
synaptic response including its efficacy, temporal dynamics, and statistical behavior is provided as well. The dynamics of synaptic response with respect to variation in its internal parameters is formulated; the synapse PSP gradient with respect to its parameters is the key tool utilized in DSNN training. The linear state space representation is extended to the neural level where a DSNN includes multiple neurons with thousands of synaptic connections. A supervised spike-in-spike-out learning rule for the synaptic adaptation of DSNNs is proposed. The learning rule is biologically plausible, and it is capable of simultaneously adjusting both LTP and STP factors of individual DSNN synapses. The learning algorithm is extended to an on-line learning algorithm called gradient-gain learning. Gradient-gain learning is a spike-domain learning rule that adjusts synaptic parameters of individual DSNN synapses to perform a spike-domain functional mapping. A comprehensive simulation analysis is performed to verify accuracy, repeatability, and scalability of the learning algorithm. The system identification test confirms the off-line learning algorithm ability to predict synaptic parameters of a DSNN containing up to 160 synaptic connections. The system identification task also confirms consistency and convergence in the synaptic adaptation process for both Type I/II DSNNs. A new similarity measure is introduced to assess similarity between two spike trains. The similarity measure identifies each spike of the test train - in comparison with the desired spike train - to three classes of: 1) similar, 2) missing, or 3) extra where the labeling process is insensitive to spiking jitter. Similarity range is the free factor in the spike similarity measure defining the maximum acceptable jitter. The similarity range plays a crucial role in the synaptic adaptation mechanism; this is because the learning objective function minimized by the supervised learning rule is an approximate measure of the similarity assessment process. The relationship between spike train similarity and neural membrane dynamics in the synaptic adaptation process is analyzed to estimate the optimum similarity range, where the simulation result estimates an optimum similarity range of about a couple of milliseconds. The synapse model’s ability in predicting PSP dynamics - and spiking activity - of CA1 pyramidal cell is demonstrated. A single synaptic connection with only 2 free parameters provides the same accuracy of second-order Volterra series approximation with 14 free parameters (Ghaderi, 2012). Gradient-gain learning with a /1 gain for missing and extra spikes has been applied for synaptic parameter adjustment where the parameters converge to the same parameters optimized using an exhaustive search. The DSNN is intrinsically a time-domain computational engine; it is applicable for processing dynamic neural/artificial signals. Gradient-gain proposes a flexible synaptic adaptation mechanism to shape DSNN temporal processing. A Type II DSNN model with only a single feed-forward synaptic connection is introduced to generate different spiking/bursting patterns of the cortical neurons. Gradient-gain with a /1 gain at missing and extra spikes shapes the model dynamics for the desired mapping. This recurrent neural structure introduces an adjustable generative model capable of performing a wide set of spike-domain mappings. A DSNN model with multiple output neurons has been applied in performing a temporal pattern classification task. The DSNN classifier performance competes with a 2-layer neural network and surpasses the GMM classifier - 96.5% classification rate -, while using a significantly fewer number of free parameters. The gain function tunes synaptic parameters of the winning neuron by injecting 1 gain resulting in a regular spike train - with time period - in response to the matching input while the other neurons are punished by 1 gain for their spiking activity in response to un-matched incoming input. It is also shown that DSNN can be trained to form an optimum decision in a minimal time, a property which is unachievable by the other competitor models. The DSNN is capable of performing the classification task with an 80% accuracy after observation only 36 milliseconds of the data - each data lasts for 400 milliseconds. Again, the gain function rewards the winning neuron with 1 gain and punishes the other neurons with 1 gain to bring winning neuron’s spiking activity to an earlier time while pushing other neurons’ spiking activity to a later time.
26
Synapse temporal dynamics is the natural substrate for neural computation. Neither of the spikedomain mappings introduced in this research are feasible for the networks with static synapses - or exclusively synaptic weight adaptation. In fact, simultaneous regulation of STP and LTP factors are required to shape the neural computations performed in this research. Temporal processing in a static synapse model is limited to the overlap time between PSPs of consecutive APs - practically less than 100 milliseconds -, whereas DSNN temporal processing is expandable to longer time periods using LTP and STP factor regulation. DSNN performs a diverse set of spike-domain processing and can be performed on a sparse spiking activity or a sequence of spike train expanded to multiple seconds. The simulation results for CA1 spike activity prediction and the temporal pattern classifier show how DSNN model performance surpasses state-of-the-art mathematical tools for neural modeling and dynamic classification tasks. The DSNN gains its computational power by: 1) augmenting synaptic dynamics, and 2) incorporating simultaneous LTP and STP adaptation. The DSNN model and the spike-generative model can be used for both modeling neural circuitries and their spiking activity prediction. VI.
CONCLUSION
This research is an attempt to show how a transition from a static synapse model to a dynamic synapse model that also includes STP and LTP adaptation is achievable. It is also shown how the computational power of the artificial neural model increases by augmenting accurate models of the corresponding neural circuitry. There are many other unknown questions about neural dynamics and the synaptic adaptation mechanism. The learning addressed here is limited to a single layer DSNN but can be extended to multiple layers or a cascading neural topology structure as observed in the brain. The spike-domain mappings performed here address a small sub-set of possible spike domain mapping observed in the brain. The next question is how the learning and neural model can be modified to perform a more diverse neural processing. The authors’ future research focuses on applying DSNN modeling in predicting spiking activity of the hippocampus cortical region where it requires a cascading neural structure.
APPENDIX A DETAILED FD MODEL This section provides the intricate model of synapse temporal dynamics. The model addresses: 1) multiple facilitation mechanisms, 2) the Ca2+-dependent vesicle recovery time constant, 3) activitydependent facilitation regulation, and 4) activity-dependent vesicle depletion. The FD model proposed in equation (1) is a simplified derivation of the following intricate model, in which by an appropriate setting of the model parameters it turns to the simplified FD model already introduced in equation (1). The calcium facilitation factor controlling calcium concentration in synapse on the pre-synaptic side changes through the time course of incoming spiking activity. In addition, multiple facilitation mechanisms with different time constants have been addressed for the calcium concentration kinetics. Equation (A.1) to (A.4) provide two facilitation mechanisms in which their Ca2+ facilitation factor is a function of impinging spiking activity to the synapse. ⁄ \ 1 a ∗ ∆a ∗ a ⁄ a 1⁄a
⁄
∗ a ∗ ⁄ 1 ∗ ∆ ∗ a
¿ ⁄ ¿ \ 1 a¿ ∗ ∆a¿ ∗ a¿ ⁄ a¿ 1⁄a
(A.1) (A.2) (A.3)
¿ ∗ a¿ ∗ ¿ ⁄ ¿ ¿ Ô ¿ 1 ¿ ∗ ∆ ¿ ∗ a
(A.4)
? ⁄ ? 1⁄? ⁄ Õ 1 ? ∗ ∆? ∗
(A.5)
Calcium-dependent vesicle recovery regulates the recovery time constant of the released vesicles in response to prolonged spiking activities. The vesicle recovery time constant decreases in response to incoming spiking activities guaranteeing a synaptic response during long periods of impinging spikes. Equation (A.5) defines the Ca2+-dependent mechanism in vesicle recovery time constant regulation.
27
The vesicle release dynamics is defined by:
⁄ 1 ∗ ? ∗ ? ∗ 1⁄ ¿ ∗ ¿ ∗ ∗
(A.6)
The other crucial mechanism in post-synaptic kinetics is the vesicle depletion process. The ability of post-synaptic side to absorb the released vesicles drops after a prolonged spiking activity. The following equation models the permeability of the post-synaptic membrane in response to pre-synaptic spiking activities. ⁄ 1 ⁄ É ±! ∗ ∆ ∗
(A.7)
The synaptic PSP response is defined by:
⁄ ⁄ ∗ ∗ ¿ ∗ ¿ ∗ ∗
(A.8)
(A.10)
⁄ ⁄ ∗ ∗ ¿ ∗ ¿ ∗ ∗
(A.9)
Equation set (A) defines the intricate model of synaptic response. The PSP response is modeled by an approximate α-synapse function, but it is possible to be defined with other PSP shapes.
APPENDIX B LINEAR OBJECTIVE FUNCTION, SIMILARITY MEASURE AND UNIQUNESS OF DSNN MEMBRANE POTENTIAL In this section, several theories are developed to support the convergence of the proposed learning algorithm. The first theory addresses the relationship between the learning objective function and the similarity measure. The second theory discusses the relationship between the neuron membrane potential and synaptic parameters in Type I/II DSNNs. A few of corollaries extends the relationship for spiking activity in Type I/II DSNN and its synaptic parameters. Theorem I: A decrement in the expected value of the learning objective function - defined in equation (10.a) - corresponds to an increase in the number of similar spikes. Proof: In a DSNN with excitatory synaptic connections, the PSP derivative in both the test and target DSNNs have the same sign - independent of ∆, parameters value - at different processing bins. With the mentioned assumption, the expected value of error at missing and extra spike times is defined by: ccÖ × c∑!∈ØÙn·p ! k ∑!∈ÚÛmmÛÜ k ! nÚÛmmÛÜ nØÙn·p ∗ c
It is assumed that ! gets a uniform distribution between E0 k . Parameter c is the expected value of membrane potential error ∎ Theorem II: The membrane potential of a synapse is a unique function of its synaptic parameters. The free parameters of each synapse are ∆, . Proof: The PSP response in the DSNN is described by equation (4.a). The BC and DC matrices are defined independent of the free synapse parameters, and they are the same for different synapses. It can be shown that B and D elements are unique functions of synapse free parameters. Thus, each parameter set generates a unique membrane potential in response to different pattern of spiking activities. (B.1)
The theorem is valid by the ApÞ AßÞ and BÞp BÞß assumption, and it is not necessarily valid for DSNN models with other free parameters including recovery and facilitation time constants ∎ Corollary I: The membrane potential of a DSNN network is a unique function of its synaptic parameters. Again, the free parameter of each synapse are ∆, . Proof: The DSNN network dynamics follow the same state space model defined for a single synapse. Thus, it is possible to consider DSNN with multiple synaptic connections as a single synapse; therefore, Theorem II is valid for the DSNN network ∎
28
The learning algorithm deals with spiking activity, but the more important question is the relationship between DSNN spiking activity and its synapse parameters. In contrast to membrane potential uniqueness, spiking activity is not a unique function of synaptic parameters. For example, DSNN networks failing to fire APs can share similar spiking activities while their synaptic parameters can be different. Thus, the spiking activity and synaptic parameters relationship can be discussed in only DSNN networks firing APs. It is also assumed that the processing resolution is fine enough that neurons fire APs at k à in which à converges to zero. DSNN networks with different numbers of synaptic connections fail to generate similar spiking activities; it can be assumed that the extra synaptic connections perturb the neuron membrane potential about à - around a spike time - changing its spiking time. Corollary II: DSNN networks with similar spiking activities in response to all possible input spike trains have similar synaptic parameters. Proof: Assume that all except one - kth synapse - synaptic connection of two DSNNs share the same parameters. Moreover, assume that there is a specific pattern of spike requiring the involvement of the kth synapse to fire an AP. Corollary I showed that the PSP of the remaining synaptic connections is the same. Thus to have a similar spike time, the kth synapse in both networks is required to have the same synaptic parameters. As a result, precise spike timing in two DSNN networks in response to all possible input spike trains requires the same parameter set and network topology ∎ The simulation result in the system identification task shows that the distance between free parameters of two DSNNs shrinks as their similarity measure increases. The convergence of parameters in the test and training network demonstrates the result of corollary II. In fact, the simulation is performed with a limited number of input spike trains and a coarse processing resolution. The coarse processing resolution relaxes the à constraint which introduces deviation in the synaptic parameter estimation. The corollary's result implies that there is a possibility of multiple parameter sets that can generate similar spiking activity whenever the number of input spike patterns is limited. The simulation result is in accordance with this assumption. As the number of input spike trains increases, the expected value of membrane potential decreases which always correlates with a higher value of similarity measure.
APPENDIX C GAIN-STDP The synaptic PSP in response a spike train can be written as: u' ∑!âã ' ' ∗ B' ∗ á' '
(C.1)
ä ä ä u' × ' ' ∗ á' ' ∗ B'
(C.2)
' is PSP response modeled by α-synapse function, B' is the vesicle release at ' spike time, and á' is the step function. ' is a decaying function at its tail. Thus the PSP response can be ä -:. approximated by its vesicle release at the last impinging spike time - ' The PSP derivative with respect to synaptic parameters is defined by: ä ⁄ ä ä u' ⁄ × ' ' ∗ á' ' ∗ B'
u' ⁄∆ × ' ä B'
ä '
∗
ä ⁄ B' ∆
ä '
∗ á'
(C.3.a)
ä '
(C.3.b)
- is a function of synaptic parameters and the history of incoming - vesicle release at ä is a linear function of , but it is a non-linear function of spikes - equation set (4) and (5). B' ä can be written as: facilitation factor and timing of previous spiking activity. B' ä ät# ät1 ä (C.4) B' j ∗ F∆, ' , ' , ' , ⋯ , ' The gradient-gain at gain injection time ' is defined by: ä ä ät# ät1 ä , ' , ' , ⋯ , ' ← ' ∗ ' ' ∗ á' ' ∗ F∆, ' ä ∆ ← ∆ ' ∗ ' ' ∗ j
∗ ∂F
ä ät# ät1 ⁄ , ' , ' , ⋯ , ' ∆ ∆, '
29
ä ∗ á' '
(C.5.a) (C.5.b)
Equation set (C.5) is a function of pre-synaptic activity and gain injection time regulated by the gain value. Equation set (C.5) defines the gain-STDP learning rule. Function defines the normalized vesicle release. The F function and its derivative are described in equation set (4) and (5).
REFERENCES
1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13.
14.
15. 16. 17. 18. 19. 20. 21. 22. 23. 24. 25. 26.
Abbott, L. F., and Wade G. Regehr. "Synaptic computation." Nature 431.7010 (2004): 796-803. Morrison, Abigail, Markus Diesmann, and Wulfram Gerstner. "Phenomenological models of synaptic plasticity based on spike timing."Biological cybernetics 98.6 (2008): 459-478. Natschlger, Thomas, Wolfgang Maass, and Anthony Zador. "Efficient temporal processing with biologically realistic dynamic synapses." Network: Computation in Neural Systems 12.1 (2001): 75-87. Liaw, Jim-Shih, and Theodore W. Berger. "Dynamic synapse: Harnessing the computing power of synaptic dynamics." Neurocomputing 26 (1999): 199-206. Tsodyks, Misha, Klaus Pawelzik, and Henry Markram. "Neural networks with dynamic synapses." Neural computation 10.4 (1998): 821-835. Dittman, Jeremy S., Anatol C. Kreitzer, and Wade G. Regehr. "Interplay between facilitation, depression, and residual calcium at three presynaptic terminals." The Journal of Neuroscience 20.4 (2000): 1374-1385. Mongillo, Gianluigi, Omri Barak, and Misha Tsodyks. "Synaptic theory of working memory." Science 319.5869 (2008): 1543-1546. Liaw, Jim‐Shih, and Theodore W. Berger. "Dynamic synapse: A new concept of neural representation and computation." Hippocampus 6.6 (1996): 591-600. Kreuz, Thomas, et al. "Measuring spike train synchrony." Journal of neuroscience methods 165.1 (2007): 151-161. Gütig, Robert, and Haim Sompolinsky. "The tempotron: a neuron that learns spike timing–based decisions." Nature neuroscience 9.3 (2006): 420-428. Ponulak, Filip, and Andrzej Kasiński. "Supervised learning in spiking neural networks with ReSuMe: sequence learning, classification, and spike shifting."Neural Computation 22.2 (2010): 467-510. Bohte, Sander Marcel. "Spiking neural networks." Unpublished doctoral dissertation, Centre for Mathematics and Computer Science, Amsterdam(2003). Ghosh-Dastidar, Samanwoy, and Hojjat Adeli. "A new supervised learning algorithm for multiple spiking neural networks with application in epilepsy and seizure detection." Neural Networks 22.10 (2009): 14191431. Yousefi, Ali, Alireza A. Dibazar, and Theodore W. Berger. "Supervised learning in a single layer Dynamic Synapses Neural Network." Neural Networks (IJCNN), The 2011 International Joint Conference on. IEEE, 2011. Izhikevich, Eugene M. "Simple model of spiking neurons." Neural Networks, IEEE Transactions on 14.6 (2003): 1569-1572. Gerstner, Wulfram, and Werner M. Kistler. Spiking neuron models: Single neurons, populations, plasticity. Cambridge university press, 2006. Mainen, Zachary F., and Terrence J. Sejnowski. "Reliability of spike timing in neocortical neurons." Science 268.5216 (1995): 1503-1506. Tiesinga, Paul HE, J-M. Fellous, and Terrence J. Sejnowski. "Attractor reliability reveals deterministic structure in neuronal spike trains." Neural computation 14.7 (2002): 1629-1650. Berry, Michael J., David K. Warland, and Markus Meister. "The structure and precision of retinal spike trains." Proceedings of the National Academy of Sciences 94.10 (1997): 5411-5416. Hatsopoulos, N., et al. "At what time scale does the nervous system operate?."Neurocomputing 52 (2003): 25-29. Zucker, Robert S., and Wade G. Regehr. "Short-term synaptic plasticity."Annual review of physiology 64.1 (2002): 355-405. Abbott, L. F., et al. "Synaptic depression and cortical gain control." Science275.5297 (1997): 221-224. Maass, Wolfgang, and Anthony M. Zador. "Dynamic stochastic synapses as computational units." Neural Computation 11.4 (1999): 903-917. Varela, Juan A., et al. "A quantitative description of short-term plasticity at excitatory synapses in layer 2/3 of rat primary visual cortex." The Journal of neuroscience 17.20 (1997): 7926-7940. Bohte, Sander M., Joost N. Kok, and Han La Poutre. "Error-backpropagation in temporally encoded networks of spiking neurons." Neurocomputing 48.1 (2002): 17-37. Booij, Olaf. "Temporal pattern classification using spiking neural networks."Unpublished master's thesis, University of Amsterdam (2004).
30
27. Davis, Bryan A., et al. "Supervised synaptic weight adaptation for a spiking neuron." Neural Networks, 2003. Proceedings of the International Joint Conference on. Vol. 4. IEEE, 2003. 28. Gütig, Robert, and Haim Sompolinsky. "The tempotron: a neuron that learns spike timing–based decisions." Nature neuroscience 9.3 (2006): 420-428. 29. Fang, Huijuan, Yongji Wang, and Jiping He. "Spiking neural networks for cortical neuronal spike train decoding." Neural computation 22.4 (2010): 1060-1085. 30. Storck, Jan, Frank Jäkel, and Gustavo Deco. "Temporal clustering with spiking neurons and dynamic synapses: towards technological applications." Neural networks 14.3 (2001): 275-285. 31. Rostro, Horacio, et al. "Back-engineering of spiking neural networks parameters." arXiv preprint arXiv:0905.4810 (2009). 32. Song, Sen, Kenneth D. Miller, and Larry F. Abbott. "Competitive Hebbian learning through spike-timingdependent synaptic plasticity." Nature neuroscience 3.9 (2000): 919-926. 33. Senn, Walter, Henry Markram, and Misha Tsodyks. "An algorithm for modifying neurotransmitter release probability based on pre-and postsynaptic spike timing." Neural Computation 13.1 (2001): 35-67. 34. Strain, Thomas J., et al. "A supervised STDP based training algorithm with dynamic threshold neurons." Neural Networks, 2006. IJCNN'06. International Joint Conference on. IEEE, 2006. 35. Bohte, S. M., and M. C. Mozer. "Reducing spike train variability: A computational theory of spike-timing dependent plasticity." Proceedings of Advances in Neural information processing (NIPS’05): 270-279. 36. Masquelier, Timothée, Rudy Guyonneau, and Simon J. Thorpe. "Spike timing dependent plasticity finds the start of repeating patterns in continuous spike trains." PLoS one 3.1 (2008): e1377. 37. Izhikevich, Eugene M., and Niraj S. Desai. "Relating stdp to bcm." Neural computation 15.7 (2003): 15111523. 38. Bush, Daniel, et al. "Reconciling the STDP and BCM models of synaptic plasticity in a spiking recurrent neural network." Neural computation 22.8 (2010): 2059-2085. 39. Burkitt, Anthony N., Hamish Meffin, and David B. Grayden. "Spike-timing-dependent plasticity: the relationship to rate-based learning for models with weight dynamics determined by a stable fixed point." Neural Computation 16.5 (2004): 885-940. 40. Pfister, Jean-Pascal, et al. "Optimal spike-timing-dependent plasticity for precise action potential firing in supervised learning." Neural computation 18.6 (2006): 1318-1348. 41. Fiete, Ila Rani. "Learning and coding in biological neural networks". Diss. Harvard University Cambridge, Massachusetts, 2003. 42. Xie, Xiaohui, and H. Sebastian Seung. "Learning in neural networks by reinforcement of irregular spiking." Physical Review E 69.4 (2004): 041909. 43. Williams, Ronald J. "Simple statistical gradient-following algorithms for connectionist reinforcement learning." Machine learning 8.3-4 (1992): 229-256. 44. Fiete, Ila R., and H. Sebastian Seung. "Gradient learning in spiking neural networks by dynamic perturbation of conductances." Physical review letters97.4 (2006): 048104. 45. Legenstein, Robert, Dejan Pecevski, and Wolfgang Maass. "A learning theory for reward-modulated spiketiming-dependent plasticity with application to biofeedback." PLoS Computational Biology 4.10 (2008): e1000180. 46. Lee, Kyoobin, and Dong-Soo Kwon. "Synaptic plasticity model of a spiking neural network for reinforcement learning." Neurocomputing 71.13 (2008): 3037-3043. 47. Hinton, Geoffrey E., Simon Osindero, and Yee-Whye Teh. "A fast learning algorithm for deep belief nets." Neural computation 18.7 (2006): 1527-1554. 48. Belatreche, Ammar, Liam P. Maguire, and Martin McGinnity. "Advances in design and application of spiking neural networks." Soft Computing 11.3 (2007): 239-248. 49. Jin, Yaochu, Ruojing Wen, and Bernhard Sendhoff. "Evolutionary multi-objective optimization of spiking neural networks." Artificial Neural Networks–ICANN 2007. Springer Berlin Heidelberg, 2007. 370-379. 50. Deco, Gustavo, and Bernd Schürmann. "Spatiotemporal coding in the cortex: information flow-based learning in spiking neural networks." Neural computation11.4 (1999): 919-934. 51. Gupta, Ankur, and Lyle N. Long. "Hebbian learning with winner take all for spiking neural networks." Neural Networks, 2009. IJCNN 2009. International Joint Conference on. IEEE, 2009. 52. Bartlett, Peter L., and Jonathan Baxter. "Infinite-horizon policy-gradient estimation." arXiv preprint arXiv:1106.0665 (2011). 53. Ghaderi, Viviane S., et al. "Analog low-power hardware implementation of a Laguerre-Volterra model of intracellular subthreshold neuronal activity."Engineering in Medicine and Biology Society (EMBC), 2012 Annual International Conference of the IEEE. IEEE, 2012.
31
S0893-6080(14)00085-9 http://dx.doi.org/10.1016/j.neunet.2014.04.001 NN 3318
To appear in:
Neural Networks
Received date: 7 January 2014 Revised date: 10 April 2014 Accepted date: 20 April 2014 Please cite this article as: Yousefi, A., Dibazar, A. A., & Berger, T. W. Synaptic dynamics: Linear model and adaptation algorithm. Neural Networks (2014), http://dx.doi.org/10.1016/j.neunet.2014.04.001 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
*Title Page (With all author details listed)
SYNAPTIC DYNAMICS: LINEAR MODEL AND ADAPTATION ALGORITHM Authors: 1. Ali Yousefi PhD Candidate, Department of Electrical Engineering, University of Southern California Address: DRB 140, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9359 Email: [email protected] 2. Alireza A. Dibazar Research Assistant Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 384, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9315 Email: [email protected] 3. Theodore W. Berger Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 140, 1042 Downey Way, Los Angeles, California, 90089-1111 Phone: 213-740-9315 Email: [email protected]
Corresponding Author: Alireza A. Dibazar Research Assistant Professor, Department of Biomedical Engineering, University of Southern California Address: DRB 384, 1042 Downey Way, Los Angeles, California, 90089 Phone: 213-740-9315 Email: [email protected]
Address: Alireza A. Dibazar 1042 Downey Way, DRB 140, Los Angeles, California, 90089
*Manuscript Click here to view linked References
SYNAPTIC DYNAMICS: LINEAR MODEL AND ADAPTATION ALGORITHM Ali Yousefi, Alireza A. Dibazar, Theodore W. Berger Abstract—In this research, temporal processing in brain neural circuitries is addressed by a dynamic model of synaptic connections in which the synapse model accounts for both pre- and post-synaptic processes determining its temporal dynamics and strength. Neurons, which are excited by the post-synaptic potentials of hundred of the synapses, build the computational engine capable of processing dynamic neural stimuli. Temporal dynamics in neural models with dynamic synapses will be analyzed, and learning algorithms for synaptic adaptation of neural networks with hundreds of synaptic connections are proposed. The paper starts by introducing a linear approximate model for the temporal dynamics of synaptic transmission. The proposed linear model substantially simplifies the analysis and training of spiking neural networks. Furthermore, it is capable of replicating the synaptic response of the non-linear facilitation-depression model with an accuracy better than 92.5%. In the second part of the paper, a supervised spike-in-spike-out learning rule for synaptic adaptation in dynamic synapse neural networks (DSNN) is proposed. The proposed learning rule is a biologically plausible process, and it is capable of simultaneously adjusting both pre- and post-synaptic components of individual synapses. The last section of the paper starts with presenting the rigorous analysis of the learning algorithm in a system identification task with hundreds of synaptic connections which confirms the learning algorithm's accuracy, repeatability and scalability. The DSNN is utilized to predict spiking activity of cortical neurons and pattern recognition tasks. The DSNN model is demonstrated to be a generative model capable of producing different cortical neuron spiking patterns and CA1 Pyramidal neurons recordings. A single-layer DSNN classifier on a benchmark pattern recognition task outperforms a 2-Layer Neural Network and GMM classifiers while having fewer numbers of free parameters and decides with a shorter observation of data. DSNN performance in the benchmark pattern recognition problem shows 96.7% accuracy in classifying three classes of spiking activity. I.
INTRODUCTION
Synaptic transmission is a dynamic process; it is regulated by a variety of short-term and long-term processes in which interplay of the processes defines the synaptic response (Liaw & Berger, 1996; Liaw & Berger, 1999; Dittman, Kreitzer, & Regehr, 2000; Zucker & Regehr, 2002; Tsodyks, Pawelzik, & Markram, 1998; Abbott, Varela, Sen, & Nelson, 1997; Mongillo, Barak, & Tsodyks, 2008; Maass & Zador, 1999). The synaptic response can change on a short time scale by an order of magnitude depending on the past pre-synaptic activity where its temporal dynamics forms the basis of neural computation (Liaw & Berger, 1996; Maass & Zador, 1999). Synaptic facilitation and depression are the core components of any mathematical model describing neurotransmitter release/recovery kinetics for the chemical synapses. In fact, synaptic models are mainly implemented through static weights and axonal delays failing to address fine spatio-temporal dynamics observed in spiking neural networks (SNNs). In this research, the facilitation depression (FD) model of the synaptic transmission is modified, and an approximate linear state space representation of the synaptic transmission is introduced. The synaptic response is modeled by a dynamic process addressing both STP and LTP characteristics where the mutual regulation of STP and LTP factors shape synapse temporal dynamics. There is biological evidence that in almost any neural system which requires a high processing-speed, the timing of a single spike is remarkably precise and reliable (Bohte, 2003; Mainen & Sejnowski, 1995; Berry, Warland, & Meister, 1997). Using simulation analysis, it was shown that neural spike trains have a fine temporal structure with up to one millisecond resolution in which the reliability of spike train is greater than the Renewal process (Hatsopoulos, Geman, Amarasingham, & Bienenstock, 2003; Tiesinga, Fellous, & Sejnowski, 2002). From the collected evidence, there is a well-established fact that precise and reliable spike time coding is a key component for neural processing and information conveyance in the brain. All of the aforementioned facts are in accordance with the
1
synaptic temporal dynamics; synapses are provoked by spiking activity and builds a spike-domain processing unit. In this research, a new family of artificial neural models with hundreds of synaptic connections is built which has the necessary components of a fast spike-domain computational engine. The synaptic plasticity has been accepted for over 60 years. Synapses with repeated use demonstrate different forms of plasticity. At some synapses the plasticity is in the form of facilitation dominance, and in others, the result is a decrement in synaptic strength and depression. Multiple short-term and long-term processes are present in synaptic plasticity, and its efficacy emerges by a combination of facilitation and depression processes which is highly dependent upon the timing of incoming action potentials (AP). In other words, synaptic plasticity is not merely the synaptic strength regulation; it is a dynamic process built upon multiple processes occurring in both the pre- and post-synaptic terminals of a synapse. Biologically realistic learning algorithms have emerged in parallel with experimental evidence of synaptic adaptation to address the learning mechanism in synapses and eventually artificial models of neural circuitries. STDP and its variants introduce the spike-time formulation of the Hebbian learning rule which are biologically plausible frameworks of synaptic adaptation (Song, Miller, & Abbott, 2000; Strain, Daid, Maguire, & Ginnity, 2006; Pfister, Toyoizumi, Barber & Gerstner, 2006; Bohte & Mozer, 2004; Masquelier, Guyonneau, & Thorpe, 2008; Senn, Markram, & Tsodyks, 2000; Morrison, Diesmann, & Gerstner, 2008; Izhikevich &. Desai, 2003; Bush, Philippides, Husbands, & Oshea, 2010; Burkitt, Meffin, & Grayden, 2004). Additionally, machine learning techniques including supervised (Bohte, 2003; Bohte, Poutr, & Kok, 2002; Boij, 2004; Davis, Erdogmus, Rao, & Principe, 2003; Gutig & Sompolinsky, 2006; Ghosh-Dastidar & Adeli, 2009; Ponulak & Kasinski, 2010; Fang, Wang, & He, 2010; Rostro, Cessac, Vasquez, & Vieville, 2009) and statistical (Fiete, 2003; Xie & Seung, 2004;Fiete & Seung, 2006; Williams, 1992; Legenstein, Pecevski, & Maass, 2008; Lee & Kwon, 2008; Deco & Schurmann,1999; Hinton, 2006) approaches are widely utilized for synaptic adaptation. Combinatory techniques plus heuristic search have also been extensively applied in learning of complex dynamic models (Belatreche, Maguire, & McGinnity, 2007; Jin, Wen, & Sendhoff, 2007). Reinforcement learning is the other biologically plausible synaptic adaptation framework which is categorized as a statistical learning approach. Except for one supervised learning algorithm which performs synaptic facilitation adaptation (Storck, Jakel, & Deco, 2001), other supervised learning algorithms, proposed thus far, fail to adjust the temporal dynamics of synapses. Though the learning rule proposed by Storck et al. demonstrates facilitation adaptation in DSNNs, its capability is limited to generating a single spike in the output neuron. Gupta & Long (Gupta & Long, 2009) argue that the STDP algorithm is the optimal learning approach for SNNs, even though it has worked well on a number of applications. One issue with STDP is noncausality of the learning process; when a post-synaptic neuron fires in a time-marching code, it is unknown whether one of the pre-synaptic neurons will fire in the near future. It is also unclear how a STDP rule operating in a time scale of milliseconds can justify most of the behavioral/learning events occurring at much higher time scale. Evidence for triplet and quadruplet STDP rules implies that paired-STDP is not a generalized learning rule for synaptic adaptation of neural circuitries. Statistical learning algorithms and specifically reinforcement learning suffer from slow convergence in large size networks. In reinforcement learning and other sampling methods, synaptic parameters of the neural system are updated with a noisy estimate of the expected reward's gradient. Variance of this estimate diverges with the increase in the network size. In fact, the correlation of reward with the fluctuations in spiking activity of a single neuron becomes weaker as the network size grows. Statistical learning maximizes the probability of a specific activity pattern rather than its precise timing. This is in contrast with the biological finding of precise and reliable spike timing. Heuristic search algorithms fail to demonstrate a learning mechanism in neural circuitries. Moreover, the processing time in finding optimum synaptic parameters increases exponentially as the network size and processing time increases.
2
In this research, a family of learning algorithms is introduced satisfying neuro-biological constraints including simultaneous pre- and post-synaptic adaptation and precise spike timing. A supervised learning algorithm capable of simultaneously adjusting STP and LTP factors of individual synapse of a DSNN is introduced to perform a spike-in-spike-out functional mapping (Yousefi, Dibazar, & Berger, 2011). The assessment of neural circuitries has generated considerable research interest. One of the frontiers of neuroscience research focuses on replacing a damaged human brain region with a neural prosthesis in which the prosthesis duplicates the functionality of the damaged region. Synaptic plasticity is the core component of any cortical neural circuitries. Reproducing the neural functionality requires a precise artificial model capable of performing the spike-domain processing task. Though functional approximation and detailed compartmental methods have been applied in reproducing neural region functionality, the learning scalability and computational complexity of both approaches limit their practicality in the prosthesis application. The Volterra series and its simplified kernel counterpart encounter an exponential growth of the model parameters in accurate modeling of the neural functionality. Though detailed conductance models demonstrate kinetics of different compartments of a neural system, these models are exhaustive and require large computational resources. Similar to the functional approximation technique, free parameters of the detailed model grow rapidly making their tuning process a complicated task. In this research, an intermediary technique for modeling the dynamics of cortical neural circuitries is proposed; the phenomenological models show a balance of biological feasibility and computation complexity. Moreover, structure of the model has the flexibility of expanding in both synaptic connections and neuron layers making the model a superior technique in reproducing the neural system functionality. It will be shown that the phenomenological model of the manuscript along with the proposed learning algorithm are able to predict the spiking activity of Hippocampal neurons significantly better than non-parametric approaches while being able to explain the possible learning mechanisms that occur in the corresponding biological system. The efficiency of the proposed model and training algorithm will be demonstrated in a dynamic pattern recognition task and its performance will be compared with a 2-layer sigmoid neural network. II.
LINEAR MODEL OF SYNAPSE DYNAMICS
The facilitation depression (FD) model of synaptic response (Tsodyks, Pawelzik, & Markram, 1998; Dittman, Kreitzer, & Regehr, 2000; Mongillo, Barak, & Tsodyks, 2008) is described by: ⁄ ⁄ 1 ∗ ∆ ∗
(1.a)
⁄ 1 ⁄ ∗ ∗
(1.b)
⁄ ⁄ ∗ ∗ ∗
(1.c)
⁄ ⁄ ∗ ∗ ∗
(1.d)
(1.e)
Variable describes facilitation dynamics in response to incoming APs, and variable represents the portion of release ready vesicles. Auxiliary variables and model an -synapse function representing the PSP - - elicited due to synaptic neurotransmitter release (Gerstner, 2006). The PSP of a single synapse in equation set (1) is a non-linear function of synapse state variables - ∗ -; and it is the interplay of facilitation and vesicle recovery/release that determines the temporal dynamics of synapse. The ∗ term defines the portion of vesicle sites that release their neurotransmitter due to an impinging AP time while parameter defines the population of neurotransmitter receptor sites. Parameter ∆ controls calcium influx in the presynaptic terminal of a synapse, and it is the main biological factor in modulating synaptic dynamics. Parameter is a post-synaptic factor representing the number of neurotransmitter receptors or simply synaptic strength. Facilitation increase and vesicle release dynamics are both AP-driven process where the AP time is defined by . The FD model proposed in equation (1) defines the fundamental biological components shaping
3
synapse temporal dynamics. Appendix A defines a more intricate synapse model including multiple pre- and post-synaptic facilitation processes in addition to the calcium dependent vesicle recovery time constant proposed by Regehr et al. The discrete counterpart state space model of synapse dynamics is built by utilizing the Euler discretization technique. The discretization resolution ∆ is defined in milliseconds with a nominal value of one millisecond. The millisecond resolution is required to accurately trace the fast kinetics of synaptic facilitation and vesicle release process. Facilitation increase and vesicle release in response to incoming AP occurs in about the one millisecond time period; thus, in the discrete model both facilitation increase and vesicle release evolve at the same time of impinging APs. Equation set (2) defines the discrete state space model of the synaptic response: !"# $1 ∆ ∗ 1 ' ( ∗ ! ∗ ∆ ∗ 1 '
% % ∆
∆
!"# $1
&
∆ %)
&
* ∗ ∆+ ∗ 1 ' ( ∗ ! *1 ∆ %&
∆ %&
(2.a)
∆+ ∗ 1 ' ∗ ! ∗ !
!"# *1 + ∗ ! ∗ * ∗ ∆+ ∗ 1 ' ∗ ! ∗ *1 ∆
%,
∆ %&
!"# -1 / ∗ ! ∗ * ∗ ∆+ ∗ 1 ' ∗ ! ∗ *1
! ! !
∆
%.
∆ %&
∆ %&
∆ %&
∆ %)
∆+ ∗ 1 ' ∗ ! ∗ !
∆+ ∗ 1 ' ∗ ! ∗ !
(2.b) (2.c) (2.d)
Variable 1 ' denotes AP-time, and it is equal to one at spike time. The state variables of equation set (2) define discrete samples of the continuous FD model state variables. Table 1 defines each parameter of the model and its nominal value. Figure (1) shows state variables of the FD model in response to a sample spike train. (2.e)
Table 1 Range of parameters in synapse FD model Parameters ∆ ∆
Facilitation time constant Resting facilitation Facilitation increment factor Vesicle recovery time constant Maximum number of release sites Rise time of PSP Rise time of PSP Processing time resolution
Value 100…500 milliseconds 0…1 0…1 100…1500 milliseconds 1...20 3…15 milliseconds 3…15 milliseconds 1 millisecond
The synapse model defined in equations (1) and (2) are a non-linear system, and it is the AP occurrence that injects non-linearity to the synaptic response (Morrison, Diesmann, & Gerstner, 2008; Yousefi, Dibazar, & Berger, 2011). Non-linearity in the synapse model appears in the ! ∗ ! term, and it becomes linear if the ! ∗ ! term is replaced by a surrogate variable. The new variable 0! ! ∗ ! is defined for the linearization process. The variable denotes the minimum vesicle release for a possible AP at time ' 1. To build the synapse linear state space model, 0!"# needs to be defined by a linear combination of the ! , ! and 0! state variables. Variable 0!"# !"# ∗ !"# can be rewritten as a linear function of ! , ! , 0! plus a new non-linear term - ! 1 ∗ ! . 0!"# !"# ∗ !"# 2# ∗ ! 21 ∗ 3# ∗ ! 31 ∗ ! ∗ ! 34
2# ∗ 3# ∗ ! ∗ ! 2# ∗ 31 ∗ ! ∗ ! ∗ ! 2# ∗ 34 ∗ ! 21 ∗ 3# ∗ ! 21 ∗ 31 ∗ ! ∗ ! 21 ∗ 34
2# ∗ 3# 21 ∗ 31 ∗ 0! 2# ∗ 31 ∗ ! ∗ ! ∗ ! 2# ∗ 34 ∗ ! 21 ∗ 3# ∗ ! 21 ∗ 34
4
(3)
0.5
100
200
300
400
0
500
1
0.5
0.5
200
300
400
0
500 G
2 0
K
100
100
200
300
400
2 0
100
200
300
400
1 0.5 0
100
200
300
400
1.5 1 0.5 0
AP
1 0.5 0
100
200 300 Time (msec)
400
500
(a)
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200
300
400
500
100
200 300 Time (msec)
400
500
1 0.5
500
200
2 0
500
100
2 0
500 K
G
0
N
1
Y /Y m a x
N
F
0.5 0
Y /Y m a x
∆F = 0.8, Nmax = 1.0 1
AP
F
∆F = 0.1, Nmax = 1.0 1
0
(b)
Figure 1 Synapse FD model response to a 500 millisecond spike train with 50 milliseconds ISI. a) Synapse with ∆F=0.1, b) Synapse with ∆F=0.8
Parameter 2# , 21 , 3# , 31 and 34 are defined using equation set (2). The linearization process can generate higher order terms of ! 6 ∗ ! and ! 6 , 8 1 whose values shrink to zero as increases. This is because both ! and ! vraibles are normalized positive terms with a maximum value of 1.0. The linearization process is completed by truncating higher-order non-linear terms of the synapse FD model, and then defining a set of new surrogate variables - along with the recursive equations for updating each new variable. For the first order linear model, the only non-linear term is ! ∗ ! and other terms of ! 6 , 8 1 is set to zero. In the second order linear model, non-linear terms include ! ∗ ! , ! 1 ∗ ! and ! 1 and other extra terms containing ! 6 , 8 2 is set to zero. The same procedure is applied for a higher order linear model in which any non-linear term consisting of ! 6 , 8 : is set to zero; : defines the order of linear model. Figure (2) shows the normalized error of vesicle release for different linearization orders - : ∈ - within the possible range of interspike interval (ISI) and ∆. Independent of linearization order and facilitation factor, the maximum vesicle release error increases with smaller ISIs. The vesicle release error decreases as the linearization order increases. Meanwhile the maximum of vesicle error is delayed to later APs of the impinging spike train. Similar to the ISI scenario, the normalized error gets larger for lower facilitation factors. Synapses with a lower facilitation factor correspond to potentiating synapses. Therefore linearization error accumulates through multiple impinging APs. At lower ISIs, the synapse vesicle pool is depleted completely after ten to fifteen
5
impinging spikes due to the lack of recovery reco time between incoming APs. Thus, hus, the maximum error is more reasonable until til a specific number of consecutive APs. The simulation results show that the maximum vesicle release error for the linear model with a third-order order approximation is below 15% following a prolonged spiking activity. For spiking activity with a maximum of ten consecutive impinging spikes, the maximum vesicle release error drops below 7.5%. The simulation resultss show that for spiking activity with ISIs of more than 100 milliseconds, the first-order first order linear approximate is an accurate model for synapse dynamics. The results result confirm that the third-order order linear approximate state space model is an accurate model of synapse synapse dynamics independent of ISI variation and the facilitation factor. In the third-order order linear model, other surrogate variables of the model - besides the 0! variable- are: ?! ! 1 ∗ ! , @! ! 1 , P! ! 4 ∗ ! and H! ! 4 . Higher igher order term terms of ! 6 ∗ ! and ! 6 , 8 3 are set to zero. Equation (4)) presents synapse linear state space model: ∗ D A!"# BC 1 ' ∗ B ∗ A! DC 1 '
! E0 0 0 0 0 0 0 1 1G ∗ A! Q ∗ A!
(4.a)
A! E! @! ! 0! ?! H! P! ! ! is the model state vector, and ! corresponds to the PSP. Elements of BC , B , DC and D are constant functions of synapse internal parameters, and 1 ' is one at AP time. Equation (4) ( introduces a generalized linear model for synaptic response response. The model also provides a compact tool for temporal and steady state analysis of the synaptic response. The following subsections analyze synapse temporal dynamics as a function of synapse internal parameters and AP variability. GR
(a)
(4.b)
(b)
Figure 2 Maximum vesicle release error in the first, second and third order linear approximate model of synapse FD model model. The error is equal to normalized vesicle release between the linear and non-linear non models defined by OOOOOOOOOOOOO IJJ ∗ |LM"I ∗ NM L M"I ∗ NM |⁄LM"I ∗ NM ; in which LM"I ∗ NM determines vesicle release at impinging spike time, and OOOOOOOOOOOOO L M"I ∗ NM is the linear estimate of vesicle release using the approximate linear model. a) Maximum vesicle error for different ISIs. For each specific ISI, the parameter space - facilitation factor and AP index - is searched for the maximum error, b) Maximum vesicle error for different facilitation tation factors. factor Number of APs is limited to 10 consecutive spikes. Synapse time constants - ST , SU - are set at IVJ, WVJ milliseconds.
6
A. Parameter Sensitivity Analysis
Each element of the state matrix - BC , B - and input vector - DC , D - is a differentiable function of the synapse internall parameters. Equation (5) (5 denotes the derivative of a state vector relative to any of the synapse parameters: XYZ[ X\
X`Z X\
-
X]^
Q∗
X\
1 ' ∗
XYZ X\
X]_ X\
/ ∗ A! BC 1 ' ∗ B ∗
XYZ X\
-
Xf^ X\
1 ' ∗
Xf_ X\
/
(5.a) (5.b)
(b)
(a)
Figure 3 Synapse temporal dynamics. a) Synaptic PSP response for ∆L J. JV, J. h and J. i plus PSP sensitivity in response to ∆L changes. Input APs has ISI of 50 milliseconds illiseconds and last for 400 milliseconds. Smaller ∆L leads to a potentiating synapse. Higher values of ∆L turns to a depressing synapse, b) Synaptic temporal selectivity in response to ISI variation for different ∆L. Smaller ∆L has a higher sensitivity for lower ISIs, ISI and the peak of synaptic efficacy move to lower ISI for higher values of ∆L. Synapse time constants - ST , SU - are set at IVJ, WVJ milliseconds.
Variable a can be any synapse parameter; in figure (3.a) .a) the PSP derivative with respect to facilitation factor - ∆ - at different values is shown. The PSP derivative with respect to ∆ is not a monotonic function, and its sign and amplitude changes depending on AP time and internal synapse parameters. The PSP derivative with respect to is equal to a normalized PSP. Thus, the PSP increases with the growth of . Equation (5)) will be used in the learning algorithm; the equation defines a recursive function to calculate the derivative of the neuron membrane potential relative to the synapse parameters. B. AP Pattern Sensitivity Analysis
Synaptic efficacy is defined by the ratio rat of vesicle release in the jk and the first impinging AP time. This ratio for different ISIs builds a quantitative measure of the synaptic temporal dynamic dynamics in response to the variability of the impinging spike train. Equation (6) denotes the efficacy ratio for a spike train with b-millisecond ISIs. Figure (3.b) shows the ratio for different ISIs and different facilitation factors. c A∗ "# 4 ⁄A "# 4 1 1 j ∗ b 1
(6)
AO!"# BC e ∗ B ∗ AO! DC e ∗ D
(7.a)
Synaptic efficacy is one possible measure to quantify synaptic temporal dynamics. s. It is also possible to identify synapse temporal dynamics using its maximum release lease in response to impinging spikes spikes. Similarly, synaptic transmission in response to the mean firing rate of incoming incoming spike trains can be defined using the linear model.. For a homogenous firing rate - e -, equation set (7) represents the synaptic response and the state variables' derivatives with respect to change in firing rate at its steady state.
7
XYOZ[ Xl
BC e ∗ B ∗
XYOZ Xl
B ∗ AO! D
(7.b)
AOmnopqr s BC e ∗ B t# ∗ DC e ∗ D
(7.c)
uOmnopqr Q ∗ AOmnopqr
(7.d)
Temporal dynamics of the synaptic response is addressed in equations (4-7). The next section extends the linear state representation to a neuron model consisting of multiple synaptic connections. C. Neuron Membrane Dynamics Model The neuron model integrates the PSP response of multiple synaptic connections plus the refractory process emerging from its previous firing activity. The neuron model is capable of receiving external stimuli in which both the external stimuli and refractory process are defined by a linear state space representation. Equation set (8) defines neuron membrane potential dynamics and its firing activity criterion: ' ∑ ' Ref z ' {o '
(8.a)
{| ' 1 1 ∆⁄ ∗ {| ' ∗ ∆ ∗ s| '
(8.b)
Ref z ' 1 1 ∆⁄ | ∗ Ref z ' Refp}~ ∗ 1 '
If ' 8 k ∩ ' ' 1 ⇒
1 '
(8.c)
1
(8.d)
Variable ' is neuron membrane potential at time ', and variable { | ' is the impact of external stimuli - s | ' - on the neuron membrane potential. Variable ' is the PSP of the jk synapse in response to incoming APs, and variable Ref z ' corresponds to the refractory effect in response to neuron spiking activity - 1 ' . Equation (8.d) defines the firing criteria of the neuron; the neuron fires an AP whenever its membrane potential passes the firing threshold. Parameter and are defined by neuron membrane characteristics, and | and Refp}~ - a negative term - model refractory dynamics. The neuron membrane dynamics follows the same generalized form proposed for synaptic response by equation (4) in which 1 ' is replaced by a matrix replicating all parent’s plus neuron's spiking activity at time '. Figure (4) shows the motif of neuron models with multiple synaptic connections. For a Type I neuron model, the network consists of feed-forward synapses projecting input spiking activities on the neuron membrane while the Type II neuron has an extra synapse projecting its firing activity to itself. The feedback synapse has different dynamics compared to the refractory process. The feedback synapse can be either excitatory or inhibitory with a dynamic response. Similarly, feed-forward synaptic connections can be excitatory or inhibitory. The dynamic synapse neural model and its morphology Type I and Type II - are defined; the next section will focus on the learning algorithm for the DSNN. Ref
Ref
(a)
(b)
Figure 4 Type I and Type II DSNN. a) Type I DSNN with multiple feed-forward synaptic connections, b) Type II DSNN with multiple feed-forward synaptic connections plus a feed-back synaptic connection. Feed-back synaptic connection has different dynamics comparing to the refractory process, and it can be either excitatory or inhibitory.
III.
SYNAPTIC ADAPTATION ALGORITHM
A biologically plausible learning rule capable of simultaneously adjusting long term plasticity - LTP and short term plasticity - STP - parameters of individual synaptic connections of a Type I/II DSNN is introduced. Facilitation factor - ∆ - and synaptic strength factor - - are adapted simultaneously
8
to shape desired neuron spiking activity. The learning rule is a supervised learning rule in which the similarity between a neuron’s spiking activity and a desired spike train is increased. Moreover, a new STDP learning rule capable of training Type I/II DSNNs with hundreds of synaptic connections is introduced. In the following sections, a similarity measure between spike trains and supervised learning rules are introduced. The new STDP and its formulation will be described in detail. A. Spike Train Similarity Measure The spike train similarity measure quantifies the coincidence between APs in the test and the desired spike trains (Kreuz, Haas, Morelli, Abarbanel, & Politi, 2007; Yousefi, Dibazar, & Berger, 2011). The proposed similarity measure quantifies the similarity between two spike trains, the test and desired spike, in two different time scales. The fine temporal scale identifies each AP of the test spike train as a: 1) similar, 2) missing, or 3) extra AP. Meanwhile, the global scale measure returns a normalized value determining the overall similarity between two spike trains. An AP is considered similar if it happens in the time vicinity of one of the APs in the desired spike train. This vicinity time interval is called the similarity range and determines acceptable jitter between spike pairs in the test and desired spike trains. It is defined to be less than half of refractory period and is generally limited to a couple of milliseconds. The overall similarity of two spike trains is defined by: C ⁄j 1, ,
C is the number of similar spikes; and denote the number of spikes in the desired and test spike trains respectively. The similarity measure - - returns a positive normalized value with a maximum of one corresponding to completely similar spike trains. The measure can be defined at different time scales, and as the time scale increases temporal resolution of the measure drops. Figure (5) shows similar, missing and extra spikes for a sample pair of test and desired spike trains in which gets value of 0.6. Similarity between two spike trains is increased by reducing the number of missing and extra spikes. To increase the similarity, the membrane potential at extra spikes needs to be dropped below firing threshold, and at the same time, the membrane potential of missing spikes needs to be increased above the firing threshold. The learning algorithm adjusts membrane potential for increasing the overall similarity between two spike trains. The following section defines the learning algorithm in DSNN. (9)
Target B'
Test D'
s
m
s
e
m
s
'
j
AP Refractory Period Similarity Range Similar Extra Missing
Figure 5 Similarity measure. Similarity measure classifies APs of a test spike train to three classes of : 1) Similar (s), 2) missing (m), and 3) extra (e) spikes. ' is the relaxed spike train defined by the comparison between B' and D' spike trains. The S value is equal to 0.6.
B. Off-line Learning Algorithm The learning algorithm modifies synaptic parameters to maximize the similarity measure between the training and desired spike trains. AP firing is a discontinuous function of synapse membrane potential, and it is approximated with a continuous function to build a tractable learning algorithm. The learning objective function approximates the similarity measure; it is defined by: E ∑ 1qp~ n ∗ Ep~ n 1 1qp~ n
∗ Eztp~ n
Ep~ n log -1 exp -α ∗ Vn Vn // log *1 exp -β ∗ Vn 1 Vn /+
Eztp~ n log1 expα ∗ Vn Vn ∗ log *1 exp -β ∗ Vn Vn 1 /+
9
(10.a) (10.b) (10.c)
Ep~ n and Eztp~ n determine the network membrane potential deviation from its desired behavior, and 1qp~ n is either zero or one depending on the desired spiking activity, or the similarity measure output, at processing time '. Both Ep~ n and Eztp~ n return a small positive, almost zero, value whenever the training network neuron membrane potential satisfies firing/non-firing criteria. The objective function E returns a small positive value at similar spike times, and its returned value grows as the number of dissimilar spike grows. It is theoretically shown - in Appendix B - that the expected value of c is proportional to the number of missing and extra spikes. Thus, reducing the objective function value yields a higher similarity measure for the training network spiking activity. In the objective function definition, free parameters α and β are set to positive large values that penalize neuron membrane potential deviation from its desired value. The optimum direction for synaptic parameter adjustment is defined by the learning objective function gradient. Synaptic parameters are adjusted in the gradient descent direction to increase the similarity between the test and the desired spike train. The learning algorithm is defined by: a a ∗ ~ c
(11.a)
P maxLower Bound, minUpper Bound, P η ∗ ~ E
(11.b)
Parameter ∆ is a positive number with a maximum value of one, and similarly is constrained to be a positive number. Excitatory or inhibitory synapses are defined independent of the sign of . Thus parameter adjustment - equation (11.b) - is defined by considering parameters' range and sign. The learning objective function returns a value close to zero for processing bins with similar spiking activities. Thus the objective function gradient is defined by gradient of either Ep~ n at missing spikes or Eztp~ n at extra spike time. In fact, the gradient of the learning objective function is defined by the sum of sparse samples of the neuron membrane potential gradient at missing and extra spike time. It is theoretically shown - in Appendix B - that the membrane potential of Type I/II DSNNs is a unique function of its synaptic parameters. The uniqueness is extended to spiking activity in Type I/II DSNNs in which the free synaptic parameters of two DSNN networks with the same number of synaptic connections and similar spiking activity converge to each other. The gradient descent learning algorithm reduces the returned value of the learning objective function, and this reduction corresponds to increases of the similarity measure. Thus, the learning algorithm converges to optimum parameter adjustment. Algorithm I describes processing steps for the learning algorithm in Type I/II DSNNs. Algorithm I Supervised Off-line Learning Algorithm in Type I/II DSNN
1. 2. 3. 4. 5.
6. 7.
Initialize synaptic parameters of the training network with random values Set the desired similarity range - Set the learning iterations – L Set j to zero – there are M set of input and output spike pairs. ¤¥' - and desired spiking activity - B' . ' ∈ For jk pair of input - 1 ¤¥, 0 a. Set § 0, D' 0 b. If § 8 @ ∪ 8 → Jump to step 6 c. ª c 0, § § 1 ¤¥' input stimuli (Equation 8.b) d. Run the network with 1 e. Calculate and build ' - ' is built by matching B' - network spiking activity- with B'
f. Calculate Ref z using ' (Equation 8.c) - Refractory compensation g. Calculate output neuron membrane potential dynamics (Equation 8.a) h. Calculate optimum parameter update (Equation 11.a) i. Update network parameters (Equation 11.b) j. Jump to step b m m 1, if j ¬ H jump to step 5 Stop learning
10
The similarity signal - ' - generates a relaxed desired spike train emerging from the spike train similarity process. Figure (5) shows the ' signal for a sample pair of training and desired spike trains. The relaxed spike train - ' - will be utilized in refractory process compensation. For a Type II DSNN, this will also determine the impinging spike train fed into the feedback synaptic connection. The refractory process of a training network initiates at each spike of the output neuron. Any mismatch in the neuron spiking activity changes the dynamics of the refractory process. Using the missing and extra spikes of similarity signal, the refractory process of a training network is compensated to replicate the desired refractory process. The feedback synaptic connection of the Type II DSNN can be broken by replacing its impinging spiking activity with the desired spiking activity. Replacing '
with the desired spiking activity - B' - makes the learning algorithm insensitive to spiking jitter, and the learning convergence is boosted as both of the desired spike train and synapse parameters are adjusted at each learning step. Practically in Type II DSNNs, the learning algorithm is repeated multiple times for each training sample. Although Algorithm I is an off-line learning algorithm, the optimum parameter adjustment is defined by observing only a sparse set of spiking activity or a partial observation of the neuron membrane potential. Regularity in the parameter update rule at missing and extra spike time suggests the possibility of building a more generalized learning rule which is online and is independent of observing the desired spiking activity. The next section will introduce a gradient-gain learning rule which is an on-line learning algorithm applicable in the synaptic adaptation of Type I/II DSNNs. C. Gradient-Gain Learning Rule The gradient of the post-synaptic membrane potential with respect to synapse parameters is positively gained at missing spike times and negatively gained at extra spike times. The gain is either a constant value, or it reflects the neuron membrane potential deviation from the firing threshold. The gain can be generalized to any function of membrane potential and its firing threshold - including a constant with alternating sign -, and this gain function - defined by ' - adapts the spiking activity of the training network. In the Algorithm I, the gain value during similar spiking activity and any silent period is set to zero, and it is equal to a constant at missing and extra spike time. It is assumed that all synaptic connections are excitatory, and desired spiking activity is generated by the same neural topology. Refractory compensation as defined in the algorithm I will be dropped if the gain function predicts spiking activity - in the next milliseconds where defines the similarity range - of the training network rather than its instantaneous spiking activity. For example, the training algorithm expects to see a spike in next milliseconds and the training network activity is impotent to fire an AP shortly. Thus, a positive gain is injected to the network boosting its membrane potential to fire an AP. If the network responds to the gain injection by firing an AP, the refractory process is already adjusted by the network response. The same scenario happens whenever the network receives a negative gain, and it reduces its membrane activity below the firing threshold. The gradient gain with a predictive behavior resembles the temporal difference - TD - learning algorithm in which the system behavior is adapted to reduce upcoming error. Algorithm II Gradient-gain Learning Algorithm in Type I/ II DSNN 1. 2. 3. 4. 5.
Initialize synaptic parameters of training network with random values Set ' 0 and ª~® 0 ( is the PSP of mth synapse, and p is the synapse free parameter) Build network state space model (Equation set 8) Run the network with ¤¥ 1' input stimuli (Equation set 8) Update network parameters with the received gain - ¯ is the forgetting factor and ' is the gain function a. ª~® ∶ ¯ ∗ ª~® ' ∗ ª~® '
p} ∶ maxLower Bound, minUpper Bound, p} η ∗ ª~® ' ' 1 Jump to step 3
b. 6. 7.
11
The generalized learning rule for a Type I/II DSNN is defined by Algorithm II. The learning rule is defined independent of refractory compensation and uses the gain function to adjust the DSNN synaptic parameters for a desired spike-domain functional mapping. The gain is any evaluative function reflecting the similarity between the training network spiking activity and a desired one. Learning algorithm II is similar to a family of approximate dynamic programming algorithms called Gradient Policy (Williams, 1992). The gain function corresponds to the reward function of the gradient policy algorithm, and the rewarding process measures the training network’s ability to perform a defined dynamic task. Similarly, the gradient-gain can be generalized to any evaluative function of spiking activity of the training network. An evaluative gain function can be defined by the instantaneous spiking activity of the training network, its average spiking activity through time, or any other combination of its temporal response. Thus, the gain function is shaping neural system response for performing a spike-domain functional task. Gradient policy learning algorithms - including REINFORCE and GOMPD (Baxter & Bartlett, 2001; Williams, 1992)- build a stochastic gradient ascent approach to adjust training network policy parameters with a guarantee of convergence. In contrast, the gradient-gain learning proposed here is a deterministic approach in which optimum parameter adjustment is determined by the variability in spiking activity of the training network. Optimum parameter adjustment of the gradient-gain learning rule is the same as in off-line learning; thus, the off-line learning convergence can be extended to this gradient-gain learning. The learning objective function is a convex function with a minimum value of zero. The learning objective function converges to zero as the network converges to its optimum parameters. Optimum parameter update is defined by two terms: 1) the local gradient of each individual synapses, and 2) the gain function. The gradient term is generally non-zero; thus, receiving a non-zero gain implies an incomplete learning process. The gain function definition requires convergence to zero as the training network learns the desired spiking activity in response to incoming stimuli. The similarity signal generated by the spike similarity process returns zero for a similar spike train; thus, the similarity signal is a viable gain function. Either defining the gain function with a stopping mechanism or dropping its time average from the instantaneous gain are possible choices to build a learning stopping criteria. The gradient-gain learning is a function of the pre- and post-synaptic spiking activity of the network, and it can be extended to a STDP learning rule. The next subsection discusses a new STDP rule for synaptic parameter adjustment in DSNNs. D. STDP Learning Rule In the gradient-gain learning, synapse parameters are updated at the time indexes - '± - at which the network receives the gain - '± . ª c'± is a function of pre-synaptic AP time and synapse internal parameters which involve the history of impinging APs for synaptic adaptation. The gain factor is a function evaluating the network post-synaptic spiking activity. Thus, the gradient-gain is a function of pre- and post-synaptic spiking activity of each individual synapse leading to a new class of STDP learning rule. In contrast to pair-wise STDP learning rules (Morrison, Diesmann, & Gerstner, 2008) in which synaptic weight is adjusted by the relative timing of a pair of pre- and post-synaptic AP time, ª c'± at gain time is determined by a history of AP times and synapse internal parameters. For the pair-wise STDP rule - plus their extension of triplet and higher order STDP rule -, synaptic weight are strengthened whenever a pre-synaptic AP precedes post-synaptic AP and weakened when the postsynaptic AP leads a pre-synaptic AP. The pair-wise STDP learning rule presents a uni-modal learning behavior. However, the STDP learning rule emerging from the gradient-gain alternates the synaptic weight strengthening/weakening mechanism - along with STP factor - depending on the gain function. The pair-wise STDP requires a matching process to find corresponding pre- and post-synaptic APs, whereas the new STDP rule is independent of an AP matching process as the gain function is defined independent of pre- and post-synaptic AP activity sequencing. In contrast to the standard STDP learning rule (Morrison, Diesmann, & Gerstner, 2008), the gradient gain STDP rule is a causal learning
12
process. The gain value can be defined efined without a matching process; therefore, the STDP learning rule emerging from gradient gain learning can be an on-line on learning rule. Figure (6)) describes the STDP learning rule defined by gradient-gain gain learning for different sets set of synapse parameters. Appendix C describes the components of the gain-STDP gain learning rule. Through this section, off-line line and on on-line learning rules for synaptic adjustment in DSNNs DSNN are introduced. The STDP plasticity learning rule addresses both STP- and LTP-based synaptic adjustment in DSNNs. The following sections will focus on the computational capability and verification of the proposed learning algorithms for DSNNs.. 1st
2nd
3rd
4th
Incoming Spike Train Gain Function
5th ∆t
Pre-synaptic neuron Post-synaptic synaptic neuron
(a)
(b) ∆ 0.05, 1
(c) ∆ 0.3, 1
(d) ∆ 0.9, 1
Figure 6 Gradient-gain gain learning rule and STDP. A positive gain corresponds to an expected AP at gain injection time, meanwhile the negative gain reflects an unexpected AP. STDP graphs shown here define STP and LTP adjustment for different facilitation factors with a fixed positive gain injected inject with different delays. ISI of impinging AP is 50 millisecondss and it includes 5 consecutive APs. a) Schematic of gradient-gain gradient learning rule, b) STP and LTP learning curve for a potentiating synapse, c) STP and LTP learning curve for a synapse with medial facilitation factor, d) STP and LTP learning curve for a depressing ressing synapse. The learning curve is a non-linear non linear function of facilitation factor and impinging AP's ISI. Although the learning curve shown here is only limited to 5 APs, practically larger number of APs are involved in synaptic adaptation. Synapse time constants - ST , SU - are set to IVJ, WVJ milliseconds.
IV.
SIMULATION RESULT
In this section, different aspects of DSNN networks and the proposed learning learning rule will be addressed through a set of benchmark problems plus replicating in vivo neural responses. First, a system identification analysis is performed on Type I/II DSNNs. Type I/II DSNNs with different number numbers of synaptic connections are trained by observing the input stimuli and the output spike train of a desired DSNN network having the same me topology. Algorithm I’s convergence and accuracy in parameter estimation in response to network size growth and similarity range variation will be studied in detail. The second set of results examines the FD model’s model capacity to predict the temporal dynamics mics of the
13
Schaffer collateral synaptic pathway which projects CA3 neuron spiking activity onto the CA1 pyramidal neurons. Free parameters of the synapse model are tuned to replicate the CA1 post-synaptic membrane potential. Furthermore, Algorithm II’s ability to estimate the optimum synaptic parameters will be demonstrated. Next, a generative model - a Type II DSNN with one feed-forward synaptic connection - with dynamic synapses is introduced to produce the different spiking and burst patterns of various cortical neurons. The gain function for the on-line learning process is defined, and capability of gradient-gain learning algorithm in generating different spiking patterns will be analyzed. For the last portion, a spatio-temporal pattern recognizer is built utilizing the DSNN model. The classifier consists of three output neurons recognizing three different classes of spiking stimuli in which each output neuron matches with one of the input classes. The network is trained by the gradient-gain learning rule to fire APs with a shorter onset delay - or a higher number of APs - for the neuron that corresponds to the class to which the incoming spike train belongs. In the system identification problem, the gain is an estimate of the membrane error between the test and target neurons. In contrast, synaptic adaptation for the generative model utilizes the spike time difference for the gain function. The gain function for the spatio-temporal pattern classifier and the CA1 neuron model switches between constant positive and negative values at the desired and undesired spiking times to regulate the firing time/rate of the output neurons. This set of problems demonstrate versatility of the gain function and how the spatio-temporal processing capability of synaptic responses can be merged to build a powerful biological computational unit. In the following subsections, the different simulation problems are discussed in detail. A. System Identification in DSNN
Algorithm I is applied for synaptic parameter estimation - Δ and - for a DSNN network by observing its input and output spiking activity. Four different properties of the proposed learning algorithm will be analyzed: 1) convergence, 2) repeatability, 3) scalability, and 4) similarity range variation. The learning properties are quantified by their convergence to unknown synaptic parameters, similarity measure improvement, and the square error in membrane potential prediction. The components and structure of simulation test are as follows: a. Both the training and target network consists of Nm excitatory synaptic connections. b. The facilitation factor for each synapse is generated independently by a uniform distribution within the range of E0.05 0.95G. c. The synaptic strength for each synapse is defined by uniform distribution in the range of E0.05,0.05 K ∗ UE0,1 G. K is dependent on the number of synaptic connections. Its minimum value is 0.2, and its nominal value is defined as 20⁄Nm . UE0,1 is a uniform distribution ranging from 0.0 to 1.0. The synaptic strength is chosen by the assumption that six synapses with the same facilitation factor of 0.5 will satisfy firing criteria at the time of their first impinging spikes. d. The firing threshold is set to be 1.0. e. The refractory process parameters τ·o¸ , Refp}~ are 2.5 msec, 5 . f. The input stimuli to the network consists of an array of Nm spike trains. Each spike train is generated independently by a Renewal process with a random rate uniformly chosen from 5 to 20 Hz. g. There are at most 1500 training spike trains and 100 test spike trains. Training progress is examined on the test data set. h. The simulation time is 400 milliseconds, and the time resolution ∆t is set to 0.2 msec. i. The nominal similarity range is 2 milliseconds unless otherwise mentioned. j. The learning rate is set to be 0.01. k. In Type II DSNNs, the learning process is repeated 5 times per each training sample. For Type I DSNNs, the learning process is applied once per training sample. l. , are set to 150,250 milliseconds.
14
1. Learning Convergence Table 2 and 3 show the learning progress for a Type I and II DSNN with 10 synaptic connections. Results in both tables show a significant improvement in the similarity measure and membrane error reduction. The similarity measure between the desired and training network grows from 0.21 - 4 similar spikes - to 0.97 - completely similar spikes - in the Type I DSNN while the membrane potential error decreases by 34.5 dB. Similarly, the same improvements are observed in the Type II DSNN. The similarity grows from 0.07 - 1 similar spike - to 0.97 - completely similar spikes - while there is more than a 31 dB improvement in the membrane potential error. The result implies an accurate estimation of the network parameters. Figure (7.a) and (7.b) show the initial guesses of the unknown parameters ∆, - and the estimated parameters after training with 1500 samples using a Type I DSNN. Figure (7.c) and (7.d) show the membrane potential of the trained and target networks in response to a sample input using the initial and trained parameters. The numerical result for the sample task verifies convergence to an optimum parameter. The next essential questions are repeatability and scalability of the learning algorithm in synaptic parameter adjustment. The following simulations will focus on Type I DSNN as the learning convergence to optimum parameter in Type II DSNN has been already observed by the convergence simulation result. Table 2 Learning progress for a Type I DSNN network with 10 synaptic connections - in average there are 19 spikes per input sample in 400 milliseconds. Training Sample 1 100 500 1000 1500
Similarity (Mean)
Similarity (Variance)
0.22 0.80 0.93 0.96 0.97
0.007 0.008 0.003 0.002 0.002
Membrane (Mean) 27.06 dB 14.35 dB 2.12 dB -2.47 dB -7.70 dB
Error
Membrane Error (Variance) 0.69 2.98 4.55 2.46 2.30
Table 3 Learning progress for a Type II DSNN network with 11 synaptic connections - in average there are 14 spikes per input sample in 400 milliseconds. Training Sample 1 100 500 1000 1500
Similarity (Mean)
Similarity (Variance)
0.07 0.85 0.95 0.95 0.98
0.005 0.007 0.003 0.003 0.001
Membrane Error (Mean) 29.86 dB 13.32 dB 3.29 dB 5.80 dB -2.36 dB
Membrane Error (Variance) 1.10 4.96 11.12 6.56 12.28
2. Repeatability Analysis For repeatability analysis, the training process is repeated 5 times each time with a new set of parameters. Table 4 shows the minimum and maximum of similarity measure plus the membrane error improvement for the trained network. The results confirm a consistent convergence to an optimum parameter set. Greater accuracy is achievable by increasing number of training samples. The average number of spikes for a trial is about 14 spikes in 400 milliseconds. The results demonstrate that in average 13 out of 14 possible spikes coincide within a 2 milliseconds accuracy to the spike train of the desired DSNN. Similarly, the membrane error measure shows a significant improvement in neuron membrane potential prediction. For the next simulation, scalability of the learning algorithm will be studied. Table 4 Learning progress for a Type I DSNN network with 10 synaptic connections, 1500 training samples for 5 different training/test set Similarity (Min) 0.907
Similarity (Max) 0.971
Membrane Error Drop (Min) 19.01 dB
15
Membrane Error (Max) 34.7 dB
(b) 3
2
2
PSP
PSP
(a) 3
1 0
0
50
100
150
200
250
300
350
1 0
400
2
0
50
100
150
200
250
300
350
400
2
AP
AP
Target Test 1
0
0
50
100
150 200 250 Time(msec)
300
350
1
0
400
0
50
(c)
100
150 200 250 Time(msec)
300
350
400
(d)
Figure 7 Parameter convergence and membrane error in system identification task. task a) Value of ∆, in the target and training networks before training, b) Value of ∆, in the target and training networks after training with 1500 samples, c) Membrane embrane and AP activity of the target and training networks network before training, d) Membrane and AP activity ty of the target and training networks network after training with 1500 samples.
3. Scalability Analysis The learning progress is analyzed for larger numbers of synaptic connections. The previous network structure was replicated but with 5, 10, 20, 40, 80 and 160 1 synaptic connections. Table 5 shows the similarity measure and membrane error for each of these networks network after training with 1500 samples. Table 5 Similarity measure and membrane error improvement for different number of synaptic connections co in a Type I DS DSNN Synaptic connections
Similarity (Mean)
Spike Count (Mean)
Membrane Error (Mean)
5
0.27→0.93
7.97
28.90→1.54 dB
10
0.37→0.91
13.13
25.62→6.60 dB
20
0.19→0.93
9.63
22.12→-3.70 dB
40
0.48→0.93
13.7
17.85→-5.3 dB
80
0.56→0.91
21.7
17.23→-0.81 dB
160
0.72→0.84
54.3
18.28→7.94 dB
Improvements to similarity measure and membrane error are observed independent of the number of synaptic connections. In spite of the fact that a larger network requires more training samples, the improvement in membrane prediction ction and similarity measure is significant for different network size sizes using the same number of training samples. Figure (8) shows the initial guess of synaptic parameters of
16
the Type I network with 80 synaptic connections connection and the trained parameter set after ter training with 1500 samples. The simulation results identify that the convergence to an optimum parameter set is repeatable and scalable. The simulation results also suggest that the neuron membrane potential is a unique function of incoming spike train and synapse parameters. Thus, the learning progress in predicting synapse parameters corresponds to the simultaneous improvement of similarity measure and membrane potential prediction. These results validate the learning capability for accurate predictions of synaptic parameters for both Type I and Type II DSNNs.
(a)
(b) Figure 8 Training progress in a Type I DSNN with 80 synaptic connections. connections a) Desired synapse parameters plus initial guess of parameters, b) Desired and trained synapse parameters trained by 1500 training samples.
The next question is the analysis of similarity range variation in the learning process. Similarity range is the core parameter in quantifying similarity between spike trains, and the following simulation simulations examine its impact on synaptic parameter estimation. 4. Similarity Range Analysis The similarity range defines the acceptable jitter in comparing the generated spiking train to the desired spike train.. Algorithm I utilizes the similarity range for matching the desired and tr training spikes for membrane prediction and eventually synaptic parameter adaptation. Larger values of the similarity range correspond to larger er errors in membrane prediction. Thus, it it is expected that the learning process fails to converge to the optimum parameter set as the similarity range increases increases.
17
30
0.95
25
0.9
0.85
0.8
1
P a ra m e te r E rr.
1
M e m b ra n e E rr (d B )
S im . M e a s u re
Spiking activity in a biological system is precise with an accuracy of about one millisecond. Thus, defining a much larger similarity range is not a viable tool for spike trains similarity assessment. In the following simulations, the similarity range for the learning algorithm changes from 2 to 20 milliseconds with one millisecond steps for a Type I DSNN with 10 synapses. The network is trained with only 500 input stimuli - Table 2 result indicates that 500 training samples is enough for an accurate estimation of the network parameters - , and its learning progress is analyzed on 100 test samples. Figure (9) shows the average similarity measure for the trained network with different similarity ranges. In addition to similarity measure, average membrane error and parameter estimation error are shown. The reference similarity range for performance analysis in the trained network is 2 milliseconds. Similarity measure shows a significant drop as the learning similarity range increases. This drop is predictable as the learning algorithm faces a larger error in membrane prediction. Membrane error grows along with similarity range growth; in fact, the membrane error grows to a level of an untrained network - refer to first row of Table 2 and 3. Similarly, the sum of square error in parameter estimation grows significantly. The simulation results imply that an accurate synaptic parameter adjustment requires a fine similarity range. It was already discussed that the maximum value for similarity range should be smaller than half of the refractory period. The similarity range only within a couple milliseconds - a range from 0 to 7 milliseconds - is a viable choice for the learning algorithm. For the similarity measure, it is possible to maintain the same level of similarity when the similarity range for both the training and the test process is kept equal. Despite maintaining a high similarity measure, the error in parameter estimation and membrane prediction diverges significantly - Figure (9.c). The result suggests that two networks with different synaptic parameters - here desired and training networks - may show a reasonable similarity measure while they have completely different temporal dynamics. Thus, similarity measure and the learning algorithm converge to optimum synaptic parameters for only a fine similarity range- a value smaller than half of refractory period.
20
15
0.8
0.6
0.4
10 0.75
0.2
5 0.7 2
4
6
8 10 12 14 Training Sim. Range
(a)
16
18
20
2
4
6
8 10 12 14 Training Sim. Range
(b)
16
18
20
0 2
4
6
8 10 12 14 Training Sim. Range
16
18
20
(c)
Figure 9 Learning progress as a function of the similarity range. a) Similarity measure versus similarity range. Similarity measure of the trained network is measured with a fixed value of 2 milliseconds, b) Membrane error versus similarity range, c) Error in parameter estimation for different similarity range. A lower similarity range shows a lower membrane and parameter estimation error, and it correlates with higher similarity measure. A longer similarity range is neither a viable choice to evaluate similarity between two spike trains nor an optimum parameter estimation rule for the learning process. Despite the fact that large similarity range maintains a reasonable level of similarity measure, it fails in accurate membrane potential and synaptic parameter prediction.
A complete set of simulations has been applied in examining off-line learning in Type I/II DSNNs. Different aspects of the learning algorithm including convergence, repeatability and scalability were examined. The simulation results verify an acceptable consistency in the system identification task independent of the network size and type. The simulation result along with the theoretical finding Appendix B - verify that the learning algorithm and its extension, gradient-gain learning, are a viable tool for simultaneous STP and LTP adaptation in DSNNs. Through the following simulations other
18
aspects of the DSNN model and the on-line gradient-gain learning will be analyzed. The first simulation examines the FD model’s capability in the prediction of CA1 pyramidal neuron spiking activity. B. FD Model and CA1 Neuron Spiking Activity Prediction A single pyramidal cell of the CA1 hippocampus area was excited by injecting APs on its presynaptic pathway with a 2-Hz Poisson spike timing. The patch-clamp technique was utilized for recording the somatic membrane potential of the CA1 neuron in response to the spiking train impinging on its Schaffer collateral synaptic pathway – the data was collected in the USC Center of Neural Engineering, Theodore W. Berger Lab -. Figure (10.a) shows 10 seconds of the CA1 pyramidal neuron membrane potential and its AP activity in response to the impinging spike train. Using an exhaustive search with a weighted least square criterion, free parameters of a single synapse model are adjusted to replicate the PSP dynamics and its spiking activity. In order to have a better estimate of the neuron's PSP response and spiking activity, the synapse model PSP response was replaced with the neuron's PSP extracted from the recording. Figure (10.b) presents the PSP prediction using the FD model plus its projected spiking activity. The average normalized error in PSP amplitude prediction is 27%. The model predicts CA1 neuron spiking activity with 6.6 milliseconds accuracy while generating only one extra spike. The parameter adaptation was applied for the whole recording time - 215 seconds - where the recorded neuron’s PSP dynamics showed other nonlinear dynamic processes including facilitation regulation and vesicle depression. For a short recording period, adjusting only ∆, parameters suffice to build a precise predictive model of the recorded neuron’s PSP response whereas a longer period was required to capture other dynamics such as facilitation regulation and vesicle depression. Appendix A describes a more intricate FD model for synapse dynamics. The on-line learning algorithm was applied in parameter prediction with a 7 milliseconds similarity range. Facilitation and recovery time constants were fixed by the outcome of the exhaustive search in which both ∆ and parameters converged to the values predicted by the exhaustive search algorithm. The membrane potential dynamics of the CA1 neuron was modeled by a single synapse with fixed time constants. However, multiple synaptic connections with different time constants - and probably synaptic delays - generate a more precise prediction of the neuron PSP dynamics. Table 6 shows the release prediction error using single and double synaptic connections. Increasing synaptic connections correlates with a better prediction of CA1 membrane potential dynamics and, correspondingly, its spiking activity prediction. It was already shown in the system identification section that the learning algorithm is a scalable process for multiple synaptic connections. Thus, the proposed learning algorithm will address synaptic parameter adjustment for multiple synaptic connections as well. Table 6 Performance of the FD model with multiple synaptic connections in predicting a CA1 neuron PSP Synapse parameters ∆T, N¼½¾ , SU , ST 0.115,10.132,500,200
0.115,10.132,500,200
0.115,10.132,500,200
Membrane error
Release error (%)
2.66
27.3%
Spike similarity (milliseconds) 6.6 + 1 extra spike
2.76
22.4%
6.3 + 1 extra spike
Through modeling the PSP response of the FD model was replaced by the realistic PSP shape of the CA1 neuron to provide a more accurate AP timing; this change was independent of the learning process or synapse temporal dynamics. Thus, it is replaced with its phenomenological approximate in the next simulations. Though it is shown that the multiple FD model with different time constants is a better approximate of the pyramidal neuron membrane potential dynamics, through the following simulations each synapse is modeled with only one FD model. However, it is possible to increase the number of synaptic connections to create a more intricate model where a set of synapses receives the same input spiking activity. The next section introduces a generative model to generate different spiking patterns of cortical neurons.
19
1
AP
0.8 0.6 0.4 0.2 0 0
1
2
3
4
5 Time(sec)
6
7
8
9
10
(a) 120
CA1 PSP/AP
100 80 60 40 20 0 -20
0
1
2
3
4
5
6
7
8
9
10
(b) 50 FD CA1 Threshold
CA1/FD PSP
40 30 20 10 0 -10
0
1
2
3
4
5
6
7
8
9
10
5
6
7
8
9
10
(c)
CA1/FD AP
2 CA1 FD
1.5 1 0.5 0 0
1
2
3
4
(d) Figure 10 CA1 PSP dynamics and synapse FD model response. a) Input spike train to a CA1 pyramidal neuron synaptic pathway, b) A single CA1 neuron PSP response and its AP activity, c) Predictive model of CA1 neuron PSP dynamics plus a synapse FD model PSP response, d) Spiking activity of the CA1 neuron plus the synapse FD model. The FD model with only one synaptic connection predicts a CA1 neuron spiking activity with an average of 6.6 milliseconds jitter accuracy plus firing an extra spike.
C. Cortical Neuron Spike Generative Model The Hodgkin-Huxley model is the most well-known biological model capable of generating the different spiking patterns of cortical neurons. Izhikevich proposed a mathematical model capable of generating different AP patterns. A Type II DSNN model can imitate the spiking patterns claimed in Izhikevich model. In the following simulations, synaptic parameters in a Type II DSNN network with only one feed-forward synaptic connection are adjusted to generate different spiking patterns of cortical neurons. The proposed DSNN model builds a spike-domain computation unit capable of mapping incoming spike trains to different spiking and bursting patterns. The network parameters - ∆ , , ∆ ¿ , ¿ - are adjusted to generate different spiking patterns in response to the input spike train. The learning process is on-line, and the model parameters are updated at the gain injection time. Input stimuli to the network have an ISI of 50 milliseconds and lasts for 300 milliseconds. The learning process is repeated 500 times, and in compliance with the on-line learning algorithm, the input stimuli are repeated with a 300 milliseconds silent gap - with this time gap, the network synapses retrieve from their previous release activities. The following gain function adjusts the synaptic parameters to generate the desired spiking pattern:
20
' b ∗ , b 0,1,2, … . 1 missing spike ∈ E' , ' , ' À1 extra spike ∈ E' , ' ,
' b ∗ , b 0,1,2, … .Å 0 Other time
(12 12)
definition is set to 3 milliseconds; thus, the network is trained to The similarity range - - in the '
generate the desired spike train within with 3 milliseconds accuracy. The training raining process for different spiking patterns start from the same initial value for the synaptic parameters. The he facilitation factor ∆ and ∆ ¿ - are set to 0.2 and 0.5, and synaptic strength - and ¿ - are set to 5.0 and 2.0. Training stops whenever the sum of absolute gain g is less than or equal to 1; the stopping ing criterion implies that the trained network is capable of generating at least all but one of the desired spikes within 3 millisecond accuracy. Figure (11) shows the learning progress and gain pattern for chattering - CH spiking activity. Chattering neurons urons can fire high-frequency bursts of 3 to 5 spikes with a relatively short inter-burst burst period in response to incoming stimuli. Figure (12) shows the neural model and the trained parameter set for different spiking activities of cortical neuron. In contrast to the system identification problem, the gain function is determined independent of membrane potential error. The gain function reflects the coincidence of the desired and training APs in a specific time interval; therefore, synaptic parameters are adjusted justed by an evaluative function of matched AP time rather than the network membrane error. Neural models with static synaptic connections fail to emulate ate the diverse spiking activities of a cortical neuron whereas dynamic synapses allow the neural model too generate different spiking and bursting patterns. The model introduced here builds ann adjustable spike-domain computational unit. The feed-forward forward and feed-back synapses of the spike generative model are able to introduce an adjustable delay in their response to the impinging APs.. By utilizing the synaptic delay, the model is capable of generating a diverse spiking pattern. The model odel response to incoming spike activity shows a diverse set of neural processing capabilities; the he processing consists of variable spike time delay, temporal selectivity of incoming activity, and neural synchrony. The next simulation utilizes a Type I network for a temporal pattern recognition task. 1 0.8 AP
0.6 0.4 0.2 0 0
50
100
150 Time (msec)
200
250
300
200
250
300
200
250
300
(c) 3
P S P /AP
2
(a)
1 0 -1 -2 -3
50
100
150 Time (msec)
(d) 1
G a in
0.5 0 -0.5 -1 0
50
100
150 Time (msec)
(e)
(b)
Figure 11 Parameter adaptation in the neuron model for the chattering - CH - spiking activity. a) Feed-forward forward synapse parameter adaptation, b) Feed-back back synapse parameter adaptation, c) Input spike train to cortical neuron, d) PSP and AP in neuron output, e) Gain pattern in response to incoming APs.
21
R S / G a in
Ref
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
50
100
150
200
250
300
T C 1 / G a in
50
100
150
200
250
300
T C 2 / G a in
50
100
150
200
250
300
R Z 1 / G a in
∆ ¿ , ¿
50
100
150
200
250
300
R Z 2 / G a in
2 1 0 -1
50
100
150
200
250
300
R Z 3 / G a in
∆ ,
50
100
150 Time (msec)
200
250
300
I B / G a in
2 1 0 -1
C H / G a in
(a)
F S / G a in
2 1 0 -1
L S / G a in
2 1 0 -1 2 1 0 -1
(b)
2 1 0 -1 2 1 0 -1 2 1 0 -1 2 1 0 -1
(c) 1
AP
0.8 0.6
2 1 0 -1
0.4 0.2 0
0
50
100
150 Time (msec)
200
250
300
(e)
(d)
Figure 12 Cortical neuron model plus different spiking patterns. patterns a) Cortical neuron model. The model is a Type II DSNN with one feed-forward synaptic connection, b) Parameters for the feed-forward forward synapse. There are 10 different parameter set, each one corresponds to one of the desired spiking activity. The RS spiking activity refers to a regular spiking, in which cortical neuron model generates an AP in response to each impinging AP. The The IB spiking activity refers to intrinsically bursting spiking activity; the cortical neuron fires multiple APs for the first incoming AP and it responds with regular spikes afterward. The CH spiking activity corresponds to a chattering neuron; the cortical cortical neuron responds to each impinging AP with multiple APs. Neuron generates more APs in responds to the first AP, and number of APs drops in the next spikes. The FS spiking activity corresponds to fast response neurons, in which neuron fires multiple spikes spike per impinging APs. The LTS refers to a low threshold spiking activity, which is similar to FS expect the neuron fires lower number of APs per incoming AP. The TC1 and TC2 spiking activity refers to Thalamo-Cortical Cortical neuron response which has a temporal selectivity. lectivity. The TC1 neuron responds to incoming AP with a delay, and the TC2 is responding to the first couple of APs. The RZs neuron define other possible temporal selectivity of cortical neurons,, each neuron is responsive to specific speci APs of incoming spike train, c) Parameters arameters for feed feed-back synapse, d) Input stimuli to cortical neuron model, e) Cortical ortical neuron model response for different spiking activity plus gain pattern. Cortical neuron training stops whenever there is an absolute gain sum less than or equal to 1.
22
B. Temporal Pattern Classifier A DSNN based neural classifier is built for recognizing different classes of spiking activity - 3 -. A Type I DSNN consisting of neurons with H synaptic connections per neuron - H 60 in this problem - is built in which synapses project incoming spiking activity to their corresponding neurons. Each neuron is responsible for encoding a specific class and it is trained to fire a higher number of APs for its corresponding class. In fact, the synaptic connections of the network bind APs in multiple spike trains to build the desired post-synaptic firing activity. The synaptic parameters of each neuron are adjusted to be sensitive to temporal features of the corresponding activity class while being insensitive to the other activity classes. The gain function of the learning algorithm is responsible for observing neurons' spiking activity to adjust network synaptic parameters by injecting the proper gain. Equation (13) defines the gain function for each output neuron at different processing time: 1 ∑ÉÊ APÆ n t k ∈ Winning Class Å Æ ' Ç ∑ÉÊ APÆ n t
k ∈ Other Classes
(13)
Parameter defines the ISI for the spiking activity in the winning class. The learning process utilizes competitive learning among the neurons. Each neuron tries to fire more APs in response to its matched incoming activity while at the same time forcing the other neurons to be silent. It is also possible to define new gain functions built upon other temporal characteristics of incoming spiking activity. The gain function can train individual neurons to fire APs sooner in response to their matched incoming spiking activity while postponing firing activity in other neurons to a later time. This second gain function is defined by: Æ ' Í
min 1,2 ∑ÉÊ APÆ n t k ≡ Winning Class, n ¬ TÐ ∑ÉÊ APÆ n t
k ∈ Other Classes, n ¬ TÐ Å 0 n TÐ
TÐ Time of Ñirst AP in winning class D,
D deÑines time gap margin
(14)
Using two different gain functions, the synaptic parameters of the network are adjusted to capture different features of the incoming spiking trains. The first gain function tends to adjust synaptic parameters of the network for a prolonged spiking activity in the winning neurons while the second gain function adjusts synaptic parameters for an abrupt AP firing. To maintain a prolonged spiking activity, the synaptic parameters of a neuron generally move toward potentiating synapses with strong synaptic connections. In contrast, synaptic parameters of a fast firing neuron move toward a depressing synapse making the synapses sensitive to isolated spiking activity. The optimum synaptic parameters are dependent on the temporal characteristics of the incoming spiking trains as well. The discriminative information of a class are encoded in the spike time of individual spike trains and the mutual spiking activity of different spiking trains. The rule of the gain function is to correlate incoming spike train features with the synaptic parameters of the neurons to generate a desired spike pattern in the network output. Thus, optimum synaptic parameters of a network become a function of the incoming spike train, the synaptic parameters, and the gain function. Network training stops whenever the gradient gain converges to zero. The gain function converges to zero when the neurons generate the desired spiking activity. Figure (13) shows the network model plus a sample spiking train for each class. The time period of the input stimuli is 400 milliseconds. There are 500 different samples of spiking activities per each activity class. The network is trained with the first 400 samples, and its performance is examined using the remaining 100 samples. The free parameters include 180 pairs of ∆, while the parameters of the inhibitory synaptic connections are fixed at 1.0,1 . The function of the inhibitory synapse is to reduce the activity of all the neurons except the winning neuron whenever there is an uncertainty about the incoming spiking activity. Whenever there is a distinct class of incoming activity, the winning neuron fires with a higher firing rate in comparison to the other neurons without the help of inhibitory
23
synaptic connection. Facilitation factor in the inhibitory synaptic connections are manually set to 1.0 with synaptic strength in order of 1.0. Therefore, inhibitory synapse has depressing dynamics; synapses with depressing dynamics are sensitive to sparse neural activity which matches the dynamics of the network spiking activity in response to indistinct input stimuli. Ref
30
25
Synapse
20 Ref
15
10 Ref
5
0 0
150 200 250 Time (msec)
300
350
400
300
350
400
30
25
25
20 Synapse
20 Synapse
100
(c)
(a) 30
15
15
10
10
5
5
0 0
50
0 50
100
150 200 250 Time (msec)
300
350
400
0
50
100
150 200 250 Time (msec)
(d)
(b)
Figure 13 Temporal pattern classifier and sample spiking activities. a) Classifier architecture. Synaptic connections to other neurons have the same structure of the first neuron, b) Sample spiking activity of class A - there are 30 input spike trains, each spike train feeds two synaptic connections, c) Sample spiking activity of class B, d) Sample spiking activity of class C.
For training purposes, input stimuli are concatenated with a 400 millisecond silent time gap between different training samples. The initial values for the synaptic parameters are fixed for all three classes, and learning rate is set to 0.001. The learning process repeats 10 times on the training data before the classification performance is examined on the test samples. Table 7 shows the network performance using both higher spike count and time of the first AP gain functions. Network performance using spike count on trained data reaches 95.6% at the end of the training process, and its performance on the test sample is 96.5%. Network performance using spike time reaches 76.6% with an average delay time of 38 milliseconds - less than 16 samples - on correctly classified samples, and the network performance in the test samples is equal to 80.8% with an average delay time of 36 milliseconds - less than 15 samples. The delay time determines the time period until the first AP fires in the wining class. The delay time defines the time window that the network requires to reach a decision about the class of activity. The performance with the spike time classifier implies that the trained network only requires 38 milliseconds - less than 16 samples as input data has been sampled with 400 samples per second of data to decide the activity class with an 80% accuracy. The spiking activity provided in this section belongs to a chain link fence intrusion activity database the fence activity database was collected by USC Neural Dynamics Laboratory in 6 different fence panels which include 300 recordings of 25 seconds - with three different classes of sequential data. The database includes 3-axis vibration signals recorded with 12 bit samples with a sampling rate of 400 Hz. In generating spike trains, the input signal is filtered by a filter bank - consisting of 5 different frequency bands -, and it is then passed to the absolute value function. Each pre-processed output is magnified and then passed to an integrate-fire neuron with a fixed firing threshold. The neuron fires an AP whenever its membrane potential passes the firing threshold. There are two gain levels for each
24
pre-processed output; thus, there are in total 10 neurons mapping input stimuli to spike trains for each axis. In all, there are 30 different spike patterns for each 400 milliseconds of a fence activity. Each spike train feeds two synaptic connections where the initial values of the synaptic parameters in these two synapses are different. To compare the network performance with well-established artificial neural networks, a 2-layer neural network with 15 input, 50 hidden, and 3 output neurons are trained to classify the same activity classes. The feature set of the artificial neural network are defined by energy of the 5 filter banks for each vibration axis. The time period of each activity is partitioned in two parts, and the energy of each partition is added to the feature set. Table 7 shows the performance of the neural network and a Gaussian mixture model - GMM - in classifying the same activity classes. The DSNN performance is comparable with GMM and 2-layer neural network classifiers; however, the number of free parameters for the DSNN are much lower than the other classifiers. The DSNN performance is almost independent of feature selection whereas the performance of other classifiers is dependent on the choice of features. Table 7 Performance of different classifiers on the test dataset of the fence intrusion database Classifier
Performance (%)
Parameters
GMM
97.4
GMM
90.7
Neural Network
97.8
Neural Network
93.9
DSNN (AP Count)
96.5
Gaussian mixture model with 2 components Energy in the first and second half of input are considered as the feature set 30 features per input sample 1862 free parameters (900+900+30+30+1+1) Gaussian mixture model with 2 components 15 features per data 482 free parameters (225+225+15+15+1+1) 2 hidden layers, each layer has 40 neurons Energy in the first and second half of input are considered as the feature set 30 features per data More than 4000 free parameters 2 hidden layers, each layer has 40 neurons 15 features per data More than 2000 free parameters 360 free parameters (60×2×3)
DSNN (First AP)
80.8
360 free parameters (60×2×3)
The DSNN network is also capable of deciding about the incoming activity class in a shorter time period than the other classifiers. Feature extraction in the other classifiers requires a pre-defined time window; otherwise, variability in the feature set significantly drops the classifiers' performance. Neural network and GMM classifiers are incapable of performing on original data; in fact, their learning algorithms are incapable of maintaining a reasonable performance on large feature sets with extensive numbers of free parameters. The DSNN classifier performance along with its superior features: 1) lower number of free parameters, 2) optimum decision in minimal time, and 3) a feature-free classifier make the model the superior choice for temporal pattern classification. V.
DISCUSSION
This research introduces an artificial model incorporating two progressive transitions to develop a more realistic neural model: 1) synaptic response is not a static process; it is a dynamic process changing in the time-course of impinging spiking activity, and 2) synaptic adaptation is not merely a synaptic weight regulation; it is a dynamic process regulating both long-term and short-term characteristics of each synapses. The novel neural model introduces a family of neuro-computational processing engines applicable in modeling a cortical neural circuitry's spike-domain functionality and a variety of spatio-temporal processing tasks. A linear state space model of synaptic transmission was introduced providing an accurate and compact mathematical representation of the synaptic dynamics. A comprehensive analytical study of the
25
synaptic response including its efficacy, temporal dynamics, and statistical behavior is provided as well. The dynamics of synaptic response with respect to variation in its internal parameters is formulated; the synapse PSP gradient with respect to its parameters is the key tool utilized in DSNN training. The linear state space representation is extended to the neural level where a DSNN includes multiple neurons with thousands of synaptic connections. A supervised spike-in-spike-out learning rule for the synaptic adaptation of DSNNs is proposed. The learning rule is biologically plausible, and it is capable of simultaneously adjusting both LTP and STP factors of individual DSNN synapses. The learning algorithm is extended to an on-line learning algorithm called gradient-gain learning. Gradient-gain learning is a spike-domain learning rule that adjusts synaptic parameters of individual DSNN synapses to perform a spike-domain functional mapping. A comprehensive simulation analysis is performed to verify accuracy, repeatability, and scalability of the learning algorithm. The system identification test confirms the off-line learning algorithm ability to predict synaptic parameters of a DSNN containing up to 160 synaptic connections. The system identification task also confirms consistency and convergence in the synaptic adaptation process for both Type I/II DSNNs. A new similarity measure is introduced to assess similarity between two spike trains. The similarity measure identifies each spike of the test train - in comparison with the desired spike train - to three classes of: 1) similar, 2) missing, or 3) extra where the labeling process is insensitive to spiking jitter. Similarity range is the free factor in the spike similarity measure defining the maximum acceptable jitter. The similarity range plays a crucial role in the synaptic adaptation mechanism; this is because the learning objective function minimized by the supervised learning rule is an approximate measure of the similarity assessment process. The relationship between spike train similarity and neural membrane dynamics in the synaptic adaptation process is analyzed to estimate the optimum similarity range, where the simulation result estimates an optimum similarity range of about a couple of milliseconds. The synapse model’s ability in predicting PSP dynamics - and spiking activity - of CA1 pyramidal cell is demonstrated. A single synaptic connection with only 2 free parameters provides the same accuracy of second-order Volterra series approximation with 14 free parameters (Ghaderi, 2012). Gradient-gain learning with a /1 gain for missing and extra spikes has been applied for synaptic parameter adjustment where the parameters converge to the same parameters optimized using an exhaustive search. The DSNN is intrinsically a time-domain computational engine; it is applicable for processing dynamic neural/artificial signals. Gradient-gain proposes a flexible synaptic adaptation mechanism to shape DSNN temporal processing. A Type II DSNN model with only a single feed-forward synaptic connection is introduced to generate different spiking/bursting patterns of the cortical neurons. Gradient-gain with a /1 gain at missing and extra spikes shapes the model dynamics for the desired mapping. This recurrent neural structure introduces an adjustable generative model capable of performing a wide set of spike-domain mappings. A DSNN model with multiple output neurons has been applied in performing a temporal pattern classification task. The DSNN classifier performance competes with a 2-layer neural network and surpasses the GMM classifier - 96.5% classification rate -, while using a significantly fewer number of free parameters. The gain function tunes synaptic parameters of the winning neuron by injecting 1 gain resulting in a regular spike train - with time period - in response to the matching input while the other neurons are punished by 1 gain for their spiking activity in response to un-matched incoming input. It is also shown that DSNN can be trained to form an optimum decision in a minimal time, a property which is unachievable by the other competitor models. The DSNN is capable of performing the classification task with an 80% accuracy after observation only 36 milliseconds of the data - each data lasts for 400 milliseconds. Again, the gain function rewards the winning neuron with 1 gain and punishes the other neurons with 1 gain to bring winning neuron’s spiking activity to an earlier time while pushing other neurons’ spiking activity to a later time.
26
Synapse temporal dynamics is the natural substrate for neural computation. Neither of the spikedomain mappings introduced in this research are feasible for the networks with static synapses - or exclusively synaptic weight adaptation. In fact, simultaneous regulation of STP and LTP factors are required to shape the neural computations performed in this research. Temporal processing in a static synapse model is limited to the overlap time between PSPs of consecutive APs - practically less than 100 milliseconds -, whereas DSNN temporal processing is expandable to longer time periods using LTP and STP factor regulation. DSNN performs a diverse set of spike-domain processing and can be performed on a sparse spiking activity or a sequence of spike train expanded to multiple seconds. The simulation results for CA1 spike activity prediction and the temporal pattern classifier show how DSNN model performance surpasses state-of-the-art mathematical tools for neural modeling and dynamic classification tasks. The DSNN gains its computational power by: 1) augmenting synaptic dynamics, and 2) incorporating simultaneous LTP and STP adaptation. The DSNN model and the spike-generative model can be used for both modeling neural circuitries and their spiking activity prediction. VI.
CONCLUSION
This research is an attempt to show how a transition from a static synapse model to a dynamic synapse model that also includes STP and LTP adaptation is achievable. It is also shown how the computational power of the artificial neural model increases by augmenting accurate models of the corresponding neural circuitry. There are many other unknown questions about neural dynamics and the synaptic adaptation mechanism. The learning addressed here is limited to a single layer DSNN but can be extended to multiple layers or a cascading neural topology structure as observed in the brain. The spike-domain mappings performed here address a small sub-set of possible spike domain mapping observed in the brain. The next question is how the learning and neural model can be modified to perform a more diverse neural processing. The authors’ future research focuses on applying DSNN modeling in predicting spiking activity of the hippocampus cortical region where it requires a cascading neural structure.
APPENDIX A DETAILED FD MODEL This section provides the intricate model of synapse temporal dynamics. The model addresses: 1) multiple facilitation mechanisms, 2) the Ca2+-dependent vesicle recovery time constant, 3) activitydependent facilitation regulation, and 4) activity-dependent vesicle depletion. The FD model proposed in equation (1) is a simplified derivation of the following intricate model, in which by an appropriate setting of the model parameters it turns to the simplified FD model already introduced in equation (1). The calcium facilitation factor controlling calcium concentration in synapse on the pre-synaptic side changes through the time course of incoming spiking activity. In addition, multiple facilitation mechanisms with different time constants have been addressed for the calcium concentration kinetics. Equation (A.1) to (A.4) provide two facilitation mechanisms in which their Ca2+ facilitation factor is a function of impinging spiking activity to the synapse. ⁄ \ 1 a ∗ ∆a ∗ a ⁄ a 1⁄a
⁄
∗ a ∗ ⁄ 1 ∗ ∆ ∗ a
¿ ⁄ ¿ \ 1 a¿ ∗ ∆a¿ ∗ a¿ ⁄ a¿ 1⁄a
(A.1) (A.2) (A.3)
¿ ∗ a¿ ∗ ¿ ⁄ ¿ ¿ Ô ¿ 1 ¿ ∗ ∆ ¿ ∗ a
(A.4)
? ⁄ ? 1⁄? ⁄ Õ 1 ? ∗ ∆? ∗
(A.5)
Calcium-dependent vesicle recovery regulates the recovery time constant of the released vesicles in response to prolonged spiking activities. The vesicle recovery time constant decreases in response to incoming spiking activities guaranteeing a synaptic response during long periods of impinging spikes. Equation (A.5) defines the Ca2+-dependent mechanism in vesicle recovery time constant regulation.
27
The vesicle release dynamics is defined by:
⁄ 1 ∗ ? ∗ ? ∗ 1⁄ ¿ ∗ ¿ ∗ ∗
(A.6)
The other crucial mechanism in post-synaptic kinetics is the vesicle depletion process. The ability of post-synaptic side to absorb the released vesicles drops after a prolonged spiking activity. The following equation models the permeability of the post-synaptic membrane in response to pre-synaptic spiking activities. ⁄ 1 ⁄ É ±! ∗ ∆ ∗
(A.7)
The synaptic PSP response is defined by:
⁄ ⁄ ∗ ∗ ¿ ∗ ¿ ∗ ∗
(A.8)
(A.10)
⁄ ⁄ ∗ ∗ ¿ ∗ ¿ ∗ ∗
(A.9)
Equation set (A) defines the intricate model of synaptic response. The PSP response is modeled by an approximate α-synapse function, but it is possible to be defined with other PSP shapes.
APPENDIX B LINEAR OBJECTIVE FUNCTION, SIMILARITY MEASURE AND UNIQUNESS OF DSNN MEMBRANE POTENTIAL In this section, several theories are developed to support the convergence of the proposed learning algorithm. The first theory addresses the relationship between the learning objective function and the similarity measure. The second theory discusses the relationship between the neuron membrane potential and synaptic parameters in Type I/II DSNNs. A few of corollaries extends the relationship for spiking activity in Type I/II DSNN and its synaptic parameters. Theorem I: A decrement in the expected value of the learning objective function - defined in equation (10.a) - corresponds to an increase in the number of similar spikes. Proof: In a DSNN with excitatory synaptic connections, the PSP derivative in both the test and target DSNNs have the same sign - independent of ∆, parameters value - at different processing bins. With the mentioned assumption, the expected value of error at missing and extra spike times is defined by: ccÖ × c∑!∈ØÙn·p ! k ∑!∈ÚÛmmÛÜ k ! nÚÛmmÛÜ nØÙn·p ∗ c
It is assumed that ! gets a uniform distribution between E0 k . Parameter c is the expected value of membrane potential error ∎ Theorem II: The membrane potential of a synapse is a unique function of its synaptic parameters. The free parameters of each synapse are ∆, . Proof: The PSP response in the DSNN is described by equation (4.a). The BC and DC matrices are defined independent of the free synapse parameters, and they are the same for different synapses. It can be shown that B and D elements are unique functions of synapse free parameters. Thus, each parameter set generates a unique membrane potential in response to different pattern of spiking activities. (B.1)
The theorem is valid by the ApÞ AßÞ and BÞp BÞß assumption, and it is not necessarily valid for DSNN models with other free parameters including recovery and facilitation time constants ∎ Corollary I: The membrane potential of a DSNN network is a unique function of its synaptic parameters. Again, the free parameter of each synapse are ∆, . Proof: The DSNN network dynamics follow the same state space model defined for a single synapse. Thus, it is possible to consider DSNN with multiple synaptic connections as a single synapse; therefore, Theorem II is valid for the DSNN network ∎
28
The learning algorithm deals with spiking activity, but the more important question is the relationship between DSNN spiking activity and its synapse parameters. In contrast to membrane potential uniqueness, spiking activity is not a unique function of synaptic parameters. For example, DSNN networks failing to fire APs can share similar spiking activities while their synaptic parameters can be different. Thus, the spiking activity and synaptic parameters relationship can be discussed in only DSNN networks firing APs. It is also assumed that the processing resolution is fine enough that neurons fire APs at k à in which à converges to zero. DSNN networks with different numbers of synaptic connections fail to generate similar spiking activities; it can be assumed that the extra synaptic connections perturb the neuron membrane potential about à - around a spike time - changing its spiking time. Corollary II: DSNN networks with similar spiking activities in response to all possible input spike trains have similar synaptic parameters. Proof: Assume that all except one - kth synapse - synaptic connection of two DSNNs share the same parameters. Moreover, assume that there is a specific pattern of spike requiring the involvement of the kth synapse to fire an AP. Corollary I showed that the PSP of the remaining synaptic connections is the same. Thus to have a similar spike time, the kth synapse in both networks is required to have the same synaptic parameters. As a result, precise spike timing in two DSNN networks in response to all possible input spike trains requires the same parameter set and network topology ∎ The simulation result in the system identification task shows that the distance between free parameters of two DSNNs shrinks as their similarity measure increases. The convergence of parameters in the test and training network demonstrates the result of corollary II. In fact, the simulation is performed with a limited number of input spike trains and a coarse processing resolution. The coarse processing resolution relaxes the à constraint which introduces deviation in the synaptic parameter estimation. The corollary's result implies that there is a possibility of multiple parameter sets that can generate similar spiking activity whenever the number of input spike patterns is limited. The simulation result is in accordance with this assumption. As the number of input spike trains increases, the expected value of membrane potential decreases which always correlates with a higher value of similarity measure.
APPENDIX C GAIN-STDP The synaptic PSP in response a spike train can be written as: u' ∑!âã ' ' ∗ B' ∗ á' '
(C.1)
ä ä ä u' × ' ' ∗ á' ' ∗ B'
(C.2)
' is PSP response modeled by α-synapse function, B' is the vesicle release at ' spike time, and á' is the step function. ' is a decaying function at its tail. Thus the PSP response can be ä -:. approximated by its vesicle release at the last impinging spike time - ' The PSP derivative with respect to synaptic parameters is defined by: ä ⁄ ä ä u' ⁄ × ' ' ∗ á' ' ∗ B'
u' ⁄∆ × ' ä B'
ä '
∗
ä ⁄ B' ∆
ä '
∗ á'
(C.3.a)
ä '
(C.3.b)
- is a function of synaptic parameters and the history of incoming - vesicle release at ä is a linear function of , but it is a non-linear function of spikes - equation set (4) and (5). B' ä can be written as: facilitation factor and timing of previous spiking activity. B' ä ät# ät1 ä (C.4) B' j ∗ F∆, ' , ' , ' , ⋯ , ' The gradient-gain at gain injection time ' is defined by: ä ä ät# ät1 ä , ' , ' , ⋯ , ' ← ' ∗ ' ' ∗ á' ' ∗ F∆, ' ä ∆ ← ∆ ' ∗ ' ' ∗ j
∗ ∂F
ä ät# ät1 ⁄ , ' , ' , ⋯ , ' ∆ ∆, '
29
ä ∗ á' '
(C.5.a) (C.5.b)
Equation set (C.5) is a function of pre-synaptic activity and gain injection time regulated by the gain value. Equation set (C.5) defines the gain-STDP learning rule. Function defines the normalized vesicle release. The F function and its derivative are described in equation set (4) and (5).
REFERENCES
1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13.
14.
15. 16. 17. 18. 19. 20. 21. 22. 23. 24. 25. 26.
Abbott, L. F., and Wade G. Regehr. "Synaptic computation." Nature 431.7010 (2004): 796-803. Morrison, Abigail, Markus Diesmann, and Wulfram Gerstner. "Phenomenological models of synaptic plasticity based on spike timing."Biological cybernetics 98.6 (2008): 459-478. Natschlger, Thomas, Wolfgang Maass, and Anthony Zador. "Efficient temporal processing with biologically realistic dynamic synapses." Network: Computation in Neural Systems 12.1 (2001): 75-87. Liaw, Jim-Shih, and Theodore W. Berger. "Dynamic synapse: Harnessing the computing power of synaptic dynamics." Neurocomputing 26 (1999): 199-206. Tsodyks, Misha, Klaus Pawelzik, and Henry Markram. "Neural networks with dynamic synapses." Neural computation 10.4 (1998): 821-835. Dittman, Jeremy S., Anatol C. Kreitzer, and Wade G. Regehr. "Interplay between facilitation, depression, and residual calcium at three presynaptic terminals." The Journal of Neuroscience 20.4 (2000): 1374-1385. Mongillo, Gianluigi, Omri Barak, and Misha Tsodyks. "Synaptic theory of working memory." Science 319.5869 (2008): 1543-1546. Liaw, Jim‐Shih, and Theodore W. Berger. "Dynamic synapse: A new concept of neural representation and computation." Hippocampus 6.6 (1996): 591-600. Kreuz, Thomas, et al. "Measuring spike train synchrony." Journal of neuroscience methods 165.1 (2007): 151-161. Gütig, Robert, and Haim Sompolinsky. "The tempotron: a neuron that learns spike timing–based decisions." Nature neuroscience 9.3 (2006): 420-428. Ponulak, Filip, and Andrzej Kasiński. "Supervised learning in spiking neural networks with ReSuMe: sequence learning, classification, and spike shifting."Neural Computation 22.2 (2010): 467-510. Bohte, Sander Marcel. "Spiking neural networks." Unpublished doctoral dissertation, Centre for Mathematics and Computer Science, Amsterdam(2003). Ghosh-Dastidar, Samanwoy, and Hojjat Adeli. "A new supervised learning algorithm for multiple spiking neural networks with application in epilepsy and seizure detection." Neural Networks 22.10 (2009): 14191431. Yousefi, Ali, Alireza A. Dibazar, and Theodore W. Berger. "Supervised learning in a single layer Dynamic Synapses Neural Network." Neural Networks (IJCNN), The 2011 International Joint Conference on. IEEE, 2011. Izhikevich, Eugene M. "Simple model of spiking neurons." Neural Networks, IEEE Transactions on 14.6 (2003): 1569-1572. Gerstner, Wulfram, and Werner M. Kistler. Spiking neuron models: Single neurons, populations, plasticity. Cambridge university press, 2006. Mainen, Zachary F., and Terrence J. Sejnowski. "Reliability of spike timing in neocortical neurons." Science 268.5216 (1995): 1503-1506. Tiesinga, Paul HE, J-M. Fellous, and Terrence J. Sejnowski. "Attractor reliability reveals deterministic structure in neuronal spike trains." Neural computation 14.7 (2002): 1629-1650. Berry, Michael J., David K. Warland, and Markus Meister. "The structure and precision of retinal spike trains." Proceedings of the National Academy of Sciences 94.10 (1997): 5411-5416. Hatsopoulos, N., et al. "At what time scale does the nervous system operate?."Neurocomputing 52 (2003): 25-29. Zucker, Robert S., and Wade G. Regehr. "Short-term synaptic plasticity."Annual review of physiology 64.1 (2002): 355-405. Abbott, L. F., et al. "Synaptic depression and cortical gain control." Science275.5297 (1997): 221-224. Maass, Wolfgang, and Anthony M. Zador. "Dynamic stochastic synapses as computational units." Neural Computation 11.4 (1999): 903-917. Varela, Juan A., et al. "A quantitative description of short-term plasticity at excitatory synapses in layer 2/3 of rat primary visual cortex." The Journal of neuroscience 17.20 (1997): 7926-7940. Bohte, Sander M., Joost N. Kok, and Han La Poutre. "Error-backpropagation in temporally encoded networks of spiking neurons." Neurocomputing 48.1 (2002): 17-37. Booij, Olaf. "Temporal pattern classification using spiking neural networks."Unpublished master's thesis, University of Amsterdam (2004).
30
27. Davis, Bryan A., et al. "Supervised synaptic weight adaptation for a spiking neuron." Neural Networks, 2003. Proceedings of the International Joint Conference on. Vol. 4. IEEE, 2003. 28. Gütig, Robert, and Haim Sompolinsky. "The tempotron: a neuron that learns spike timing–based decisions." Nature neuroscience 9.3 (2006): 420-428. 29. Fang, Huijuan, Yongji Wang, and Jiping He. "Spiking neural networks for cortical neuronal spike train decoding." Neural computation 22.4 (2010): 1060-1085. 30. Storck, Jan, Frank Jäkel, and Gustavo Deco. "Temporal clustering with spiking neurons and dynamic synapses: towards technological applications." Neural networks 14.3 (2001): 275-285. 31. Rostro, Horacio, et al. "Back-engineering of spiking neural networks parameters." arXiv preprint arXiv:0905.4810 (2009). 32. Song, Sen, Kenneth D. Miller, and Larry F. Abbott. "Competitive Hebbian learning through spike-timingdependent synaptic plasticity." Nature neuroscience 3.9 (2000): 919-926. 33. Senn, Walter, Henry Markram, and Misha Tsodyks. "An algorithm for modifying neurotransmitter release probability based on pre-and postsynaptic spike timing." Neural Computation 13.1 (2001): 35-67. 34. Strain, Thomas J., et al. "A supervised STDP based training algorithm with dynamic threshold neurons." Neural Networks, 2006. IJCNN'06. International Joint Conference on. IEEE, 2006. 35. Bohte, S. M., and M. C. Mozer. "Reducing spike train variability: A computational theory of spike-timing dependent plasticity." Proceedings of Advances in Neural information processing (NIPS’05): 270-279. 36. Masquelier, Timothée, Rudy Guyonneau, and Simon J. Thorpe. "Spike timing dependent plasticity finds the start of repeating patterns in continuous spike trains." PLoS one 3.1 (2008): e1377. 37. Izhikevich, Eugene M., and Niraj S. Desai. "Relating stdp to bcm." Neural computation 15.7 (2003): 15111523. 38. Bush, Daniel, et al. "Reconciling the STDP and BCM models of synaptic plasticity in a spiking recurrent neural network." Neural computation 22.8 (2010): 2059-2085. 39. Burkitt, Anthony N., Hamish Meffin, and David B. Grayden. "Spike-timing-dependent plasticity: the relationship to rate-based learning for models with weight dynamics determined by a stable fixed point." Neural Computation 16.5 (2004): 885-940. 40. Pfister, Jean-Pascal, et al. "Optimal spike-timing-dependent plasticity for precise action potential firing in supervised learning." Neural computation 18.6 (2006): 1318-1348. 41. Fiete, Ila Rani. "Learning and coding in biological neural networks". Diss. Harvard University Cambridge, Massachusetts, 2003. 42. Xie, Xiaohui, and H. Sebastian Seung. "Learning in neural networks by reinforcement of irregular spiking." Physical Review E 69.4 (2004): 041909. 43. Williams, Ronald J. "Simple statistical gradient-following algorithms for connectionist reinforcement learning." Machine learning 8.3-4 (1992): 229-256. 44. Fiete, Ila R., and H. Sebastian Seung. "Gradient learning in spiking neural networks by dynamic perturbation of conductances." Physical review letters97.4 (2006): 048104. 45. Legenstein, Robert, Dejan Pecevski, and Wolfgang Maass. "A learning theory for reward-modulated spiketiming-dependent plasticity with application to biofeedback." PLoS Computational Biology 4.10 (2008): e1000180. 46. Lee, Kyoobin, and Dong-Soo Kwon. "Synaptic plasticity model of a spiking neural network for reinforcement learning." Neurocomputing 71.13 (2008): 3037-3043. 47. Hinton, Geoffrey E., Simon Osindero, and Yee-Whye Teh. "A fast learning algorithm for deep belief nets." Neural computation 18.7 (2006): 1527-1554. 48. Belatreche, Ammar, Liam P. Maguire, and Martin McGinnity. "Advances in design and application of spiking neural networks." Soft Computing 11.3 (2007): 239-248. 49. Jin, Yaochu, Ruojing Wen, and Bernhard Sendhoff. "Evolutionary multi-objective optimization of spiking neural networks." Artificial Neural Networks–ICANN 2007. Springer Berlin Heidelberg, 2007. 370-379. 50. Deco, Gustavo, and Bernd Schürmann. "Spatiotemporal coding in the cortex: information flow-based learning in spiking neural networks." Neural computation11.4 (1999): 919-934. 51. Gupta, Ankur, and Lyle N. Long. "Hebbian learning with winner take all for spiking neural networks." Neural Networks, 2009. IJCNN 2009. International Joint Conference on. IEEE, 2009. 52. Bartlett, Peter L., and Jonathan Baxter. "Infinite-horizon policy-gradient estimation." arXiv preprint arXiv:1106.0665 (2011). 53. Ghaderi, Viviane S., et al. "Analog low-power hardware implementation of a Laguerre-Volterra model of intracellular subthreshold neuronal activity."Engineering in Medicine and Biology Society (EMBC), 2012 Annual International Conference of the IEEE. IEEE, 2012.
31
