CSP算法研究资料包.7z

# 分析mne中CSP的计算，为自己实现CSP算法做准备
# 协方差矩阵计算，这里有几个点需要注意下，这里使用的是 Optimizing Spatial Filters for Robust EEG Single-Trial Analysis 文章中的实现流程。
import pickle
import numpy as np
import scipy.io as io
from mne.decoding import CSP
import matplotlib.pyplot as plt

import mne

mne.set_log_level（False）

with open（r‘F:\BaiduNetdiskDownload\BCICompetition\BCICIV_2b_gdf\Train\CueLeftRight\trainData.pkl‘‘rb‘） as f:
    epochsTrain = pickle.load（f）
trainData = epochsTrain.get_data（[‘EEG:C3‘ ‘EEG:C4‘]）
trainLabels = epochsTrain.events[: -1]



csp = CSP（n_components=2 reg=None cov_est=‘epoch‘ norm_trace=False）
dataCSP = csp.fit_transform（trainData trainLabels）
print（dataCSP.shape）
plt.figure（0）
plt.plot（dataCSP）



EEG_Channels = trainData.shape[1]
EEG_Trials = trainData.shape[0]
classLabels = np.unique（trainLabels）
EEG_Classes = len（classLabels）
covMatrix = list（range（EEG_Classes））
trialCov = np.zeros（[EEG_Channels EEG_Channels EEG_Trials]）
for i in range（EEG_Trials）:
    E = trainData[i]
    EE = np.dot（E E.T）
    trialCov[: : i] = EE/np.trace（EE）

for i in range（EEG_Classes）:
    covMatrix[i] = np.mean（trialCov[: : trainLabels == classLabels[i]] 2）

print（‘a‘ covMatrix[0]）
print（covMatrix[1]）
covTotal = covMatrix[0] + covMatrix[1]

#########特征值降序
eigenvalues featurevectors = np.linalg.eig（covTotal）
egIndex = np.argsort（eigenvalues）[::-1]
eigenvalues = np.sort（eigenvalues）[::-1]
Ut = featurevectors[: egIndex]
# 矩阵白化
P = np.dot（np.diag（np.sqrt（1/eigenvalues）） Ut.T）
# 矩阵P作用求公共特征向量transformedCov1
transformedCov1 = P.dot（covMatrix[0]）.dot（P.T）
# 计算公共特征向量transformedCov1的特征向量和特征矩阵
eigenvalues featurevectors = np.linalg.eig（transformedCov1）
egIndex = np.argsort（eigenvalues）[::-1]
eigenvalues = np.sort（eigenvalues）[::-1]
U1 = featurevectors[: egIndex]
# 计算投影矩阵W
CSPMatrix = U1.T.dot（P）

# 计算特征矩阵
FilterPairs = 2
feature_train = np.zeros（[EEG_Trials 2*FilterPairs+1]）
feature_test = np.zeros（[EEG_Trials 2*FilterPairs+1]）
Filter = np.r_[CSPMatrix[:FilterPairs] CSPMatrix[-FilterPairs:]]

#从每一次实验中提取CSP特征
for t in range（EEG_Trials）:
    # 将数据投影到CSP空间
    projectedTrial_train = Filter.dot（trainData[t]）
    projectedTrial_test = Filter.dot（trainData[t]）
    # 投影后信号的方差作为特征
    variances_train = np.var（projectedTrial_train 1）
    variances_test = np.var（projectedTrial_test 1）
    for f in range（len（variances_train））:
        feature_train[t f] = np.log（variances_train[f]）

    for f in range（len（variances_test））:
        feature_test[t f] = np.log（variances_test[f]）



CSP_Train_faetures = feature_train[: 0:4]
CSP_Test_features = feature_test[: 0:4]



plt.figure（1）
plt.plot（CSP_Train_faetures[: :2]）

plt.show（）