Brain machine interface decoding method based on spiking neural network
The present invention discloses a brain machine interface decoding method based on spiking neural network, comprising: (1) constructing a liquid state machine model based on a spiking neural network, the liquid state machine model consists of an input layer, an middle layer and an output layer, wherein, a connection weight from the input layer to the middle layer is W hh , a loop connection weight inside the middle layer is W hh , a readout weight from the middle layer to the output layer is W yh ; (2) Inputting a neuron spike train signal, and training each weight with the following strategy: (2-1) Using STDP without supervision to train the connection weight W hh from the input layer to the middle layer; (2-2) Setting the loop connection weight W hh inside the middle layer by means of distance model and random connection, and obtaining a middle layer liquid information R(t); (2-3) Using ridge regression with supervision to train the readout weight W yh from the middle layer to the output layer, and establishing a mapping between the middle layer liquid information R(t) and the output motion information, and finally outputting a predicted motion trajectory. The present invention can quickly train a model in a relatively short time, predict the arm motion trajectory in real time, and achieve an improvement in efficiency and accuracy.
1. A brain machine interface decoding method based on a spiking neural network, wherein the method comprises the following steps:
(1) constructing a liquid state machine model based on a spiking neural network, the liquid state machine model consists of an input layer, a middle layer and an output layer, wherein, a connection weight from the input layer to the middle layer is W hx , a loop connection weight inside the middle layer is W hh , a readout weight from the middle layer to the output layer is W yh ;
(2) inputting a neuron spike train signal, and training each weight with the following strategy:
(2-1) using Spile-timing-dependent plasticity (STDP) without supervision to train the connection weight W hx from the input layer to the middle layer; wherein in the step (2-1), using the trace of the pre- and post-neuron synapses to achieve STDP optimization the connection weight W hx from the input layer to the middle layer, the trace is the accumulation of the decaying pulse signal, and the formula is:
T
pre
(
t
)
=
U
(
t
)
+
T
pre
(
t
-
1
)
e
-
1
τ
T
post
(
t
)
=
X
(
t
)
+
T
post
(
t
-
1
)
e
-
1
τ
wherein, U and X are pulses fired by the input layer and middle layer, respectively, T pre (t) and T post (t) are the traces accumulated by the pre-synapse and the post-synapse at time t due to pulse firing, respectively, and T parameter is the decay factor, controlling the speed of trace decay when the trace is accumulated, the STDP algorithm using trace optimization is:
W hx ( i,j )= W hx ( i,j )−lr×Δ W hx ( i,j )
wherein, lr is the learning rate, the weight update amount ΔW hx (i,j) between neurons i and j is:
Δ
W
hx
(
i
,
j
)
=
{
W
hx
(
i
,
j
)
-
W
min
,
W
hx
(
i
,
j
)
>
0
and
T
pre
<
T
u
p
W
hx
(
i
,
j
)
+
W
min
,
W
hx
(
i
,
j
)
<
0
and
T
pre
<
T
u
p
W
hx
(
i
,
j
)
-
W
max
,
W
hx
(
i
,
j
)
>
0
and
T
pre
>
T
d
o
w
n
W
hx
(
i
,
j
)
+
W
max
,
W
hx
(
i
,
j
)
<
0
and
T
pre
>
T
d
o
w
n
wherein, T up and T down are the upper and lower target values of the traces of the pre-synapses, respectively, W max and W min are the maximum/minimum value of the weight change caused by STDP optimization, and taking the maximum value of W hx and 0.0;
when the post-synapse neuron fires pulses, the traces of the pre-synapses are detected, when the current accumulated traces of the pre-synapses T pre (t) is greater than T up , it is considered that there is a significant causal relationship between the pre-synapses pulse train and the post-synapse pulse firing in the current period of time, and enhance its corresponding weight, when the T pre (t) is lesser than T up it is considered that the pre-synapses is not directly related to the post-synapse pulse firing in the current period of time, and weaking the weight or disconnect;
(2-2) setting the loop connection weight W hh inside the middle layer by means of distance model and random connection, and obtaining a middle layer liquid information R(t); wherein, the specific process of step (2-2) is:
initializing the connection weight W hh as a two-dimensional matrix with a standard normal distribution:
W hh ˜N (0,1)
scaling the normalized W hh ;
setting the pool in the middle layer as a three-dimensional structure, the pool is composed of multiple cubes with a side length of 1, each neuron is distributed at the vertex of the cube, then judging whether there is a connection based on the Euclidean distance between two points of the neuron, and breaking connection between distant neurons; wherein, the specific way to judge whether there is a connection based on the Euclidean distance between two points of the neuron is:
W hh is a square matrix in [N, N], wherein N is the number of neurons in the middle layer, and the element w ij defines the size of the connection between neuron j and neuron i, whether the connection between neurons i and j obeys probability:
p
(
i
,
j
)
=
Ce
-
D
(
i
,
j
)
λ
wherein, the parameter λ defaults to 2, C defaults to 1.0, D(i,j) function is a distance function, which measures the distance between neurons i and j, using the square of the Euclidean distance, the formula is:
D ( i,j )=(pos( i )−pos( j )) 2 ;
(2-3) using the middle layer liquid information as an input, and using ridge regression with supervision to train the readout weight W yh from the middle layer to the output layer, and establishing a mapping between the middle layer liquid information R(t) and the output motion information, and finally outputting a predicted motion trajectory; wherein, the specific process of step (2-3) is:
when training the connection weight W yh from the middle layer to the output layer, fixing the generated input layer W hx and middle layer W hh connection; when the data arrives, calculating the real-time update of the middle layer neuron membrane potential:
V ( t )= V ( t− 1)+ƒ( W hx U ( t ))+ƒ( W hh X ( t− 1))
wherein, ƒ(W hx U(t)) is the current contribution of the input layer pulse to the post-synapses neuron, and ƒ(W hh X(t−1)) is the current contribution of the pulses fired in the past in the middle layer loop synapse to the current moment;
calculating the output pulse X i (t) of the ith neuron in the middle layer:
X
i
(
t
)
=
{
1
,
V
i
(
t
)
≥
V
thres
0
,
V
i
(
t
)
<
V
thres
,
wherein, V thres is the membrane potential firing threshold, V i (t) is the membrane potential of the i-th neuron in the middle layer at time t, when the neuron membrane potential exceeds the threshold, a pulse is fired, then V i (t)=V is set back to the resting voltage, waiting for the next accumulation-distribution event;
the expression formula of the middle layer liquid information R(t) is:
R
(
t
)
=
X
(
t
)
+
R
(
t
-
1
)
e
-
1
τ
wherein, τ parameter is the decay factor, controlling the speed of trace decay, the smaller the τ, the smaller the impact of the pulses fired in the past on the current moment;
the optimization objective function using ridge regression with supervision training is:
J
(
W
y
h
)
=
∑
t
(
Y
(
t
)
—
W
y
h
R
(
t
)
)
+
λ
R
W
y
h
2
2
wherein, Σ t (Y(t)−W yh R(t)) is used to reduce the deviation between the task objective and the prediction, and λ R ∥W yh ∥ 2 2 is a penalty item, responsible for minimizing the variance of the parameter W yh , λ R is the weight coefficient, which is used to control the proportion of penalty items and is responsible for balancing the prediction deviation and variance;
the trained connection weight W yh is fixed for real-time motion signal prediction:
Y ( t )= W yh R ( t ).
2. The brain machine interface decoding method based on a spiking neural network according to claim 1 , wherein, the value of the decay factor τ is 4 times the time window value.
3. The brain machine interface decoding method based on a spiking neural network according to claim 1 , wherein, the membrane potential firing threshold V thres =0.15, the weight coefficient λ R =20.