CN105809124A - DWT- and Parametric t-SNE-based characteristic extracting method of motor imagery EEG(Electroencephalogram) signals - Google Patents

Abstract

The invention provides a DWT- and Parametric t-SNE-based characteristic extracting method of motor imagery EEG(Electroencephalogram) signals. First, effective time and frequency ranges of EEG characteristics are determined by using a Wigner-Ville distribution and power spectrum; the EEG signals in a specific time and frequency segment is subjected to three-layer discrete wavelet decomposition and statistical characteristic quantity including the average value, the energy average value, the mean square error and the like are calculated and are taken as the time frequency characteristic of the EEG signals; at the same time, a parameterization t-SNE algorithm is utilized for performing non-linear characteristic mapping on said wavelet coefficients and embedded coordinates corresponding to a low-dimensional space are taken as the non-linear characteristic; the two characteristics are standardized and a characteristic vector including both the time frequency information and the non-linear information of the EEG signals in the specific time frequency segment is obtained. According to the invention, EEG characteristics of compactness and completeness are obtained and a method for solving a problem of poor generalization performance of a traditional manifold learning algorithm in pattern classification application through fitting a multilayer forward propagation neural network to nonlinear mapping is proposed, so that accuracy of pattern classification of MI-EEG signals is improved further.

Description

The feature extracting method of Mental imagery EEG signals based on DWT and Parametric t-SNE

Technical field

The present invention is a kind of EEG Processing technology, is applied particularly to brain-computer interface (Brain-Computer Interface, BCI) Extraction to Mental imagery EEG signals feature in system, uses wavelet transform (Discrete Wavelet Transform, DWT) With parametrization t-distribution adjacent embedding (Parametric t-Distributed Stochastic Neighbor Embedding, Parametric at random T-SNE) method combined carries out feature extraction and fusion to Mental imagery EEG signals.

Background technology

Motion imagination EEG signals (Motor Imagery Electroencephalography, MI-EEG) is contained and is considered as the fortune of person Dynamic wish and abundant nervous physiology information, receive much concern in research fields such as brain cognition, brain application, and correctly understands and accurate The characteristic information really extracting MI-EEG is the key of its successful Application.

There is the features such as individual difference, non-linear, non-stationary and time-varying sensitivity for MI-EEG signal, owing to small echo becomes Change (Wavelet Transform, WT) on multiple dimensioned, signal to be decomposed, it is achieved localize while time and frequency zone, enter And the time-frequency transient state information effectively obtained under different scale different frequency in non-stationary signal, wavelet transformation is as a kind of classical Time-Frequency Analysis Method be widely applied in terms of MI-EEG signal time-frequency characteristics extraction.

But, human brain is a structure, the biosystem of function high complexity, is a typical nonlinear system.MI-EEG Signal not only comprises information and the rhythmicity feature of abundant Mental imagery consciousness, and has obvious nonlinear organization feature. Traditional based on lineary system theory time-feature extracting method frequently, original signal information can be inevitably resulted in and lose, also Contained nonlinear structure characteristics cannot be excavated from the eeg data of higher-dimension.Manifold learning (Manifold Learning, ML) recovery low dimensional manifold structure data can be gathered from higher-dimension, and obtain corresponding embedded coordinate, keep data message foot Dimensionality Reduction is realized, it is thus achieved that potential manifold structure under enough complete meanings.At present, with Local Liner Prediction (Locally Linear Embedding, LLE) and ML algorithm that Isometric Maps (Isometric Mapping, Isomap) is representative, Necessarily applied in EEG signal nonlinear characteristic maps.

But in actual applications, it has been found that existing epidemiology learning method has following deficiency: (1) is quick to data noise Sense, for having the most non-stationary, randomness and there is the EEG signal of noise and use manifold learning to carry out feature and carry Take or during Feature Dimension Reduction, be easily destroyed low-dimensional embedded structure, thus effect characteristics quality；(2) EEG signal has significantly Time-frequency distributions feature and non-linear behavior so that be difficult to obtain its substitutive characteristics merely with ML algorithm, it is impossible to ensure spy comprehensively Levy compactness and the completeness of vector, the most also can cause characteristic information redundancy and feature mismatch problem；(3) traditional stream Shape learning method can only carry out Data Dimensionality Reduction to given data set, it is impossible to produce one explicit from higher-dimension observation space to low The mapping relations of dimension embedded space, this is unfavorable for the generalization ability of manifold learning data outer for sample.Therefore cause at present with Manifold learning is that the Method of Data with Adding Windows of representative pattern classification accuracy rate in terms of MI-EEG signal is the highest.

Summary of the invention

In EEG feature extraction field, apply the deficiency of existence for existing manifold learning, the present invention propose a kind of based on Wavelet transform (Discrete Wavelet Transform, DWT) and parametrization t-distribution are adjacent at random embeds (Parametric t-SNE) The method of Mental imagery EEG feature extraction.Utilize the method to be possible not only to acquisition and can characterize MI-EEG signal comprehensively The brain electrical feature of time-frequency and nonlinear transformations, moreover it is possible to by multilayer feedforward neural network empty from high-dimensional data space to low dimensional feature Between nonlinear mapping relation record get off, and the feature extraction of the new samples being applied to outside training sample, this not only solves The extensive problem concerning study of traditional manifold learning, also improves the classification accuracy rate of EEG signals.

The technical solution used in the present invention is: first, and EEG signals is carried out pretreatment；Then, wavelet transform pair is used Often lead EEG signals and carry out multi-resolution decomposition；The statistics of the wavelet coefficient then calculating the sub-band relevant to Mental imagery task is special Levy the time-frequency characteristics as this EEG signals, use the parameterized t-SNE wavelet coefficient number to special frequency channel EEG signals According to dimensionality reduction, using the low-dimensional embedded coordinate of this EEG signals as its nonlinear characteristic；Finally, the characteristic vector being made up of the two It is input in support vector machine classifier classify, and with classification accuracy for joining according to in Parametric t-SNE algorithm Being optimized of number, finally extracts the brain electrical feature under optimized parameter.

The physiology of EEG feature extraction is according to being: when cerebral cortex region is activated, the metabolism in this region and blood flow Increase, cause the brain wave Alpha rhythm and pace of moving things (8～13Hz composition) and the reduction of the Beta rhythm and pace of moving things (14～30Hz composition) amplitude, be referred to as Event-related desynchronization；Under brain tranquillization or inert condition, the Alpha rhythm and pace of moving things and the Beta rhythm and pace of moving things show amplitude and substantially increase, It is referred to as event-related design.Therefore, people when the unilateral hands movement of the imagination, the brain of its offside respective primary sensorimotor cortex The electricity Alpha rhythm and pace of moving things and Beta rhythm and pace of moving things amplitude can reduce；And the EEG amplitude of the homonymy correspondence rhythm and pace of moving things can raise.This is follow-up EEG The determination of signal extraction time-frequency characteristics particular sub-band, it is provided that theoretical basis.This Event-related desynchronization of EEG signals (Event-Related Desynchronization, ERD) and event-related design (Event-Related Synchronization, ERS) Phenomenon becomes the fundamental basis analyzing and judging right-hand man's Mental imagery EEG signals.

Analyze based on above, the inventive method to implement process as follows:

Step 1, Signal Pretreatment.First, the feature produced MI-EEG signal from neuro physiology angle dissects, And when using Wigner-Ville distribution and average power spectrumanalysis to determine brain electrical feature effective, scope frequently, due to ERS/ERD Phenomenon shows more obvious in C3 and C4 leads EEG signals, is averaged so C3 and C4 two is only led signal by the present invention The analysis of power spectrum, its mean power P (j) is calculated by following formula:

$P\left(j\right)=\frac{1}{N}{\mathrm{\Σ}}_{i=1}^N{x}^2(i,j)---\left(1\right)$

Wherein, (i, j) represents that certain leads the jth data of MI-EEG signal i & lt experiment to x, and N is experiment number.Mean power Analysis of spectrum, respectively as in figure 2 it is shown, combine the time-frequency scope of the above-mentioned feature extraction determined, uses finite impulse response (FIR) (Finite Impulse Response, FIR) wave filter carries out 8-30Hz bandpass filtering to EEG signals, preliminary obtain EEG signals ERD and ERS physiological phenomenon shows significant target data segment；

Step 2, the EEG signals obtaining step 1 carries out the discrete wavelet transformation of three layers.L layer scattering small echo for signal f (t) Conversion can be expressed as:

$f\left(t\right)=A_L+{\mathrm{\Σ}}_{j=1}^L{D}_j---\left(2\right)$

In formula, L is Decomposition order, A_LRepresent low pass approximation component, D_jFor details coefficients under j yardstick, j=1,2 ..., L.Thus, Signal f (t) is divided into multiple sub-band.If the sample frequency of signal f (t) is f_s, then A_L, D_L, D_L-1…D₁Each component divides Not corresponding frequency band range is followed successively by: [0, f_s/2^L+1],[f_s/2^L+1,f_s/2^L],[f_s/2^L,f_s/2^L-1],…,[f_s/2²,f_s/2]。

Step 3, time-frequency characteristics extracts.Owing to wavelet coefficient have expressed signal in time domain and the Energy distribution of frequency domain, and Alpha joint Rule (8～13Hz) is closer to the frequency range (8～16Hz) of D3, and the Beta rhythm and pace of moving things (14～30Hz) and the frequency range of D2 (16～30Hz) are closer to, the most obvious two wave bands of ERD/ERS phenomenon of this MI-EEG signal just.Therefore, herein Feature extraction will be carried out based on D2 and D3 details coefficients.In order to portray brain electrical feature from energy point of view, herein will be in step 3 Average, average energy value and the mean square deviation of wavelet details component D2 and D3 obtained is as the time-frequency characteristics of EEG signals.

OrderRepresent C_iLead the jth layer details coefficients coefficient D of MI-EEG signal_j(i=3,4j=2,3), K=1,2 ..., n, its average is defined as:

${\stackrel{\‾}{\mathrm{\μ}}}_j^i=\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{m}_{j,k}^i---\left(3\right)$

Average energy value calculating formula is:

${\stackrel{\‾}{E}}_j^i=\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{\left(m_{j,k}^i\right)}^2---\left(4\right)$

Mean square deviation calculates according to following formula:

$S_j^i=\sqrt{\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{(m_{j,k}^i-{\stackrel{\‾}{\mathrm{\μ}}}_j^i)}^2}---\left(5\right)$

Comprehensive ERD and ERS phenomenon leads the performance on signal at C3, C4, more significantly may be used to make extracted feature have Indexing, during definition-feature F frequently₁∈R^6×1For

$F_1={\[\frac{{\stackrel{\‾}{\mathrm{\μ}}}_2^3-{\stackrel{\‾}{\mathrm{\μ}}}_2^4}{\left|\right|{\stackrel{\‾}{\mathrm{\μ}}}_2^3-{\stackrel{\‾}{\mathrm{\μ}}}_2^4\left|\right|},\frac{{\stackrel{\‾}{\mathrm{\μ}}}_3^3-{\stackrel{\‾}{\mathrm{\μ}}}_3^4}{\left|\right|{\stackrel{\‾}{\mathrm{\μ}}}_3^3-{\stackrel{\‾}{\mathrm{\μ}}}_3^4\left|\right|},\frac{{\stackrel{\‾}{E}}_2^3-{\stackrel{\‾}{E}}_2^4}{\left|\right|{\stackrel{\‾}{E}}_2^3-{\stackrel{\‾}{E}}_2^4\left|\right|},\frac{{\stackrel{\‾}{E}}_3^3-{\stackrel{\‾}{E}}_3^4}{\left|\right|{\stackrel{\‾}{E}}_3^3-{\stackrel{\‾}{E}}_3^4\left|\right|},\frac{S_2^3-S_2^4}{\left|\right|S_2^3-S_2^4\left|\right|},\frac{S_3^3-S_3^4}{\left|\right|S_3^3-S_3^4\left|\right|}\]}^T---\left(6\right)$

Wherein, ‖. ‖ represents and seeks 2-norm.

Step 4, Nonlinear feature extraction.Due to wavelet details component D2 and D3 have expressed ERD and ERS phenomenon time The Energy distribution situation of frequency domain, so using parameterized t-SNE algorithm to carry out Data Dimensionality Reduction for D2 and D3 after concatenation, Excavate and intrinsic geometry structure potential in data set for reconstruction.

For one group of high dimensional data X=[x after series connection₁,x₂,…,x_n]∈R^D×n, data dimension before wherein D is dimensionality reduction, N is number of samples, and the data after dimensionality reduction are Y=[y₁,y₂,…,y_n]∈R^d×n, d is the data dimension after dimensionality reduction.This algorithm By by the range information between Data In High-dimensional Spaces, the similarity being converted between data point, and attempt in lower dimensional space Reduce this similarity.

Step 4.1, for the high-dimensional data space before dimensionality reduction, the similarity between data is defined as conditional probability p_j|i, i.e. data x_jPhase For data x_iSimilarity be proportional to x_iCentered by the probability density of Gauss distribution:

$p_{j|i}=\frac{\exp (-||x_i-x_j|{|}^2/2{\mathrm{\σ}}_i^2)}{{\mathrm{\Σ}}_{k\≠i}\exp (-||x_i-x_k|{|}^2/2{\mathrm{\σ}}_i^2)}---\left(7\right)$

Wherein, σ_iRepresent with data x_iCentered by Gauss variance.

P is understood by the symmetry between data_j|i=p_i|j, so

$p_{ij}=\frac{p_{j|i}+p_{i|j}}{2}---\left(8\right)$

In order to overcome " congested problem " existed between data point in lower dimensional space, t-distribution is used to replace height in lower dimensional space This distribution, similarity definition is similar with higher dimensional space, is denoted as

$q_{ij}=\frac{{(1+|\left|f\right(x_i|W)-f(x_j|W\left)\right|{|}^2/\mathrm{\α})}^{-\frac{\mathrm{\α}+1}{2}}}{{\mathrm{\Σ}}_{k\≠i}{(1+|\left|f\right(x_k|W)-f(x_i|W\left)\right|{|}^2/\mathrm{\α})}^{-\frac{\mathrm{\α}+1}{2}}}---\left(9\right)$

Wherein, f represents the dimensionality reduction from higher dimensional space to lower dimensional space and maps, i.e. f (X)=Y, this mapping is by with W as weight Multilamellar of based on restricted Boltzmann machine (Restricted Boltzmann Machines, RBMs) propagated forward neutral net institute Definition, α represents the value of t-distribution degree of freedom.

Step 4.2, is denoted as P by the joint probability distribution between Data In High-dimensional Spaces point, connection between data point in lower dimensional space Close probability distribution to be denoted as the data that the core of Q, Parametric t-SNE is so that in lower dimensional space and maintain higher dimensional space as far as possible Similarity relation.Difference between joint probability distribution P and Q uses Kullback-Leibler divergence to measure, this algorithm Object function is:

$C=\underset{i}{\Σ}KL\left(P_i\right|\left|Q_i\right)=\underset{i}{\Σ}\underset{j}{\Σ}{p}_{ij}log\frac{p_{ij}}{q_{ij}}---\left(10\right)$

Following iterative of concrete optimizations process use

$Y^{\left(t\right)}=Y^{(t-1)}+\mathrm{\η}\frac{\mathrm{\δ}C}{\mathrm{\δ}y}+\mathrm{\α}\left(t\right)(Y^{(t-1)}-Y^{(t-2)})---\left(11\right)$

Wherein, t is iterations, and η is learning rate, and α (t) is momentum term during the t time iteration.

Just Y is can get after t iteration^(t)∈R^d×n。

Use high dimensional data X=[x in the method₁,x₂,…,x_n]∈R^D×nEmbedded coordinate Y in lower dimensional space is as correspondence Nonlinear characteristic F of data₂∈R^d×n。

Step 5, in order to obtain the characteristic vector comprising EEG signals Time-Frequency Information and nonlinear transformations, from the unified order of magnitude Brain electrical feature is carried out data classification, to the above brain electrical feature F extracted in this algorithm₁And F₂It is standardized and carries out serial spy Levy fusion, characteristic vector F of this algorithm can be obtained, be denoted as

$F=\[\frac{F_1}{\left|\right|F_1\left|\right|},\frac{F_2}{\left|\right|F_2\left|\right|}\]---\left(12\right)$

Wherein, ‖. ‖ represents second order norm.

Compared with prior art, the invention have the advantages that

(1) present invention compares tradition Mental imagery brain electrical feature method based on discrete small wave converting method, owing to comprising MI-EEG Time-frequency and nonlinear transformations, the characteristic vector of acquirement can more fully characterize the feature of EEG signals.Use support to The brain electrical feature using the present invention to extract is carried out in the experiment of pattern classification as grader by amount machine, compare traditional linear or Method of Nonlinear Dimensionality Reduction, this method obtains the classification accuracy rate of the highest 94.1%；Find simultaneously, compare other manifold learning arithmetic, Use Parametric t-SNE algorithm that higher-dimension brain electric information can be made to have, at three-dimensional feature visualization, become apparent from poly- Class characteristic distributions, for data visualization technique is applied to the preferred of brain electrical feature, and is provided for reference frame.

(2) the present invention is directed to tradition manifold learning and lack generalization ability not in MI-EEG signal mode identification neighborhood is applied Foot, i.e. cannot provide dominant mapping relations while high dimensional data is carried out dimensionality reduction, proposes to use based on restricted Boltzmann The multilamellar propagated forward neutral net of machine records the data nonlinear mapping relation from higher dimensional space to low dimensional feature space, for Data outside training sample use this multilayer neural network to complete Data Dimensionality Reduction, improve the real-time of algorithm application and extensive energy Power.

Accompanying drawing explanation

Fig. 1 .1 leads mean power situation of change for C3, C4 during imagination left hand motion.

Fig. 1 .2 leads mean power situation of change for C3, C4 during imagination right hand motion.

Fig. 2 is experiment sequential chart.

Fig. 3 is electrode position distribution.

Fig. 4 is the implementing procedure figure of this method.

Detailed description of the invention

In the present invention, specific experiment is to carry out under Windows 7 (32) system uses the simulated environment of Matlab2011a.

The MI-EEG data set that the present invention uses derives from " BCI Competition2003 " standard database, by Austria GRAZ University BCI research center provides.Therapy lasted 9s every time, concrete sequential is as shown in Figure 2.When t=0～2s, experimenter protects Hold resting state；When t=2s, display providing of short duration prompting sound while a continuously display tracking cross, experiment is opened Begin；The arrow to the left or to the right that tracking cross is randomly generated by one when t=3s replaces, and requires that experimenter is according to arrow simultaneously The motion of guided imagery right-hand man.Whole experiment is made up of 280 experiments, wherein is used for training for 140 times, is used for surveying for 140 times Examination, use AgCl as electrode, sample frequency is 128Hz, data by differential electrode from the 10-20 lead system of international standard Tri-passages of C3, CZ and C4 obtain, electrode is placed as shown in Figure 3.

In conjunction with concrete eeg data storehouse, as shown in Figure 4, it is as described below that inventive algorithm is embodied as step:

Step 1, Signal Pretreatment.First, in conjunction with mechanism of production and the signal own characteristic of EEG signal, based on Wigner-Ville Distribution carries out time frequency analysis to MI-EEG signal data collection, by the Energy distribution situation of signal in time frequency space, rationally determines Effective frequency window during feature extraction.The nervous physiology information contained by this MI-EEG data set knowable to time frequency analysis has substantially Time-frequency distributions feature, i.e. in 8～30Hz frequency ranges, signal energy is more concentrated, and this signal is mainly at 10H and 20Hz There is at frequency significant zonal distribution feature.Therefore use finite impulse response filter that EEG signals carries out 8-30Hz band logical Filtering；Filtered EEG signals is calculated to the mean power of MI-EEG data set according to formula 1.Result such as Fig. 1 .1-1.2 institute Show, as seen from the figure: right-hand man's MI-EEG signal that C3 and C4 leads in the mean power of 8～30Hz frequency ranges at corresponding brain Motor sensory area can show obvious ERD and ERS phenomenon, and shows the most prominent within 3.5～7s periods, and this is true The time range determining feature extraction provides foundation.

Step 2, the frequency band range obtaining step 1 is 8～30Hz, and time range is that the EEG signals in 3.5～7s periods is carried out The discrete wavelet transformation of three layers.After decomposition, the dimension of each layer wavelet coefficient A3, D1, D2, D3 is 76 respectively, 280, 144、76.

Step 3, time-frequency characteristics extracts.Owing to wavelet coefficient have expressed signal in time domain and the Energy distribution of frequency domain, and wavelet systems The most obvious two wave bands of ERD/ERS phenomenon of number frequency range MI-EEG signal just corresponding to D2 and D3.Therefore, Feature extraction will be carried out herein based on D2 and D3 details coefficients.Herein by the wavelet details component D2 for obtaining in step 3 The statistics such as average, average energy value and the mean square deviation with D3 as the time-frequency characteristics of EEG signals, computing formula be respectively formula 3, 4 and 5.

Step 4, Nonlinear feature extraction.In order to fully excavate the nonlinear transformations in MI-EEG signal, use parameterized T-SNE algorithm carries out Data Dimensionality Reduction for D2 and D3 after concatenation, excavates and intrinsic geometry structure potential in data set for reconstruction. For one group of high dimensional data X=[x after series connection₁,x₂,…,x_n]∈R^220×280, data dimension before wherein D is dimensionality reduction, n For number of samples, the data after dimensionality reduction are Y=[y₁,y₂,…,y_n]∈R^d×n, d is the data dimension after dimensionality reduction.

Step 4.1, in order to range information between data point in higher-dimension observation space is converted to the similarity between two data, and Keep its neighborhood relationships constant in low dimensional feature space, conditional probability p is defined as the similarity between High dimensional space data_j|i, This conditional probability can calculate according to formula 7 and formula 8.

In order to overcome " congested problem " existed between data point in lower dimensional space, t-distribution is used to replace height in lower dimensional space This distribution, similarity definition is similar with higher dimensional space, in like manner calculates data point conditional probability q in lower dimensional space according to formula 9_ij, Wherein the value of t-distribution degree of freedom α is 3.In this step, used in the present invention is multilamellar based on restricted Boltzmann machine Propagated forward neutral net, its network structure is 220-600-600-2500-d, i.e. the input layer of this neutral net, three hidden layers It is respectively 220,600,600,2500 and d, the excitation function between three first layers neutral net with the neuron number of output layer For sigmoid function, the excitation function between hidden layer and output layer is linear function, the learning rate of above two excitation function It is respectively 0.12 and 0.014.

Step 4.2, is denoted as P by the joint probability distribution between Data In High-dimensional Spaces point, connection between data point in lower dimensional space Close probability distribution to be denoted as the data that the core concept of Q, Parametric t-SNE algorithm is so that in lower dimensional space and remain high as far as possible The similarity relation of dimension space, i.e. data point are the most similar to joint probability distribution P in lower dimensional space and Q at higher dimensional space. In this algorithm, the difference between P and Q uses Kullback-Leibler divergence to measure.The object function of this algorithm such as formula 10 Shown in, iterative shown in use formula 11 is optimized.

In concrete iterative optimization procedure, being set to of parameter: total iterations is 70, dynamic in iterative 10 Quantifier, is set to 0.5 in front ten iteration, and remaining is set to the 0.7. fine setting stage in neutral net, uses conjugation ladder Degree algorithm carries out 35 iteration optimization with the method for batch processing.In the present invention, the low-dimensional Embedded dimensions of optimum is obtained by optimizing The value of d is 11.

Use high dimensional data X=[x in the method₁,x₂,…,x_n]∈R^220×280Embedded coordinate Y conduct in lower dimensional space Nonlinear characteristic F of corresponding data₂∈R^11×280。

Step 5, owing to the brain electrical feature using different characteristic extracting method to extract is the letter carried out MI-EEG from different perspectives Breath is understood, in order to obtain the characteristic vector comprising EEG signals Time-Frequency Information and nonlinear transformations and right from the unified order of magnitude Brain electrical feature carries out data classification, to the above brain electrical feature F extracted in this algorithm₁And F₂Use formula 12 is standardized, and can obtain Characteristic vector F of this algorithm.

For the Mental imagery EEG signals based on wavelet transform with the adjacent embedding at random of parametrization t-distribution that the present invention is proposed The effectiveness of feature extracting method is verified, utilize support vector machine (from LibSVM software kit,Http:// csie.ntu.edu.tw/～cjlin/) feature using the present invention to extract is carried out pattern classification experiment.Data base as using Correlational study contrast understands, and the feature extracting method that the present invention proposes, in the case of intrinsic dimensionality is less, obtains the highest by 94.10% Classification accuracy rate, the most as shown in table 1.

Table 1 and the contrast of correlational study result

Owing to characteristic vector F not only comprising feature F of the Time-Frequency Information of EEG signals₁, further comprises and there is linear separability feature F₂, Therefore the characteristic vector that the method obtains can more fully characterize the feature of EEG signals；In addition Parametric is being used During t-SNE algorithm carries out Data Dimensionality Reduction, use multilamellar propagated forward neutral net will be empty from data observation space to feature Between nonlinear mapping relation record get off, so for training sample outside new data can directly use neutral net number According to dimensionality reduction, not only increase the Generalization Capability of tradition manifold dimension reduction method, more improve EEG Processing real-time.

Claims (1)

1. the feature extracting method of Mental imagery EEG signals based on DWT and parametrization t-SNE, it is characterised in that: include Following steps,

Step 1, Signal Pretreatment；First, the feature produced MI-EEG signal from neuro physiology angle dissects, And when using Wigner-Ville distribution and average power spectrumanalysis to determine brain electrical feature effective, scope frequently, due to ERS/ERD Phenomenon shows more obvious in C3 and C4 leads EEG signals, is averaged so C3 and C4 two is only led signal by the present invention The analysis of power spectrum, its mean power P (j) is calculated by following formula:

$P\left(j\right)=\frac{1}{N}{\mathrm{\Σ}}_{i=1}^N{x}^2(i,j)---\left(1\right)$

Wherein, (i, j) represents that certain leads the jth data of MI-EEG signal i & lt experiment to x, and N is experiment number；Mean power Analysis of spectrum, respectively as in figure 2 it is shown, combine the time-frequency scope of the above-mentioned feature extraction determined, uses finite impulse response (FIR) (Finite Impulse Response, FIR) wave filter carries out 8-30Hz bandpass filtering to EEG signals, preliminary obtain EEG signals ERD and ERS physiological phenomenon shows significant target data segment；

Step 2, the EEG signals obtaining step 1 carries out the discrete wavelet transformation of three layers；L layer scattering small echo for signal f (t) Conversion can be expressed as:

$f\left(t\right)=A_L+{\mathrm{\Σ}}_{j=1}^L{D}_j---\left(2\right)$

In formula, L is Decomposition order, A_LRepresent low pass approximation component, D_jFor details coefficients under j yardstick, j=1,2 ..., L；Thus, Signal f (t) is divided into multiple sub-band；If the sample frequency of signal f (t) is f_s, then A_L, D_L, D_L-1…D₁Each component divides Not corresponding frequency band range is followed successively by: [0, f_s/2^L+1],[f_s/2^L+1,f_s/2^L],[f_s/2^L,f_s/2^L-1],…,[f_s/2²,f_s/2]；

Step 3, time-frequency characteristics extracts；Owing to wavelet coefficient have expressed signal in time domain and the Energy distribution of frequency domain, and Alpha joint Rule (8～13Hz) is closer to the frequency range (8～16Hz) of D3, and the Beta rhythm and pace of moving things (14～30Hz) and the frequency range of D2 (16～30Hz) are closer to, the most obvious two wave bands of ERD/ERS phenomenon of this MI-EEG signal just；Therefore, herein Feature extraction will be carried out based on D2 and D3 details coefficients；In order to portray brain electrical feature from energy point of view, herein will be in step 3 Average, average energy value and the mean square deviation of wavelet details component D2 and D3 obtained is as the time-frequency characteristics of EEG signals；

OrderRepresent C_iLead the jth layer details coefficients coefficient D of MI-EEG signal_j(i=3,4 j=2,3), K=1,2 ..., n, its average is defined as:

${\stackrel{\‾}{\mathrm{\μ}}}_j^i=\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{m}_{j,k}^i---\left(3\right)$

Average energy value calculating formula is:

${\stackrel{\‾}{E}}_j^i=\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{\left(m_{j,k}^i\right)}^2---\left(4\right)$

Mean square deviation calculates according to following formula:

$S_j^i=\sqrt{\frac{1}{n}\underset{k=1}{\overset{n}{\Σ}}{(m_{j,k}^i-{\stackrel{\‾}{\mathrm{\μ}}}_j^i)}^2}---\left(5\right)$

Comprehensive ERD and ERS phenomenon leads the performance on signal at C3, C4, more significantly may be used to make extracted feature have Indexing, during definition-feature F frequently₁∈R^6×1For

$F_1={\[\frac{{\stackrel{\‾}{\mathrm{\μ}}}_2^3-{\stackrel{\‾}{\mathrm{\μ}}}_2^4}{\left|\right|{\stackrel{\‾}{\mathrm{\μ}}}_2^3-{\stackrel{\‾}{\mathrm{\μ}}}_2^4\left|\right|},\frac{{\stackrel{\‾}{\mathrm{\μ}}}_3^3-{\stackrel{\‾}{\mathrm{\μ}}}_3^4}{\left|\right|{\stackrel{\‾}{\mathrm{\μ}}}_3^3-{\stackrel{\‾}{\mathrm{\μ}}}_3^4\left|\right|},\frac{{\stackrel{\‾}{E}}_2^3-{\stackrel{\‾}{E}}_2^4}{\left|\right|{\stackrel{\‾}{E}}_2^3-{\stackrel{\‾}{E}}_2^4\left|\right|},\frac{{\stackrel{\‾}{E}}_3^3-{\stackrel{\‾}{E}}_3^4}{\left|\right|{\stackrel{\‾}{E}}_3^3-{\stackrel{\‾}{E}}_3^4\left|\right|},\frac{S_2^3-S_2^4}{\left|\right|S_2^3-S_2^4\left|\right|},\frac{S_3^3-S_3^4}{\left|\right|S_3^3-S_3^4\left|\right|}\]}^T---\left(6\right)$

Wherein, ‖. ‖ represents and seeks 2-norm；

Step 4, Nonlinear feature extraction；Due to wavelet details component D2 and D3 have expressed ERD and ERS phenomenon time The Energy distribution situation of frequency domain, so using parameterized t-SNE algorithm to carry out Data Dimensionality Reduction for D2 and D3 after concatenation, Excavate and intrinsic geometry structure potential in data set for reconstruction；

For one group of high dimensional data X=[x after series connection₁,x₂,…,x_n]∈R^D×n, data dimension before wherein D is dimensionality reduction, N is number of samples, and the data after dimensionality reduction are Y=[y₁,y₂,…,y_n]∈R^d×n, d is the data dimension after dimensionality reduction；This algorithm By by the range information between Data In High-dimensional Spaces, the similarity being converted between data point, and attempt in lower dimensional space Reduce this similarity；

Step 4.1, for the high-dimensional data space before dimensionality reduction, the similarity between data is defined as conditional probability p_j|i, i.e. data x_jPhase For data x_iSimilarity be proportional to x_iCentered by the probability density of Gauss distribution:

$p_{j|i}=\frac{\exp (-||x_i-x_j|{|}^2/2{\mathrm{\σ}}_i^2)}{{\mathrm{\Σ}}_{k\≠i}\exp (-||x_i-x_k|{|}^2/2{\mathrm{\σ}}_i^2)}---\left(7\right)$

Wherein, σ_iRepresent with data x_iCentered by Gauss variance；

P is understood by the symmetry between data_j|i=p_i|j, so

$p_{ij}=\frac{p_{j|i}+p_{i|j}}{2}---\left(8\right)$

In order to overcome " congested problem " existed between data point in lower dimensional space, t-distribution is used to replace height in lower dimensional space This distribution, similarity definition is similar with higher dimensional space, is denoted as

$q_{ij}=\frac{{(1+|\left|f\right(x_i\left|W\right)-f\left(x_j\right|W\left)\right|{|}^2/\mathrm{\α})}^{-\frac{\mathrm{\α}+1}{2}}}{{\mathrm{\Σ}}_{k\≠i}{(1+|\left|f\right(x_k\left|W\right)-f\left(x_i\right|W\left)\right|{|}^2/\mathrm{\α})}^{-\frac{\mathrm{\α}+1}{2}}}---\left(9\right)$

Wherein, f represents the dimensionality reduction from higher dimensional space to lower dimensional space and maps, i.e. f (X)=Y, this mapping is by with W as weight Multilamellar of based on restricted Boltzmann machine (Restricted Boltzmann Machines, RBMs) propagated forward neutral net institute Definition, α represents the value of t-distribution degree of freedom；

Step 4.2, is denoted as P by the joint probability distribution between Data In High-dimensional Spaces point, connection between data point in lower dimensional space Close probability distribution to be denoted as the data that the core of Q, Parametric t-SNE is so that in lower dimensional space and maintain higher dimensional space as far as possible Similarity relation；Difference between joint probability distribution P and Q uses Kullback-Leibler divergence to measure, this algorithm Object function is:

$C=\underset{i}{\Σ}KL\left(P_i\right|\left|Q_i\right)=\underset{i}{\Σ}\underset{j}{\Σ}{p}_{ij}log\frac{p_{ij}}{q_{ij}}---\left(10\right)$

Following iterative of concrete optimizations process use

$Y^{\left(t\right)}=Y^{(t-1)}+\mathrm{\η}\frac{\mathrm{\δ}C}{\mathrm{\δ}y}+\mathrm{\α}\left(t\right)(Y^{(t-1)}-Y^{(t-2)})---\left(11\right)$

Wherein, t is iterations, and η is learning rate, and α (t) is momentum term during the t time iteration；

Just Y is can get after t iteration^(t)∈R^d×n；

Use high dimensional data X=[x in the method₁,x₂,…,x_n]∈R^D×nEmbedded coordinate Y in lower dimensional space is as correspondence Nonlinear characteristic F of data₂∈R^d×n；

Step 5, in order to obtain the characteristic vector comprising EEG signals Time-Frequency Information and nonlinear transformations, from the unified order of magnitude Brain electrical feature is carried out data classification, to the above brain electrical feature F extracted in this algorithm₁And F₂It is standardized and carries out serial spy Levy fusion, characteristic vector F of this algorithm can be obtained, be denoted as

$F=\[\frac{F_1}{\left|\right|F_1\left|\right|},\frac{F_2}{\left|\right|F_2\left|\right|}\]---\left(12\right)$

Wherein, ‖. ‖ represents second order norm.