Applied Soft Computing 122 (2022) 108811

Contents lists available at ScienceDirect

Applied Soft Computing
journal homepage: www.elsevier.com/locate/asoc

Data augmentation for Convolutional LSTM based brain computerinterface system∗∗

Kahoko Takahashi a , Zhe Sun b , , Jordi Solé-Casals c,d,e , , Andrzej Cichocki f , Anh Huy Phan f ,∗Qibin Zhao g , Hui-Hai Zhao h , Shangkun Deng i , Ruggero Micheletto a ,aCognitive Information Science Lab., Yokohama City University, 22-2 Seto, Kanazawa Ward, Yokohama, Japan
Computational Engineering Applications Unit, Head Office for Information Systems and Cybersecurity, RIKEN, 2-1 Hirosawa,Wako-shi, Saitama, Japanc
Data and Signal Processing Research Group, Department of Engineering, University of Vic–Central University of Catalonia, Vic
08500, Catalonia, Spaind
Department of Psychiatry, University of Cambridge, Cambridge CB2 3EB, United KingdomeCollege of Artificial Intelligence, Nankai University, Tianjin 300071, Chinaf
Laboratory Tensor networks and deep learning for applications in data mining, Skolkovo Institute of Science and Technology, Moscow, RussiagTensor Learning Team, RIKEN Center for Advanced Intelligence Project, Tokyo, Japanh
Alibaba Quantum Laboratory, Alibaba Group, Beijing 100102, Chinai
College of Economics and Management, China Three Gorges University, Yichang, Chinab

articleinfo

Article history:
Received 20 July 2020
Received in revised form 15 March 2022
Accepted 27 March 2022
Available online 8 April 2022
Keywords:
Convolutional LSTM neural network
Data augmentation
Brain computer interface
Empirical Mode Decompositiona b s t r a c t
Electroencephalogram (EEG) is a noninvasive method to detect spatio-temporal electric signals inhuman brain, actively used in the recent development of Brain Computer Interfaces (BCI). EEG’s
patterns are affected by the task, but also other variable factors influence the subject focus on thetask and result in noisy EEG signals difficult to decipher. To surpass these limitations methods basedon artificial neural networks (ANNs) are used, they are inherently robust to noise and do not require
models. However, they learn from examples and require lots of training data-sets. This will increase
costs, need research time and subjects effort. To reduce the number of experiments necessary for
network training, we devised a methodology to provide artificial data from a limited number oftraining data-sets. This was done by applying Empirical Mode Decomposition (EMD) on the EEG frames
and intermixing their Intrinsic Mode Function (IMFs). We experimented on motor imagery (MI) tests
where participants were asked to imagine movement of the left (or right) arm while under EEG
recording. The EEG data were firstly transformed using the Morlet wavelet and then fed to an originallydesigned Convolutional Neural Network (CNN) with long short term memory blocks (LSTM-RNN). The
introduction of artificial frames improved performances when compared with standard algorithms.
The artificial frames become advantageous even when the number of available real frames was onlyof 7 or 8. In a test with two subjects (200 recordings for each subject), we reached an accuracy better
than 88% for both subjects. Improvements due to the artificial data were especially noticeable for theunder-performing subject, whose EEG had lower accuracy. Imagination recognition accuracy was about
89% with 360 training frames, in which 300 were artificially created starting from 60 real ones. We
believe this methodology of synthesizing artificial data may contribute to the development of noveland more efficient ways to train neural networks for brain computer interfaces.
© 2022 Elsevier B.V. All rights reserved.

1. Introduction
Brain computer interface (BCI) devices provide means to sendcommands of various complexity to computers without using
any verbal or tactile actions. These systems, like many others∗ Corresponding authors.
E-mail addresses: zhe.sun.vk@riken.jp (Z. Sun), jordi.sole@uvic.cat(J. Solé-Casals), ruggero@yokohama-cu.ac.jp (R. Micheletto).
https://doi.org/10.1016/j.asoc.2022.108811
1568-4946/© 2022 Elsevier B.V. All rights reserved.

(e.g., fNIRS, ECoG, etc.) are able to read brain signals and recognize specific commands. In simple words, a BCI device connects
between the user and the actual machinery with using brain
signals. A BCI apparatus helps people with handicaps, it may also
contribute to the realization of devices with applications in a
multitude of fields [1–3].
However, brain signals are generated by the electric activityof billions of neurons firing in a random, uncorrelated manner.
The signals are tiny in intensity (few microvolts) and of low frequency, thus they are difficult to analyze and classify. Moreover,K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

to assign meaning to specific EEG events, we have to understandsubjects differences and variability within same subjects [4].
In this sense, it is very important to realize a reliable androbust algorithm able to recognize and classify properly EEG
patterns. Most used in the past are linear classifiers, like linear
discriminant analysis (LDA), regularized LDAs and support vector machines (SVMs) that generally outperforms other classifiers[5–7]. Recently, artificial neural networks (ANNs) emerged as
useful classifier, because they are characterized by robustness tonoise. When properly trained, they can overcome subject variability and can adapt to extremely non-linear behavior [8]. However,to be used for brain–computer interfaces, ANN should be trainedlabeling each EEG pattern to the correct mental state in orderto obtain the optimal ANN response [9,10]. The necessity of agreat amount of calibration data makes the training of ANNs
inconvenient, especially when dealing with subjects of age orwith disabilities. Generating artificial data by randomly sampling
the available dataset is a possible strategy [11].
In this study we propose a new method that uses a Convolutional LSTM network to analyze EEG patterns with a reduced
number of training sets. The subjects are asked to imagine movement of the left or right hand (motor imagery tests) while EEG
are being recorded for two seconds. Each of these training set islabeled accordingly and used for the network learning phase. We
developed an original method to synthesize artificial training data
starting from a limited number of real tests on the subjects andthen using empirical mode decomposition [12,13] (EMD), mixing
in an opportune way the data intrinsic mode functions (IMF) as
described in details in Section 2.3.
The imagination of a motor action can induce recognizable
brain signals. If we analyze the spectral content, we notice theemergence of recognizable frequencies in correspondence to particular mental states [14]. For example, θ rhythm (4–7 Hz) isfound during sleep and low brain activities and β rhythm (13–35 Hz) in alertness and motor activity. Brain waves show not only
rhythmical patterns, but have also more structured behaviors.
For instance, in preparation of a motor action, frequencies slowly
reduce in amplitude. These frequencies show dynamical properties. The recall of cognitive event, emotional events or motor
actions can accelerate/decelerate or suppress/enhance α and βrhythms [15] in a complex way, the phenomenon is called event
related de-synchronization. Experiments confirm that even thesimple perception of a motor action performed by a third person,
can result as suppression of β waves [16] in the EEG signal.
Electroencephalogram events are mapped to specific areas of thebrain. For example, different body movements result to signals in
localized in regions of various shape and extensions on the skull.
The area related to hand movements is especially large and easy
to map with EEG. For this reason, many imagery EEG experiments
are done asking subjects to imagine motor actions involving theirhands or arms. One of the main difficulties in detecting correctly
EEG signals is noise filtering. An approach is to firstly filter each
EEG channel with a band-pass, then a spatial filter is applied,the cluster obtained are projected on positions on the skull tobe identified. This method, or those derived from it, are called
common spatial pattern (CSP) technique [9]. They are widely used
to reduce noise, but noticeably the band-pass filtering used are
obtained by trial and error in a case-by-case fashion.

Fig. 1. The sketch of the structure of the CNN used. Inml represents the image(or frame) extracted from the EEG potential data, n and m are coordinates on
the scalp and l is the frame number. The array W is the convolution filter,
it is smaller and have three components x, y and l. The first two are spatial
components running over the wider range n and m. The latter l represents thetime (frame). This operation extracts the features associated to the filter Wxyl .

Fig. 2. A batch of time dependent EEG data are filtered by a convolution layer(see Fig. 1), the resulting frame is reduced in size by pooling the maximum
value for each four adjacent data as in this simple sketch.

structured in an N × M grid where N and M represent, loosely,
the coordinates along the scalp. In this way, we can interpret theEEG data at a certain instant of time as an image, or frame, ofsize N × M. Since our data are time series coming from multiple
EEG channels, we have a multidimensional input that is of size
N × M × L, where L is the frame number (trials), in a certain EEG
data batch.
The convolution process filters each of these L images with afilter W that runs over the N and M coordinates. The size of thefilter is X by Y by the frame number L. The operation results in a
new array Oij given by the expression:
Oij = Σx Σy Σl Wxyl I(i+x−1)(j+y−1)l + b(1)

The array Oij represents the EEG signal map, convoluted withthe filter W and integrated along all the frames available in
the particular data batch. The data size and CNN structure are
introduced in Table 1. This can be called the feature map due tothe filter W . See Fig. 1 for a simple schematic representation of
the convolution process (the filter W is optimized in order toobtain the best classification performances, in other words thefilter is actually a variable).
The Oij frame is as large as the original data, a pooling layer in
the network is introduced to reduce the size of it. We used themax pooling procedure, each four adjacent data in the frame issubstituted by one single data corresponding to the maximum of
the four (see Fig. 2).
For better optimization with only positive numbers, a rectified
linear unit (ReLU) layer is added. This operator simply sets any

2. Methods
As mentioned above, in this study we want to improve previous methodologies using highly adaptable, robust to noise andtrainable convolutional neural networks (CNN) [17], structured
with convolution layers, pooling layers and rectification layers.
The input of the net at a certain instant is the EEG voltage,2

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Fig. 3. In panel (a) the scheme of the LSTM block. A memory cell passes its value Ct to the next time step depending on calculations involving the values of thegates f , i and z. In (b) the time explicit structure of the network.

negative values to zero and leave untouched all the other data.
With this configuration the network can learn to classify a special
set of training data. These are data that are EEG recordings fromsubjects, with the correct left or right imagery label manuallyinserted by a human.
2.1. The long short-term memory cell
The so called Long short-term memory (LSTM) [18] is a particular recurrent neural network (RNN) that is used to analyze
time series data. An RNN can memorize previous data in the timeseries and use them for the optimization process. Compared to asimple RNN, a LSTM network has longer memorization time andthis results in successful applications to natural language processing and speech recognition. It is widely used for time series andsignals processing. In this study we use a LSTM network with theforget gate introduced by Gers [19].
Calling the input signal at time t as Xt , the corresponding LSTM
block output is ht . The latter depends not only on Xt but also on
the previous output ht −1 . This is valid for all point in time Xt ,Xt −1 , Xt −2 etc. The LSTM block involves memory gates that havea function to selectively pass information to neurons. During thelearning phase, the gates actually extend or shorten the memory
length to maximize the fit with the desired output.
The scheme of the LSTM block we used is as in Fig. 3(a).
The first output gate acts as a standard RNN network, so its
output is:ot = σ (Wo xt + Ro ht −1 + bo )Fig. 4. The core structure of an artificial neural network. Layer (x − 1) and layer
x are fully connected. The output of one layer becomes the input of the next.
The output O(x−2) of layer (x − 2) is the input of neurons in layer (x − 1) and soon. Eq. (8) transfers values from one layer to the next through the weights Wij .

The function of the memory cell is to keep the ct value or not.
Depending on this decision ct −1 calculation is executed and finallythe output is:zt = ot × tanh(ct )As in any RNN the LSTM block acts on each Xt , as schematicallyshown in Fig. 3(b).
2.2. Optimization
We use here a fully connected artificial neural network withfeatures extraction filters and memory blocks. The specific typeof memory we use is called long short-term memory (LSTM).
Hereafter we describe briefly some details about the structure of
our network.
The fully connected network is structured in layers. Withineach layer, neurons are not connected to each other. Howevereach neuron of the layer is feeding and receiving information to
all the neurons of the adjacent layers (Fig. 4). The output of eachlayer is transferred to the next layer through the weights W ij. In
this way the output for the neuron i of layer x is:(2)

where σ is a sigmoid function, Wo and Ro are weight arrays andbo is a bias constant. The forget gate ft and the input gates it andzt have similar structure:
ft = σ (Wf xt + Rf ht −1 + bf )(3)

it = σ (Wi xt + Ri ht −1 + bi )

(4)

and
zt = tanh(Wi xt + Ri ht −1 + bi )Oxi = L (Iix )
(5)

(8)

−θ
θ is athreshold and N the number of neurons in the layer. The logistic
function L is 1+1e−x .
The optimization of the weights Wij is obtained by the backpropagation method [20]. The cost function that estimates theerror between a target Yi and the actual output Oi is the crosswhere Iix is the neuron activation ΣiN Wijx−1 Oxi −1

These three elements are connected to the memory cell ct throughelement-wise addition or multiplication modules as in Fig. 3. Thatis:ct = it × zt + ft × ct −1

(7)

(6)
3

x−1
. Herei

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

entropy derived from the logistic regression theory [21]. For eachlayer x we define the cost as:
E = −Σi [Yix lnOxi + (1 − Yix )ln(1 − Oxi )](9)The new weights Wij are optimized by adding a small term
proportional to the gradient:dE
dWijx−1=

dE dOxi
dOxi

dIix

dIix
dWijx−1(10)

From Eq. (9) we can calculate these partial derivatives :
dE
dOxi
dOxi
Iix
dIix
Wijx

Yx

1 − Yix

Oi

1 − Oxi= − ix +(11)= Oxi (1 − Oxi )= Oxi −1

Fig. 5. The positions and channel naming of the 62 electrodes used in theexperiment.

multiplying these together few terms go away and finally we
have:
dE
(12)
= Oxi −1 (Oxi − Yix )dWijx−1

2.4. Synthesis of artificial frames
We use the data obtained by Phan [25] on two subjects (S1and S2) using an EEG of 62 electrodes placed like in Fig. 5. The
measure was made for two seconds, 200 times per subject. The
EEG sampling rate was set to 500 Hz, so each dataset was an arrayof size 1000 × 62 (500 Hz for 2 s, 62 channels).
Of the 200 total datasets, we used 60 of them as training, andthe remaining 140 for validation. Since we intend to generate
artificial training data with our original method, we call these 60
datasets the real training set and the other that we will create
artificial training sets.
The artificial frames are constructed using an iterative procedure. As introduced above, firstly we decompose the original
signal EEG(t) in its empirical modes using the intrinsic mode function IMF (t) and the trend component rn (t) [12]. E is an operator
that does the spline interpolation of all the maxima and all theminima of the signal and the apex T represents the transpose ofa vector, the decomposition is done accordingly to the following
steps:

We multiply this term by the logistic derivative (second in
Eq. (11)), so ∆Wijx−1 will be bigger around the center and minimalon the extremes of it:∆Wijx−1 = −ε Oxi −1 (Oxi − Yix )Oxi (1 − Oxi )(13)

where ε is a variable learning rate given by:ε = ε01
1 + β Ne

(14)

Here ε0 is the initial learning rate and β the learning attenuation.
Ne are the epochs, the number of times the learning algorithm
have to pass through the entire training dataset to update its
weights (an epoch can be divided in batches of data extracted
at random from the set).
To speed up the convergence we introduce a momentum [22,23] term α . Overall, considering that

∆Wijx−1 = Wijx−1 (t + 1) − Wijx−1 (t)at time step t the new Wijx (t + 1) will be:

1. A function s(t) = ri−1 (t) is defined, i begins with i = 1 andr0 (t) = EEG(t).
2. Calculate the upper and lower envelopes (E u [s(t)], E l [s(t)])
3. We define the average of the upper and lower envelope asm(t) = (E u [s(t)] − E l [s(t)])/2.
4. The intrinsic mode function (IMF) candidate h(t) is definedas the difference between the average m(t) and the signal:h(t) = s(t) − m(t).
5. s(t) is reassigned as s(t) = h(t) and the procedure returns
s(t)T s(t)to point (2) until the ratio m(t)T m(t) reaches 40 dB.

Wijx−1 (t + 1) = Wijx−1 (t) −ε Oxi −1 (Oxi − Yix )Oxi (1 − Oxi ) +α Wijx−1 (t) (15)2.3. Empirical mode decomposition
To improve network training efficiency and reduce the number of tedious and time-consuming tests with humans, we created artificial EEG frames from a number of real motor imaginary
experiments. Empirical mode decomposition has been used for
EEG signal processing [13]. We decomposed the EEG frames using
Empirical Mode Decomposition and created new ones intermixing their intrinsic mode functions (IMF) [24]. A number of labeledhomogeneous EEG data are decomposed and randomly mixed,
creating new artificial data, each different from one another, but
still compatible to the same label because obtained from different
sets of the same motor-imagery contents.
The time series EEG(t) is decomposed in the sum of n intrinsic
mode functions and a trend component rn (t) as
EEG(t) = Σkn=1 IMFk (t) + rn (t)6. Then the first intrinsic mode function IMFi (t) = h(t) isfound.
7. The residual is defined as the difference between the original signal and the current IMF ri (t) = ri−1 (t) − IMFi (t)8. The procedure is repeated from (2) until the ratio
ri−1 (t)T ri−1 (t)ri (t)T ri (t)reaches 50 dB.

In this way the EEG(t) pattern is decomposed in n intrinsic mode
functions plus a residual, as in Eq. (16). To create new artificial
EEG frames, we choose random IMFs and use Eq. (16) to build
new artificial frames. We have, of course, to mix frames belonging
to the same group of imagery data (left hand or right one). In
Fig. 6 we show a simplified scheme of the procedure.

(16)

The main characteristics of the IMFk (t) function are that (1) the
difference between the number of peaks and the number of zerocrossing points should be one or zero, (2) the envelope of thepeaks and minimum points on average should be zero.
4

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Fig. 6. A sketch on how the IMF decomposition works. In (a) is shown the EEG time series decomposed to N intrinsic mode functions. In (b) an artificial frame iscreated composing various intrinsic mode functions chosen from random EEGs of the same imagery dataset (in the figure IMF i − j means the jth intrinsic mode
function from the ith EEG).

2.5. Time resolved frequency analysis
EEG signals can be analyzed by their frequency and theirposition in space. Since a simple Fourier analysis would criticallydepend on the time-window length, we used a Morlet based
wavelet analysis. The so-called mother wavelet is a function thatis variable in shape, enlarging at low frequencies and getting
sharper at higher. In this way the analysis is balanced and moresuitable for wide frequency range signals as our imagery EEG timeseries. The general wavelet form we used is the following:
1Wψ (a, b) = √a∫ ∞s(τ )ψ ∗ (

τ −b
a

−∞

)dτ

(17)

where b is the time reference and a is the scaling parameter. As
the mother wavelet function we use the Morlet function below:
1Fig. 7. The shape of the mother wavelet Morlet, Eq. (18), with fc = 1 Hz andβ = 1 s2 .

2

− tβ

· ei2π fc t · eψMorlet (t) = √
πβ

(18)

here fc is the center frequency and β is the band width parameter.
The smooth shape of this function is typical of biological andphysiological reactions, often used in brain research (see a plot
in Fig. 7).

of them (Fig. 9). Each time series is a 1000 points vector (2 sof EEG recording at 500 Hz sampling rate). There are 62 EEGchannels so the network input is a bi-dimensional array of size
62 × 1000. After the wavelet/Fourier transformation it becomesa 23 × 50 × 62 three dimensional array that is fed to the CNN
network. The structure of the CNN network is given in Table 1.
In our network the kernel filter (that is the filter operatingon the input bidimensional array of data) is of size 3 × 5 with
1 × 1 stride size (the strides are the number of pixels by which
the kernel filter moves over the data, 1 × 1 means that thefilter moves one step vertically and one horizontally) and has
100 filters. The pooling kernel is instead 3 × 3 with 2 × 2
stride size. The network used here has fewer layers than thatof Kwon [28], and a reduced number of layers improves computational efficiency and this has no detrimental effects on ourstudy since we have a lower number of training data. After thesecond pooling, the output is transformed to a one-dimensional
vector of length 3600. This is fed to a fully connected networkwith two output. The two labels values are the vector [1, 0] or [0,1], corresponding to the two motor imagery thought of moving
the right or left arm respectively. In the network training phase,

2.6. The convolutional neural network
We used two different network configurations in this study[17,26]. The first is a CNN network that categorizes the resultof EEG wavelet analysis, the second it is a Convolutional LSTM
network classifying the spectrum of the same EEG analysis [27](see Fig. 8).
The CNN used is of the same structure as the one used by
Kwon for the recognition and classification of emotions by electroencephalographic analysis [28]. Before the signal is fed to
the network it is converted to its time–frequency distribution(TFD) described above. The conversion is done through Fourier orwavelet transformation where an 8–30 Hz spectrum is shown for
an interval of 2 s. Given our setup constraint of 500 Hz sampling
rate (corresponding to 2 ms integration time), we have chosen
20 samples to calculate each spectrum. Thus, each wavelet timewindow has a duration of 40 ms and there are a total of 50
5K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Fig. 8. A simplified sketch of the CNN used in this study. The input is a three dimensional array in which the first dimension is the wavelet transformation output
is a vector of 23 values (8–30 Hz interval divided in 23 bins), the second is the number of time-windows (of 40 ms duration for a total of 2000 ms) and the third
the EEG channels on the scalp (Nch = 62 channels).

Fig. 9. The wavelet (top panel Fourier) transformation results in images of size 50 × 23. The normalized time-series is time-windowed in chunks of 40 ms duration,with 20 samples in each of them. The integration time is 2 ms (500 Hz sampling rate) and since the total recording time is of 2 s, there are a total of 50 time-windows,
each corresponds to pixels on the horizontal axis. The transformation spectrum results in 23 points ranging from 8 to 30 Hz on the vertical axis.
Table 1The CNN structure and data sizes as implemented with the KERAS library [29].
The bi-dimensional convolutional network output is filtered by a rectified linear
unit (ReLU) and fed to a pooling layer. The procedure is repeated and the result
is transformed in a flat vector with the Flatten operator. This vector finally is fed
to a regular ANN network composed by two densely connected layers rectified
by a ReLU unit. Weights are optimized using the Adam optimizer. The layer typeis a mnemonic and in parenthesis is the actual library function used.
Layer (type)Output shapeconvolution_1 (Conv2D)
activation_1 (ReLU)
pooling_1 (MaxPooling2D)
convolution_2 (Conv2D)
activation_2 (ReLU)
pooling_2 (MaxPooling2D)
flatten_1 (Flatten)
layer_1 (Dense)
activation_3 (ReLU)
layer_2 (Dense)(23, 50, 100)
(23, 50, 100)
(11, 24, 100)
(9, 20, 100)
(9, 20, 100)
(4, 9, 100)
(3600)
(62)
(62)
(2)

a function proportional to the output distance from these labels
will be used to optimize all the CNN network variable parameters.
To train the CNN network we used 60 of the 200 EEG recording
to create artificial training data. The remaining 140 measurements were used for the network testing and validation. In a first
test, using 30 EEG measurements, we created two sets of 20 and
30 artificial training data for a total of 50 (30+20) or 60 (30+30)
training sets. In another test, the set of 60 EEG recording wasused to create sets of 120, 300 and 600 artificial data. The CNN in
this case was trained with a total of 180, 360 and 660 sets. Oncethe training was complete, the CNN performance was tested andvalidated for subject 1 and 2 using the 140 untrained recording.
2.7. CNN and LSTM network
With the single CNN structure it is not possible to enhance thetime dependence of the EEG signal, thus here we constructed anenhanced CNN-LSTM hybrid network as shown in Fig. 10. For theexperiment we used the 45 EEG channels shown in Fig. 11. Like in6

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Fig. 12. An example of the normalized EEG input, organized in a (5 × 9) array.
Channel labels correspond to Fig. 11.
Table 2
Parameters we are using in our CNN+LSTM neural network.
Fig. 10. A sketch of the Convolutional LSTM used in this study.

Layer (type)Output shape

time_distributed_1 (Time Distributed)
time_distributed_2 (Time Distributed)
lstm_1 (LSTM)
dense_1 (Dense)
activation_1 (Activation)
dense_2 (Dense)
activation_2 (Activation)(50, 7, 3, 128)
(50, 2688)
(256)
(16)
(16)
(2)
(2)

Table 3
The accuracy for the CNN with a short time-window Fourier transformation(STFT) and wavelet transformations with different fb . In gray the best results,
the three trials are performed with 60 training sets and 140 evaluation sets
(subject 1).
Accuracy for each trial

fb = 1
fb = 2
fb = 3fb = 4fb = 5fb = 6
STFT

Average

#1

#2

#3

87.9%
88.6%
88.6%
90.0%
90.0%
88.6%
81.4%

88.6%
89.3%
86.4%
88.6%
89.3%
87.1%
83.6%

88.6%
88.6%
88.6%
85.0%
87.1%
88.6%
86.4%

88.4%
88.8%
87.9%
87.9%
88.8%
88.1%
83.8%

Fig. 11. The 45 channels used in the CNN-LSTM experiments are shown withthe filled dots.

3. Results

the previous CNN analysis, after the wavelet transformation, we
have a tensor of dimension (23 × 50 × 45). The 45 channels are
spatially organized in a (5 × 9) array (see Fig. 11). We have one ofthis array for each of the 23 resulting frequencies and each of the
50 time-windows. The data tensor is of size (9 × 5 × 23 × 50).
In Fig. 12 is shown an example of one of the 5 × 9 EEG frames. In
the CNN+LSTM configuration, the first 40 ms chunk of input data
is given to the CNN and subsequently fed to the LSTM network.
The output of it is given to a standard fully connected ANN. Likein the previous case, the network final layer is composed by two
neurons, one corresponds to the ‘‘left arm’’ movement label, theother for ‘‘right arm’’ one. The optimization process proceeds
like in any CNN, updating the network weights accordingly toan optimization algorithm. Once the network is optimized, we
evaluate its performances with untrained input data and compare
the network output with the subject actual imagery label (seeTable 2).

Firstly we compare the network accuracy for the CNN withthe time windowed Fourier transformation (short time Fouriertransformation, STFT) against the one with the wavelet with fb =1 up to fb = 6. Training is done with 200 sets of imagery data, 140
are used for training and 60 are reserved for accuracy evaluation.
Results are shown in Table 3 for subject 1 and Table 4 for subject
2.
Noticeably, using wavelet transformation gives better accuracy
up to 5% for subject 1 and 1% for subject 2 (grayed in the tables).
We did not notice strong dependence on the parameter fb , sohereafter, we use a network with wavelet spectral pre-processing
with fb = 6.
To test our CNN+LSTM network we firstly train the systemwith 60 real EEG labeled data. Then we create 120, 300 and 600
artificial data with the methodology explained above and use
them again as training sets. In total, the network is trained with180, 360 and 660 sets.
Then, we use from 6, 12 and 30 real EEG data to create artificial
frames in order to obtain a total of 60 training sets. We train the7

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Table 4
The accuracy for the CNN with a short time-window Fourier transformation for
subject 2.
Accuracy for each trial

fb = 1
fb = 2
fb = 3fb = 4fb = 5fb = 6
STFT

Average

#1

#2

#3

85.0%
85.7%
87.1%
87.9%
86.4%
85.7%
86.4%

88.6%
87.1%
85.0%
88.6%
88.6%
87.1%
85.7%

85.7%
85.0%
87.9%
85.7%
85.7%
87.1%
87.1%

86.4%
85.9%
86.7%
87.4%
86.9%
86.7%
86.4%

CNN+LSTM network again with those. This procedure is done forsubject 1 and subject 2.
In the next four tables we compare the accuracy of a CNN
network with a CNN with a long short-term memory cell (LSTM).
In Tables 5 and 6 are shown the results for the bare CNN network
for subject one and two respectively. Then the same results areshown for the CNN+LSTM network again for the two subjects,
Tables 7 and 8.
4. Discussion
We compared our results with the work of Cichocki and Phan(2010) where similar motor-imagery brain computer interface
datasets were classified using the CSP method, an advanced tensor decomposition approach [25]. We notice that the CNN+LSTM
network shows an increase of 5.18% for subject one and a decrease of 3.57% for the second subject. However, as it is evident
from the results in Tables 5 and 6, when we add 600 artificial training data with our method the accuracy raises for bothsubjects. In particular, the accuracy of 88.21% for subject 2, and89.46% obtained with subject1, are excellent values similar to previous motor-imagery results [25]. Compared to the CNN network
performance, our proposed CNN+LSTM system increases accuracyof 3.21% for subject one and 0.36% for subject two ( Table 9).
Considering both our tests with CNN and CNN+LSTM, the additional EMD artificial data in the network training produce betteror, depending on subjects, slightly better accuracy compared topast literature accuracy on EEG motor-imagery tests.

Fig. 13. The effect of EMD frames on the CNN (light color) and CNN+LSTM(dark color) networks. Accuracy compared for both subjects (see Tables 5–8). The introduction of the LSTM produces improved precision in almost allcases. On the vertical axis the accuracy is shown in percentage, averaged overfour experiments with identical conditions and different random seeds. On
the horizontal axis the numbers of total training frames are shown, below inparentheses the numbers of artificial frames are included in the count.

4.1. CNN and CNN+LSTM compared
In Fig. 13 we show the effect of our artificial frames on the
training process, for both networks CNN and CNN+LSTM compared. For subject one, the first 30 and 60 real EEG frames (zeroartificial training sets) result in an increase for both CNN and
CNN+LSTM models. However, with 30 frames CNN ha slightly
less accuracy than CNN+LSTM and the opposite happens with 60
training frames. A similar behavior is observed for subject two:for 30 training sets CNN+LSTM has better performances, for 60
the CNN network improves and the two models have almost thesame accuracy. We can say that for a low number of training data,
the CNN+LSTM performs better than the CNN, whereas when thedata are abundant the two networks appear to be equivalent. The
LSTM memory appears to be critical when the training data arescarce, whereas this become unimportant when the training isplenty.

Fig. 14. The accuracy of the CNN+LSTM network response in function of real
training data. The accuracy improves with the number of training sets, forsubject 1 this effect is especially dominant. Data are extracted from Tables 7
and 8 (no artificial data).

For both subjects we observe a lower accuracy when the network training is scarce. Increasing the number of training set, as
expected, improves the overall accuracy, nevertheless it remains
noticeably lower for subject 1. Since for equal number of training
frames, the network does not perform equally with the two
subjects, we can speculate that subject 1 was distracted or lessfocused on the task. In other words, he/her had other thoughts

4.2. The dependence on training sets number
In Fig. 14 we show the CNN+LSTM accuracy response whentrained with an increasing number of real EEG frames.
Clearly the performance gets better with training, however the
improvement is reduced for higher number of training frames.
8K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Table 5Classification accuracy with a CNN network (subject 1). The simulation is repeated 4 times using different random
seeds. The number of real and artificial data are changed as in the table. The validation tests are done with 140
EEG recording. The grayed cell emphasizes the best accuracy.
Number of training data

Accuracy of each validation

Real data

Artificial data

Total

#1

#2

#3

#4

30
30
30
60
60
60
60

0
20
30
0
120
300
600

30
50
60
60
180
360
660

80.00%
83.57%
84.29%
88.57%
87.86%
87.86%
90.00%

80.71%
78.57%
85.71%
87.14%
89.29%
87.86%
90.00%

82.14%
82.14%
80.71%
88.57%
90.00%
87.86%
90.00%

78.57%
80.71%
83.57%
87.86%
88.57%
85.00%
87.86%

AverageSESD

80.36%
81.25%
83.57%
88.04%
88.93%
87.14%
89.46%

1.29%
1.85%
1.82%
0.59%
0.80%
1.24%
0.93%

0.64%
0.92%
0.91%
0.30%
0.40%
0.62%
0.46%

Table 6Classification accuracy with a CNN network (subject 2). Data are shown in the same fashion as Table 5.
Number of training data

Accuracy of each validation

Real data

Artificial data

Total

#1

#2

#3

#4

30
30
30
60
60
60
60

0
20
30
0
120
300
600

30
50
60
60
180
360
660

72.14%
76.43%
75.00%
85.71%
85.00%
87.14%
87.86%

70.71%
75.51%
72.86%
87.14%
87.14%
88.57%
88.57%

70.00%
75.71%
74.29%
87.14%
86.43%
87.86%
87.86%

70.00%
75.00%
72.86%
85.71%
85.71%
88.57%
88.57%

AverageSESD70.71%
75.66%
73.75%
86.43%
86.07%
88.04%
88.21% (a)

0.87%
0.51%
0.93%
0.71%
0.80%
0.59%
0.36%

0.44%
0.26%
0.46%
0.36%
0.40%
0.30%
0.18%

Table 7
Classification accuracy with Convolutional LSTM neural network (subject 1). The number of validation sets is 140.
The grayed cell emphasizes the best parameters. The data in bold are those used in the plot of Fig. 13.
Number of training data

Accuracy of each validation

AverageSESD

64.29%
77.14%
78.57%

64.82%
76.25%
81.25%

1.85%
1.70%
1.63%

0.92%
0.85%
0.81%

81.43%

82.14%

80.89%

1.37%

0.69%

85.00%
65.71%
76.43%

85.00%
56.43%
77.14%

84.29%
64.11%
76.96%

0.71%
4.58%
1.06%

0.36%
2.29%
0.53%

83.57%

82.86%

83.75%

0.78%

0.39%

57.14%
78.57%
84.29%

59.29%
80.71%
80.71%

61.79%
78.57%
83.04%

3.86%
1.52%
1.46%

1.93%
0.76%
0.73%

86.43%

85.71%

85.89%

0.31%

0.15%

85.71%

86.43%

86.25%

0.31%

0.15%

85.71%
87.14%
88.57%
88.57%

85.71%
86.43%
87.14%
88.57%

86.07%
86.61%
88.93%
87.32%

0.36%
1.06%
1.29%
1.28%

0.18%
0.53%
0.64%
0.64%

Real data

Artificial data

Total

#1

#2

#3

#4

6
12
20

0
0
0

6
12
20

65.71%
75.00%
81.43%

67.14%
74.29%
82.14%

62.14%
78.57%
82.85%

30

030

81.43%

78.57%

40
6
12

0
44
38

40
50
50

83.57%
65.71%
75.71%

83.57%
68.57%
78.57%

30

20

50

85.00%

83.57%

6
12
20

56
48
40

60
60
60

67.14%
78.57%
84.29%

63.57%
76.43%
82.86%

30

30

60

85.71%

85.71%

40

20

60

86.43%

86.43%

60
60
60
60

0
120
300
600

60
180
360
660

86.43%
85.00%
90.71%
86.43%

86.43%
87.86%
89.29%
85.71%

that were perturbing the motor-imagery EEG patterns, this wasreflected on the network ability to recognize them (Fig. 14).
In Fig. 15 we plot the reciprocal of the training sets number against the accuracy. For subject 1, there is a clear linear
correlation, with regression coefficient of R2 = 0.97. In this reciprocal space, the negative inclination and position of the intercept
suggests the theoretical accuracy limit for an infinite number oftraining sets, for subject 1 this is about 87%. In the case of thesecond subject, we observe the same negative trend, however
linearity is less prominent, especially for the last two points ofthe curve. Those correspond to tests with less training data where
the network seems to produce much better accuracy than withsubject 1. However, considering only the first four points (morethan 20 training sets) the linear trend is clear, it shows a fit of
R2 = 0.98 and an intercept of about 93%. This is suggesting abetter accuracy for this subject in the theoretical case of a great
number of training sets.
Using the same reciprocal plot, we can visualize the effect of
artificial frames on the network accuracy (see Fig. 16). Here weuse an increasing number of real EEG frames used for training
the network and add a variable number of artificially generated
frames in order to reach a total of 60 training sets. The plot shows
the accuracy using only real frames with filled dots, and theaccuracy with the real and artificial ones with asterisks. Clearly,
the artificial frames enhance the overall accuracy of the network,
the only exception is when the number of real frames is verysmall. The inversion point is indicated by the arrows (0.13–0.14
corresponds to about 7–8 frames). This inversion is expected
since the artificial frames are created from real EEG data. Whenthe number of those is limited, the new frames are synthesized
by a small number of intrinsic functions and the random mixing
these in the empirical mode decomposition presumably ends up
in reproducing a bad copy of the original EEG frames, actually degrading the network ability to classify properly the EEG patterns.
On the other hand, when the number of available real EEG frames
surpasses the threshold approximately indicated by the arrows,
the effect of the intrinsic mode functions intermixing begins to
introduce positive effects on the CNN’s classification ability.
9K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Table 8Classification accuracy with Convolutional LSTM neural network (subject 2). The grayed cell emphasizes the bestparameters.
Number of training data

Accuracy of each validation

CSP
CNN
CNN + LSTM

SE

SD

86.43%
83.57%
82.14%

85.00%
85.18%
83.75%

1.43%
1.06%
1.78%

0.71%
0.53%
0.89%

85.71%

86.43%

86.61%

0.78%

0.39%

87.86%
90.00%

88.57%
90.00%

88.93%
89.46%

1.07%
0.59%

0.54%
0.30%

87.14%

90.00%

88.04%

1.17%

0.59%

85.71%
87.14%
85.00%

85.00%
87.14%
83.57%

83.21%
87.68%
85.36%

2.79%
1.37%
1.29%

1.39%
0.69%
0.64%

86.43%

88.57%

87.14%

1.52%

0.76%

90.00%

89.29%

89.64%

0.36%

0.18%

91.43%
90.71%
90.71%
86.43%

90.00%
90.71%
92.14%
86.43%

90.36%
90.36%
90.54%
86.25%

0.80%
0.62%
1.06%
1.28%

0.40%
0.31%
0.53%
0.64%

Artificial data

Total

#1

#2

#3

#4

6
12
20

0
0
0

6
12
20

83.57%
85.71%
84.29%

83.57%
85.00%
82.14%

86.43%
86.43%
86.43%

30

030

87.86%

86.43%

36
40

0
0

36
40

90.71%
89.29%

88.57%
88.57%

30

20

50

87.14%

87.86%

6
12
20

54
48
40

60
60
60

83.57%
90.00%
85.71%

78.57%
86.43%
87.14%

30

30

60

88.57%

85.00%

40

20

60

89.29%

90.00%

60
60
60
60

0
120
300
600

60
180
360
660

89.29%
90.71%
90.00%
87.86%

90.71%
89.29%
89.29%
84.29%

Table 9The comparison of our results with CSP method. The
CNN+LSTM approach improved the accuracy for bothsubjects.
MethodAverageReal data

to synthesize artificial EEG frames from subject’s real recording,and evaluated the effect of them on the network performances.
The main findings of this study are the following:

Accuracy• The CNN and CNN+LSTM networks we proposed performed

Subject 1

Subject 2

82.86%
88.04%
86.07%

90.00%
86.43%
90.36%

similarly or better than previous method used in literature
for the same motor-imagery tasks. In particular, our results compare well with CSP, having, depending on subjects,
similar or better results (Table 9).
• The convolutional network with long short-term memory
blocks (CNN+LSTM) improves results of a similar network
without the LSTM block. This is apparent especially when
the number of available training data is less (Fig. 13).
• The addition of artificial frames tends to improve the accuracy of a CNN+LSTM network, and this improvement isprone to reach a plateau for greater number of artificial
frames added (Fig. 13).
• The (accuracy) vs (reciprocal of training sets) space shows aquasi linear behavior that can be used to infer the network
accuracy limits (Figs. 15 and 16).
Our study demonstrates for the first time that with a CNN+LSTM
network, it is possible to improve the EEG pattern recognition byintroducing artificial frames realized with intermixing intrinsic
mode functions in an empirical mode decomposition process.
This helps in better training neural networks for classification of
motor-imagery tasks, with less training sets and better accuracy.
Our work may pave the way for the realization of next generation
brain computer interfaces or other future devices that rely uponbrain-wave analysis.

Fig. 15. The network classification accuracy against the reciprocal of the numberof training sets. The intercept with the vertical axis indicates the theoreticalaccuracy limit for an ideal infinite number of training sets. Subject1 seems toconverge to a maximum accuracy of about 87%–88%. However, the other subject
data show linearity only when the network is well trained. In fact the lasttwo points corresponding to 12 and 6 training sets only seem to depart fromlinearity. The regression line excluding these two suggests a theoretical limit ofabout 93%.

CRediT authorship contribution statement
Kahoko Takahashi: Conceptualization, Methodology, Investigation, Software, Writing. Zhe Sun: Methodology, Investigation, Supervision, Software, Reviewing and editing.
Jordi Solé-Casals: Methodology, Investigation, Supervision,
Software, Reviewing and editing. Andrzej Cichocki: Methodology, Supervision. Anh Huy Phan: Methodology, Software.
Qibin Zhao: Methodology, Data collection. Hui-Hai Zhao:
Methodology, Investigation. Shangkun Deng: Methodology.
Ruggero Micheletto: Methodology, Investigation, Supervision,
Software, Reviewing and editing.

5. Conclusions
In this study we compared a convolutional neural network(CNN) and one with a long short-term memory block (CNN+LSTM)
in the task of recognizing brain wave spatial–temporal patterns
in motor-imagery experiments. Because of the scarce availabilityof real EEG recordings on human subjects, we proposed a method
10

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811

Fig. 16. A comparison of the CNN+LSTM network accuracy with and without artificial frames. The dotted line and ‘‘*’’ symbol display the accuracy when artificial
frames are used (S1+ and S2+). These are advantageous for both subjects, however improvement is especially noticeable for the subject whose EEG shows lower
accuracy (subject 1). In that case the artificial frames bring the network performance to the same level of subject 2 (>85% for 60 or more frames). The arrows indicate
the point where the use of artificial frames becomes advantageous, about 0.13–0.14, corresponding to more than 7–8 frames. On the horizontal axis we have thereciprocal of the training set number from real data(1/6, 1/12, 1/20, 1/30, 1/40 and 1/60). The number of artificial frames added to those is enough to reach 60
frames in total with the exception of 60 real frames where we added 300 artificial ones (6+54, 12+48, 20+40, 30+30, 40+20 and 60+300). Data are extracted from
Tables 7 and 8.

Declaration of competing interest[10] K. Shiu, S. Alok, M. Kabir, T. Tatsuhiko, A deep learning approach formotor imagery EEG signal classification, in: 2016 3rd Asia-Pacific World
Congress on Computer Science and Engineering, APWC on CSE, IEEE, 2016,
pp. 34–39.
[11] S. Soman, S. Saxena, et al., Eigensample: A non-iterative technique foradding samples to small datasets, Appl. Soft Comput. 70 (2018) 1064–1077.
[12] C. Josep, O. Rupert, G. Christoph, S. Jordi, A new method to generate
artificial frames using the empirical mode decomposition for an EEG-based
motor imagery BCI, Front. Neurosci. 12 (2018) 308.
[13] S.J.D.J. Z. Zhiwen, C. Andrzej, Y. Zhenglu, S. Zhe, A novel deep learning
approach with data augmentation to classify motor imagery signals, IEEE
Access 7 (2019) 15945–15954.
[14] P. Nunez, E. Nunez, R. Srinivasan, A. Srinivasan, O.U. Press, Electric Fieldsof the Brain: The Neurophysics of EEG, Oxford University Press, 2006.
[15] P.J. A., The functional significance of mu rhythms: translating ‘‘seeing’’ and‘‘hearing’’ into ‘‘doing’’, Brain Res. Rev. 50 (1) (2005) 57–68.
[16] G. Pfurtscheller, F.L. da Silva, Event-related EEG/MEG synchronization anddesynchronization: basic principles, Clin. Neurophysiol. 110 (11) (1999)
1842–1857.
[17] S.R. Tibor, S.J. Tobias, F.L.D. Josef, G. Martin, E. Katharina, T. Michael,
H. Frank, B. Wolfram, B. Tonio, Deep learning with convolutional neuralnetworks for EEG decoding and visualization, Hum. Brain Map. 38 (11)(2017) 5391–5420.
[18] S. Hochreiter, J. Schmidhuber, Long short-term memory, Neural Comput.
9 (8) (1997) 1735–1780.
[19] F.A. Gers, J. Schmidhuber, Recurrent nets that time and count, in:Proceedings of the IEEE-INNS-ENNS International Joint Conference on
Neural Networks. IJCNN 2000. Neural Computing: New Challenges andPerspectives for the New Millennium, Vol. 3, IEEE, 2000, pp. 189–194.
[20] D.E. Rumelhart, G.E. Hinton, R.J. Williams, Learning representations byback-propagating errors, Nature 323 (1986) 533–536.
[21] Y. LeCun, Y. Bengio, G. Hinton, Deep learning, NATURE 521 (7553) (2015)
436–444.
[22] X. Yu, N.K. Loh, W.C. Miller, New acceleration technique for the backpropagation algorithm, in: IEEE International Conference on Neural
Networks, San Francisco, California, 28 March-1 April 1993, IEEE, 1993,
pp. 1157–1161.
[23] C.C. Yu, B.D. Liu, A backpropagation algorithm with adaptive learning rate
and momentum coefficient, in: 2002 International Joint Conference on
Neural Networks, IJCNN ’02, Honolulu, HI; United States; 12 May 2002,
IEEE, 2002, pp. 1218–1223.
[24] N. Huang, Z. Shen, S. Long, M. Wu, H. Shih, Q. Zheng, N. Yen, C. Tung,
H. Liu, The empirical mode decomposition and the Hilbert spectrum for
nonlinear and non-stationary time series analysis, Proc. Royal Soc. A-Math.
Phys. Eng. Sci. 454 (1971) (1998) 903–995.
[25] A.H. Phan, A. Cichocki, Tensor decompositions for feature extraction andclassification of high dimensional datasets, Nonlinear Theory Appl., IEICE
1 (1) (2010) 37–68.
[26] V.J. Lawhern, A.J. Solon, N.R. Waytowich, S.M. Gordon, C.P. Hung, B.J.
Lance, EEGNet: a compact convolutional neural network for EEG-based
brain-computer interfaces, J. Neural Eng. 15 (5) (2018).

The authors declare that they have no known competing financial interests or personal relationships that could have appeared
to influence the work reported in this paper.
Acknowledgments
J.S-C. contribution was partially supported by the Universityof Vic - Central University of Catalonia, Spain (No. R0947) andis also based upon COST Action CA18106, supported by COST
(European Cooperation in Science and Technology). The work of
A. C and A.H P. was partly supported by Ministry of Scienceand Higher Education grant No. 075-10-2021-068, and the joint
projects: Artificial Intelligence for Life between SKOLTECH andthe University of Sharjah.
References[1] J.R. Wolpaw, N. Birbaumer, D.J. McFarland, G. Pfurtscheller, T.M.
Vaughan, Brain–computer interfaces for communication and control, Clin.
Neurophysiol. 113 (6) (2002) 767–791.
[2] B.J. Lance, S.E. Kerick, A.J. Ries, K.S. Oie, K. McDowell, Brain–computerinterface technologies in the coming decades, Proc. IEEE 100 (Special
Centennial Issue) (2012) 1585–1599.
[3] J. V., F. Lotte, M. Tangermann, Brain-computer interfaces: Beyond medical
applications, Computer 45 (4) (2012) 26–34.
[4] G.J. N. Luis Fernando, Brain computer interfaces, a review, Sensors 12 (2)(2012) 1211–1279.
[5] F. Lotte, L. Bougrain, A. Cichocki, M. Clerc, M. Congedo, A. Rakotomamonjy,
F. Yger, A review of classification algorithms for EEG-based brain–computerinterfaces: a 10 year update, J. Neural Eng. 15 (3) (2018) 031005.
[6] K. Muller, C.W. Anderson, G.E. Birch, Linear and nonlinear methods forbrain-computer interfaces, IEEE Trans. Neural Syst. Rehabil. Eng. 11 (2)(2003) 165–169.
[7] S. Soman, et al., High performance EEG signal classification usingclassifiability and the twin SVM, Appl. Soft Comput. 30 (2015) 305–318.
[8] Z. Dalin, Y. Lina, Z. Xiang, W. Chen, S. Wang, B. Robert, B. Boualem,
Cascade and parallel convolutional recurrent neural networks on EEGbased intention recognition for brain computer interface, in: Proceedingsof the Thirty-Second AAAI Conference on Artificial Intelligence (AAAI18,
Humans and AI), February 2-7, 2018, New Orleans, la, USA., AAAI, 2018,
pp. 1703–1710.
[9] H. Ramoser, J. Muller-Gerking, G. Pfurtscheller, Optimal spatial filtering of
single trial EEG during imagined hand movement, IEEE Trans. Rehabil. Eng.
8 (4) (2000) 441–446.
11

K. Takahashi, Z. Sun, J. Solé-Casals et al.

Applied Soft Computing 122 (2022) 108811[27] P. Wang, A. Jiang, X. Liu, J. Shang, L. Zhang, LSTM-based EEG classificationin motor imagery tasks, IEEE Trans. Neural Syst. Rehabil. Eng. 26 (11)(2018) 2086–2095.

[28] Y.-H. Kwon, S.-B. Shin, S.-D. Kim, Electroencephalography based fusion two-dimensional (2D)-convolution neural networks (CNN) model for
emotion recognition system, Sensor 18 (5) (2018) 1383.
[29] F. Chollet, et al., Keras, 2015, https://keras.io.

12