Journal of Robotics, Networking and Artificial Life, Vol. 4, No. 1 (June 2017) 97–101___________________________________________________________________________________________________________A parameter optimization method for Digital Spiking Silicon Neuron model
Takuya Nanami

Institute of Industrial Science, University of Tokyo
Tokyo, Japan
Filippo Grassia

LTI Lab., University of Picardie Jules VerneSaint-Quentin, FranceTakashi Kohno

Institute of Industrial Science, University of Tokyo
Tokyo, Japannanami@sat.t.u-tokyo.ac.jp, kohno@sat.t.u-tokyo.ac.jp

AbstractDSSN model is a qualitative neuronal model designed for efficient implementation in a digital arithmetic circuit. In
our previous studies, we extended this model to support a wide variety of neuronal classes. Parameters of the DSSN
model were hand-fitted to reproduce neuronal activity precisely. In this work, we studied automatic parameter
fitting procedure for the DSSN model. We optimized parameters of the model by the differential evolution
algorithm in order to reproduce waveforms of the ionic-conductance models and reduce necessary circuit resources
for the implementation.
Keywords: Silicon neuronal network, Spiking neuron model, Differential evolution, FPGA

bursting, low-threshold spiking, elliptic bursting, andparabolic bursting.
To reproduce a variety of spiking properties in neurons,appropriate parameter sets for neuronal models have to
be found. The relationship between the values of theparameters and the behavior of the model is generally
complex, and it is difficult and time-consuming to find
an appropriate parameter set. Automatic parameter
fitting procedures have been studied as a solution to thisproblem. Parameters for Hodgkin-Huxley-type models
have been tuned automatically in [6] and [7]. Pospischil
et al [8] reproduced regular and fast spiking neuron
classes by simulated annealing based optimizationmethod using a cost function composed of the firing rate
and the adaptation time. Buhry et al. [9][10] compared
several heuristic algorithms for automatic parameter
estimation on a Hodgkin-Huxley-type model, where thedifferential evolution algorithm had a best performance.
They also applied the algorithm to a neuromimetic
analog integrated circuit.
In this work, we applied the differential evolution
algorithm to parameter fitting of the DSSN model inregular and fast spiking settings, because it shares its

1. Introduction
A wide variety of neuronal models have been used insilicon neuronal networks because of the trade-off
between reproducibility of neuronal activities andcomputational efficiency. For example, ionicconductance models can reproduce a neuronal activity
accurately but demands enormous computational
resources in large-scale implementations. In contrast,
integrate-and-fire based models such as the LIF andIzhikevich models [1] can be implemented with lessresource, because they approximate a spiking process by
the resetting of the state variables. However, they have
reduced reproducibility of complex neuronal activities.
For example, these models assume fixed maximum
membrane potentials during the spike process, whereas
the spike intensity is nonuniform in the nervous system[2]. The DSSN model [3] is a qualitative neuronalmodel designed for efficient implementation in a digitalarithmetic circuit. In our previous studies [4][5], weextended the DSSN models to support various neuronal
activities; regular spiking, fast spiking, intrinsically

Copyright © 2017, the Authors. Published by Atlantis Press.
This is an open access article under the CC BY-NC license (http://creativecommons.org/licenses/by-nc/4.0/).
97

Journal of Robotics, Networking and Artificial Life, Vol. 4, No. 1 (June 2017) 97–101___________________________________________________________________________________________________________and constants in this qualitative model are purely withno physical units. A cubic shaped 𝑣-nullcline is a key toreplicate the spiking dynamics in qualitative models.
Because multipliers are resource-consuming in a digitalarithmetic circuit, the DSSN model adopts a piecewise
quadratic function for the 𝑣 -nullcline so that its
numerical integration step requires only one
multiplication between variables. The architecture of theDSSN circuit was explained in our previous study [11].

mathematical structures with Hodgkin-Huxley-typemodels. Regular-spiking (RS) is a most typical classesof cortical neurons and is characterized by spikefrequency adaptation; the spike frequency decreases
over time in response to a constant stimulus input.
Conversely, Fast-spiking (FS) neurons maintain firing at
a constant frequency. Parameters were optimized toreproduce waveforms of the ionic-conductance models
in [8] and reduce the circuit resource requirements forimplementation.
The remainder of this paper is organized as follows.
Section 2 introduces our neuron model, differential
evolutional algorithm, details of our parameter
optimization procedure, and GPU-based implementation.
The result is shown in Section 3. Section 4 summarizes
the work and suggests ideas for future.

2.2. Differential evolution algorithm
Differential evolution (DE) algorithm [12] is a popularheuristic method to solve optimization problems using areal number function. It is characterized by its mutation
process that uses geometrical location of whole
population. In the following, we explain each step in theDE algorithm. We define 𝑥𝑘𝑖 (𝑗) as the gene j of the i th
individual of the k th generation.
 Initialization
We generate M-dimensional parameter vectors 𝑥0𝑖(𝑖 = 0, … , N − 1), where N and M are the number of theindividuals and parameters of the DSSN model.
Parameter vectors are chosen randomly within theboundary constraints of the DSSN model.
 Mutation𝑖We generate new parameter vector 𝑥𝑚𝑚𝑚
as follows:𝑖∀𝑖 = 1, . . . , N,𝑥𝑚𝑚𝑚= 𝑥𝑘𝑟1 + 𝐹 ∙ (𝑥𝑘𝑟2 − 𝑥𝑘𝑟3 ),where r1, r2, and r3 are integer randomly chosen from[0, N − 1]. F is usually set to 0.5.
 Crossover𝑖
We generate new parameter vectors 𝑥𝑐𝑐𝑐𝑐𝑐
that inherit𝑖genes of 𝑥𝑚𝑚𝑚 with a probability CR.
∀𝑖 = 1, … , N,
∀𝑗 = 1, . . . , M,𝑖(𝑗)𝑥𝑖𝑖(𝑟𝑟𝑟𝑟 ≤ 𝐶𝐶)
𝑚𝑚𝑚𝑖(𝑗) = � 𝑖
𝑥𝑐𝑐𝑐𝑐𝑐𝑜𝑜ℎ𝑒𝑒𝑒𝑒𝑒𝑒,𝑥𝑘 (𝑗)
where CR is usually set to 0.5.
 Selection
A selection is conducted by comparing the cost function
𝑖values of 𝑥𝑘𝑖 and 𝑥𝑚𝑚𝑚
, respectively, as follows:∀𝑖 = 1, . . . , N,𝑖𝑖(𝑗) 𝑖𝑖 𝑓(𝑥𝑐𝑐𝑐𝑐𝑐𝑥𝑐𝑐𝑐𝑐𝑐) ≤ 𝑓(𝑥𝑘𝑖 )𝑖𝑥𝑘+1 = � 𝑖𝑜𝑜ℎ𝑒𝑒𝑒𝑒𝑒𝑒,𝑥𝑘
where 𝑓(𝑥) returns the value of the cost for parameter
vector 𝑥.
A set of procedures including mutation, crossover, andselection is repeated for 10000 times, and the parameter

2. Method
2.1. Digital Spiking Silicon Neuron model
The 3-variable DSSN model is a qualitative neuronmodel that can simulate several classes of neuronal
activities by Euler's method with fixed point operation
including RS and FS [4]. Equations for regular and fast
spiking are given by
𝑑𝑑 ∅= (𝑓(𝑣) − 𝑛 − 𝑞 + 𝐼0 + 𝐼𝑠𝑠𝑠𝑠 ),𝑑𝑑 𝜏
𝑑𝑑 1= (𝑔(𝑣) − 𝑛),
𝑑𝑑 𝜏
𝑑𝑑 𝜀
= (ℎ(𝑣) − 𝑞),𝑑𝑑 𝜏2𝑎𝑓𝑓 �𝑣 − 𝑏𝑓𝑓 � + 𝑐𝑓𝑓 (𝑣 < 0)𝑓(𝑣) ≡ �2(𝑣 ≥ 0),𝑎 �𝑣 − 𝑏 � + 𝑐𝑔(𝑣) ≡ �ℎ(𝑣) ≡ �

𝑓𝑓

𝑓𝑓

2

𝑓𝑓

𝑎𝑔𝑔 �𝑣 − 𝑏𝑔𝑔 � + cgn
2

𝑎𝑔𝑔 �𝑣 − 𝑏𝑔𝑔 � + 𝑐𝑔𝑔
𝑎ℎ𝑛 (𝑣 − 𝑏ℎ𝑛 )2 + chn
2

�𝑣 < 𝑟𝑔 ��𝑣 ≥ 𝑟𝑔 �,(𝑣 < 𝑟ℎ )(𝑣 ≥ 𝑟ℎ ),𝑎ℎ𝑝 �𝑣 − 𝑏ℎ𝑝 � + 𝑐ℎ𝑝
where 𝑣 corresponds to the membrane potential, and 𝑛and 𝑞 are the fast and slow variables, respectively, thatabstractly describe the activity of the ion channels.
Parameter 𝐼0 is a bias constant and 𝐼𝑠𝑠𝑠𝑠 represents theinput stimulus. Parameters ∅, ε, and τ control the time
constants of the variables. Parameters, 𝑎𝑥 , 𝑏𝑥 , and 𝑐𝑥 ,where 𝑥 is 𝑓𝑓, 𝑓𝑓, 𝑔𝑔, 𝑔𝑔, ℎ𝑛, or ℎ𝑝, are constants that
adjust the nullclines of the variables. All of the variables

98

Journal of Robotics, Networking and Artificial Life, Vol. 4, No. 1 (June 2017) 97–101___________________________________________________________________________________________________________vector that produces the lowest value in the costfunction is finally adopted.
2.3. Parameter optimization procedureThe DSSN model for RS and FS classes has 24parameters, exploring whose spanning space is notrealistic with average computing power available tothose who are not specialized in high performancecomputing (HPC). To reproduce the dimension of theparameter space to be explored by the DE algorithm, wefirst explored the parameters that control the dynamics
of the fast state variables 𝑣 and 𝑛 to fit the spike
generation mechanism. Then, the parameters related tothe slow state variable 𝑞 were determined to fit theadaptation of the spike-frequency.
In the DSSN model’s equation, some parameters are acoefficient of a variable. In its circuit implementation, amultiplication between a coefficient and a state variable
was realized by shifters and adders, and we need alarger number of couples of an adder and a shifter as thenumber of digits with value 1 in the fixed point
representation of the coefficient increases. Parameters
were optimized for not only reproducing waveforms of
the ionic-conductance model but also reducing thecircuit size.
In the procedure to determine the parameters that
control the dynamics of the two fast variables, we fixed𝑣 -nullcline to facilitate parameter exploration whileconsidering the dynamics of the spike generation
process. It does not severely restrict the dynamical
property of the model, because the relation between the𝑣- and 𝑛- nullclines is one of the major factors that rulethe dynamical structures in our neuronal model. We
defined the cost function as follows:
𝑓1 (𝑥) = 𝑓w (𝑥) + 𝑘1 𝑓b1 (𝑥),where 𝑓w (𝑥) is a square error between the waveforms of
the DSSN model with parameter vector 𝑥 and a target
waveform. The target waveform was calculated by afast subsystem of the DSSN model in our previous
study [4], which is tuned to reproduce mathematicalstructure of the ionic conductance model in [8]. Notethat minimum and maximum values of waveforms arenormalized to 0 and 1 for comparison. Function 𝑓b1 is anumber of digits with value 1 in the binary fixed-point
expression of the coefficients in 𝑥 . Constant 𝑘1 is tobalance this couple of cost functions. Parameters that
control the dynamics of the fast variables were

Figure 1: Waveforms of the ionic-conductance model (red)and the DSSN model (blue) in response to a weak (left) andstrong (right) step stimulus input. (a)-(c) corresponds toregular spiking neuron cells. (d) and (e) corresponds to fast
spiking neuron cells.

determined to minimize the cost function 𝑓1 by the DE
algorithm.
The remaining parameters related to slow variable 𝑞
were determined to reproduce the spike frequency
adaptation. The cost function is
𝑓2 (𝑥) = 𝑓s (𝑥) + 𝑘2 𝑓b2 (𝑥),where 𝑓s (𝑥) is a square error between spike timings ofthe DSSN model with a parameter vector 𝑥 and that of
the target ionic-conductance model in response to a step
stimulus. Function 𝑓b2 is a number of digits with value 1
in the binary fixed-point expression of the coefficientsin 𝑥, and 𝑘2 is a constant to balance this couple of costfunctions.
The DE algorithm is known to be compatible withparallel computation schemes because of its intrinsic
parallel structure. We accelerated the calculation by theGraphic processing unit (GPU) to obtain betterparameter solution by executing a larger number of DE
steps in a time period. The calculation of the costfunction, in which the differential equations are solvedby the Euler’s method for thousands of times, wasparallelly processed on a GPU. GPU (NVIDIA TeslaK40) –based implementation for calculating costfunctions of 214 parameter vectors consumed 1.578[s],and it was about six times faster than a single-threadprogram on Xeon E5-2670 CPU, which consumed
9.51[s].

99

Journal of Robotics, Networking and Artificial Life, Vol. 4, No. 1 (June 2017) 97–101___________________________________________________________________________________________________________We compiled the DSSN models for Virtex-7
XC7VX690T FPGA using Xilinx Vivado Design Suite.
Device utilization is listed in Table 1. In the column
labeled (f) and (g), the resource usage for RS and FS
settings in [4] is shown. The requirement for LUTs wasreduced down to less than half in all the settings.
4. Conclusion

Figure 2: 𝐶𝑉 − 𝐿𝑉 characteristics of the ionic-conductancemodel (red) and the DSSN model (blue).

In this work, we developed a parameter optimization
method for the DSSN model on the basis of the DE
algorithm. By splitting the optimization process in 24dimentional parameter space into two steps, we could
find parameter sets with which the DSSN model
reproduces the characteristics activities of RS and FScells without using special HPC systems. Cost functions
𝑓𝑏𝑏 (𝑥 = 1, 2) were introduced so that the circuit size isreduced. We also confirmed qualitative similarity
between the ionic-conductance model and the DSSN
model by measuring statistics for the spike timing. In
our future work, we will improve this optimization
method to cover above other neuron classes.

Table. 1 : Device Utilization. The labels (a)-(e) correspond to
those in Fig.1.

3. Result
We found five parameter sets with which the DSSN
model reproduces the activities of the ionic-conductancemodels for three RS and two FS cells whosecharacteristics are slightly different from each other (Fig.
1).
We evaluated the similarity of the spiking patterns
between the ionic-conductance model and the DSSN
model by calculating 𝐶𝑉 and 𝐿𝑣 [13] that are statistics
for the spike sequence while changing the stimulus
intensity (Fig. 2). They are defined as follows:

Acknowledgement
This work was partially supported by JSPS SAKURA
Program and JST PRESTO and CREST.

References
1. Eugene M. Izhikevich, Simple model of spiking neurons,
IEEE Transactions on neural networks, 14(6), pp. 15691572, 2003
2. Andrew A.V. et al., A Model of a Segmental Oscillator in
the Leech Heartbeat Neuronal Network., Journal ofComputational Neuroscience, Vol. 10, pp. 281-302,
2001.
3. T. Kohno and K. Aihara. Digital spiking silicon neuron:
Concept and behaviors in gj-coupled network.
Proceedings of International Symposium on Artificial
Life and Robotics 2007, 2007.
4. T. Nanami and T. Kohno. Simple cortical and thalamic
neuron models for digital arithmetic circuit
implementation, Frontiers in Neuroscience, Vol. 10,No.181, pp. 1—12, 2016
5. T. Nanami, K. Aihara, and T. Kohno, Elliptic andparabolic bursting in a digital silicon neuron model, 2016
International Symposium on NOLTA, Japan, November
27th–30th, 2016
6. Tawfik B, Durand DM, Nonlinear parameter-estimation
by linear association : Application to a 5-parameter𝑛1
𝐶𝑉 = �
�(𝑇 − 𝑇�)2 /𝑇�𝑛−1

𝐿𝑉 =𝑖=1𝑛−113(𝑇𝑖 − 𝑇𝑖+1 )2
�(𝑇𝑖 + 𝑇𝑖+1 )2𝑛−1
𝑖=1

where 𝑇𝑖 denotes the 𝑖 th interspike interval, 𝑇�𝑖 is theaverage of 𝑇𝑖 , and n is the number of spikes in thesequence. Both models' simulation results draw similar
curves on the 𝐶𝑉 − 𝐿𝑉 plane for each cell. The similarityof these characteristics will be estimated based on therequirements specified by the circuits' applications.. Theresults in Figs. 1 and 2 were obtained by numericalintegration with the Euler's method (dt = 2−13 s) where
each variable was expressed by 18-bit fixed point with
14-bit fraction part. The same results are expected to beobserved in the following FPGA implementation,
because this condition is the same as our VHDL codes.

100Journal of Robotics, Networking and Artificial Life, Vol. 4, No. 1 (June 2017) 97–101___________________________________________________________________________________________________________passive neuron model. IEEE Trans Biomed Eng 41:461–469, 1994.
7. Haufler D, Morinc F, Lacaille JC, Skinner FK, Parameter
estimation in single-compartment neuron models using
asynchronization-based
method.
Neurocomputing
70:1605–1610, 2007.
8. M. Pospischil et al.. Minimal Hodgkin-Huxley typemodels for different classes of cortical and thalamic
neurons. Biological Cybernetics, 99(4-5):427–441, 2008.
9. Buhry, L., Sa¨ıghi, S., Giremus, A., Grivel, E., andRenaud, S. (2008). Parameter estimation of the Hodgkin–
Huxley model using metaheuristics: Application toneuromimetic analog integrated circuits. In Proceedingsof IEEE International Conference on Biomedical Circuits
and Systems (pp. 173–176). Piscataway, NJ: IEEE.

10. L. Buhry et al., Automated parameter estimation of the
Hodgkin-Huxley model using the differential evolution
algorithm: application to neuromimetic analog integrated
circuits, Neural computation 23 (10), 2599-2625, 2011.
11. Takuya Nanami and Takashi Kohno, "An FPGA-based
cortical and thalamic silicon neuronal network," Journalof Robotics, Networking and Artificial Life, Vol. 2, No. 4,
pp. 238-242, Mar., 2016.
12. R. Storn, K. Price, A Simple and Efficient Heuristic for
global Optimization over Continuous Spaces, Journal of
Global Optimization, Vol. 11, pp. 341-359, 1997.
13. Shinomoto, S., Shima, K., and Tanji, J. Differences inspiking patterns among cortical neurons. Neural Comput.
15, 2823–2842. 2003

101