认知科学实践：openBMI运动想象(二分类)

【Cheng】

已于 2024-02-01 22:44:24 修改

阅读量1.1k

点赞数 23

文章标签：分类数据挖掘人工智能机器学习算法

于 2024-02-01 17:26:27 首次发布

本文链接：https://blog.csdn.net/qq_62861601/article/details/135973103

版权

1、数据集介绍：

(一）数据集背景：

(二)数据集下载以及OpenBMI-master.zip工具箱介绍下载，matlab下载：

1、数据集介绍：

(一）数据集背景：

项目邀请54名健康受试者(年龄24-35岁;25名女性)参加了实验。其中38名被试者为初始脑机接口使用者，其他都有脑机接口实验的经验。所有被试都没有神经学、精神病学和其他相关疾病的病史。受试者坐在60(±5)cm的椅子上面对21英寸液晶显示器(刷新率: 60 Hz；分辨率: 1600 × 1200)。水平和垂直的近似视角分别为 37.7度和28.1度。在本项目只涉及到MI范例，故只说明MI范例的过程：对于每一个被试，每个实验的前3秒在监视器中心出现一个黑色的固定十字，之后出现左右箭头视觉提示，被试进行用合适的手抓握得想象任务4秒。每个任务结束后，屏幕空白时间为6s(±1.5 s)(见图1)。实验分为train和test两个阶段（注意这里并不是训练集和测试集）；每个阶段有100个实验，同时平衡左手和右手的图像任务。本实验分为两天，分别为session01和session02，每个被试者在两天分别做100个实验。

Motor-imagery 任务过程

(二)数据集下载以及OpenBMI-master.zip工具箱介绍下载，matlab下载：

数据集下载

数据集下载之前先了解一下数据集细节：

基于脑电图（EEG）的脑机接口（BCI）系统主要分为三大范式：运动想象（MI）、事件相关电位（ERP）和稳态视觉诱发电位（SSVEP）。我们记录了脑电图信号

(1)大量受试者(54名受试者）

(2)多个会议（两个不同天数[session1，session2]）

(3)使用多种范式（MI、ERP和SSVEP）

本次实验只需用到MI范式，下载的时候只需要下载MI范式的数据集即可。这个数据的文件名是sess01_subj1_EEG_MI.mat的格式数据集及工具箱下载链接：https://ftp.cngb.org/pub/gigadb/pub/10.5524/100001_101000/100542/

两个session里面分别是54个被试的实验数据，看要求你们需要多少个被试的数据，自己下载。然后每个被试都有三种范式的数据，只需要下载MI范式即可。

OpenBMI工具箱

OpenBMI工具箱提供了多种功能和工具，包括：

脑电图（EEG）信号处理：包括数据采集、预处理、特征提取等。
实时脑信号处理：支持实时监测和处理脑部信号，用于实时控制应用程序或设备。
机器学习算法：提供多种机器学习算法，用于分析和解释脑部信号，实现脑控制接口。
开放数据集：包括多个开放的脑电数据集，用于研究和开发BMI算法。

下载链接在上面。

matlab下载

参考博客：MATLAB下载+安装教程

参考

参考论文：如果想详细了解数据集可以仔细研读一下这篇论文，会有很大收获。

EEG dataset and OpenBMI toolbox for three BCI paradigms: an investigation into BCI illiteracy

2、数据预处理：

数据集结果以及详细介绍：

每一个文件，也就是单个session中的某一个subject，例如sess01_subj1_EEG_MI.mat，这个文件中其实包含了两个文件，通过Matlab读取的时候也能发现，一个是train，一个是test，但其实并不是真正意义上的训练集和测试集。如下图所示：

在Matlab打开之后可以看到数据的结构：

读取数据后，得到这么一个struct数据（不管是train还是test），因为脑电信号是随时间变化的图像，以下为部分理解

X：输入，可知维度

t：取样时间，这里1X100，表示有100个取样点

fs：取样频率

y_dec：y标签，1表示左，2是右。（这里也知道是一个二分类问题了）

y_logic：y标签的逻辑值

y_class：y标签的类

class：类别

chan：通道数，1X62表示是62个通道

time_interval：采样持续时间，这里是4000ms

这个时候，我们能知道特征和标签，之后会经常提到trial（这个概念对后续数据划分设计很重要），实际上对应就是y标签的数量，这里实际上是进行了100次间歇性实验，所以有100个结果，那么，一个train或者test对应就是100个trial

特征提取

看到x的维度非常奇怪，脑电信号有很多通道，自然是非常奇怪。我们要使用一些算法，把某些频率的信号进行增强，然后提取出其中重要的特征，然后作为处理后的输入。

常见的算法有CSP，FCSP等等，这些东西都已经集成在工具箱内，你需要看懂那些函数的输入和输出是什么，然后就可以对x进行滤波操作

代码（matlab）

clear;
sum_loss = double(0);
for i = 1:54
    if i < 10
        idx = strcat('0', string(i));
    else
        idx = string(i);
    end
    % 读取数据集
    str = sprintf('E:/pre_process_data/session1/s%d/sess01_subj%s_EEG_MI.mat', i, idx);
    [CNT_tr,CNT_te] = Load_MAT(str);
    % 训练集预处理
    CNT_tr = prep_selectChannels(CNT_tr,{'Name', {'FC5', 'FC3', 'FC1', 'FC2', 'FC4', 'FC6', 'C5', 'C3', 'C1', 'Cz', 'C2', 'C4', 'C6', 'CP5', 'CP3', 'CP1', 'CPz', 'CP2', 'CP4', 'CP6'}});
    CNT_tr = prep_filter(CNT_tr,{'frequency',[8,30]});划分频率8到30赫兹
    SMT_tr = prep_segmentation(CNT_tr,{'interval',[1000,3500]});%数据分段，10ms分一段，为250X20X100（trials）
%得到SMT下一步数据
%测试集预处理
    CNT_te = prep_selectChannels(CNT_te,{'Index',1:20});
    CNT_te = prep_filter(CNT_te,{'frequency',[8,30]});
    SMT_te = prep_segmentation(CNT_te,{'interval',[1000,3500]});
%CSP训练
   [CSP_tr,CSP_W] = func_csp(SMT_tr,{'nPatterns',[2]});%csp空间滤波，得到滤波后的训练集、滤波参数，滤波评分
   FT_tr=func_featureExtraction(CSP_tr,{'feature','logvar'});%对滤波后的训练集进行特征提取
   CF_PARAW=func_train(FT_tr,{'classifier','LDA'});%训练LDA模型，得到分类器参数
%CSP测试
    CSP_te = func_projection(SMT_te,CSP_W);
    FT_te = func_featureExtraction(CSP_te,{'feature','logvar'});
    cf_out = func_predict(FT_te,CF_PARAW);
    loss = eval_calLoss(FT_te.y_dec,cf_out);
%计算平均loss
    sum_loss = sum_loss + double(loss);
    loss
    FT_tr_x = FT_tr.x;
    FT_tr_y = FT_tr.y_dec;
    FT_te_x = FT_te.x;
    FT_te_y = FT_te.y_dec;
    strx_tr = sprintf('D:/rzkx_lab/data1/train_sess1/sess01_subj%s_FT_tr_x.mat', idx);
    stry_tr = sprintf('D:/rzkx_lab/data1/train_sess1/sess01_subj%s_FT_tr_y.mat', idx);
    strx_te = sprintf('D:/rzkx_lab/data1/test_sess1/sess01_subj%s_FT_te_x.mat', idx);
    stry_te = sprintf('D:/rzkx_lab/data1/test_sess1/sess01_subj%s_FT_te_y.mat', idx);
%保存处理好的结果
    save(strx_tr,"FT_tr_x");
    save(stry_tr,"FT_tr_y");
    save(strx_te,"FT_te_x");
    save(stry_te,"FT_te_y");
    sprintf('第%d被试的loss是', i)
end

avg_loss = sum_loss / 54

然后你会发现在Matlab左侧出现了一堆变量

Matlab会存储你的每一行代码运行结果，这点我觉得比较可以的，然后，FT变量就是提取出来的特征变量（就是feature的缩写），比如FT_te

进入这个变量，你会发现某些变量变化了，比如x，现在变成了4X100大小，也就是说，CSP滤波提取出了4个特征。

你需要做的就是把这里的x和y_dec作为数据的特征和标签导出来，也就是代码部分里的save。

3、跨被试和不跨被试介绍：

(一）划分方式

假如我们把所有的数据集都已经提取过，然后我们需要对数据进行拼接和划分，确定哪些是模型的训练集和测试集.

这里介绍一下概念，不跨被试的意思就是，仅仅对单个被试者的数据进行建模和预测，不把其它的被试者数据拼接。跨被试的意思就是，将你手上被试者的数据拼接在一起，进行建模和预测

第一种不跨被试

把每个subj的train和test拼接在一起（也就是200个trial），作为不跨被试的数据集，然后对数据集进行交叉验证，至于是几折，则依据情况确定（这个概念不懂的自己去查）

第二种不跨被试

把每个subj的train和test拼接在一起（200个trial），然后以session01作为训练集，session02作为测试集

跨被试

1个分类器，把每个subj的train和test拼接在一起（200个trial），然后9个subj拼接在一起作为训练集，剩余一个作为测试集（其实是留一验证法）

所有算法针对跨被试和不跨被试的一般处理流程：

1、跨被试：针对某一个具体的session，每次采用留一法选取其中一个subject的数据作为测试集，将其余所有subject信息合并作为训练集，放入所设计的分类器中，最后以准确率的平均值作为最终结果呈现。

2、不跨被试(不跨session)：针对每一个subject，放入分类算法中利用交差验证训练subject，以交叉验证得到的准确率的平均值作为最终呈现结果。

3、不跨被试(跨session)：使用session01的数据作为训练集，选取session02中对应subject的数据作为测试集，放入分类算法中进行预测得出准确率，以平均值作为最终呈现结果。

4、分类器的选择：

我这里用到的有SVM、逻辑回归、随机森林、朴素贝叶斯、KNN，你们也可以用其他更多的分类器模型。

这个是跨被试的参考代码。

import statistics

from scipy.io import loadmat
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score
from sklearn import svm
from sklearn.neighbors import KNeighborsClassifier
from sklearn.naive_bayes import GaussianNB
from sklearn.ensemble import RandomForestClassifier
import pyecharts.options as opts
from pyecharts.charts import Bar
import pandas as pd
import matplotlib.pyplot as plt

plt.rcParams["font.sans-serif"] = ["SimHei"]
plt.rcParams["axes.unicode_minus"] = False

# 将mat文件转化为dataframe，其中包含转置操作，返回值是dataframe数据集
def preprocess(data):
    for i in data.keys():
        if i == "FT_tr_x":
            data = data["FT_tr_x"]
        elif i == "FT_tr_y":
            data = data["FT_tr_y"]
        elif i == "FT_te_x":
            data = data["FT_te_x"]
        elif i == "FT_te_y":
            data = data["FT_te_y"]
        else:
            continue
    data = pd.DataFrame(data)
    data = data.T
    return data


# 拼接每个subj的train和test成为200个trial，返回值是200个trial
def one_mat(filename1, filename2):
    x1_train = loadmat(filename1)
    x1_train = preprocess(x1_train)
    x2_train = loadmat(filename2)
    x2_train = preprocess(x2_train)
    x_train = pd.concat([x1_train, x2_train])
    return x_train


# 每次选择一个i当作test，其余当作train数据，把中间数据放入一个list中返回
def left_one(all_data: list):
    all_train_test_data = []
    for i in range(len(all_data)):
        test_data = all_data[i]
        train_data = pd.DataFrame([])
        for j in range(len(all_data)):
            if i != j:
                train_data = pd.concat([train_data, all_data[j]])
        all_train_test_data.append((train_data, test_data))
    return all_train_test_data


# 通过svm计算 accuracy
def accuracy_model_svm(x_train, x_test, y_train, y_test):
    model = svm.SVC(C=10, kernel='linear', probability=True)
    model.fit(x_train, pd.DataFrame(y_train).values.ravel())
    y_pre = model.predict(x_test)
    accuracy = accuracy_score(y_test, y_pre)
    return accuracy


# 通过logic计算 accuracy
def accuracy_model_logic(x_train, x_test, y_train, y_test):
    log = LogisticRegression(max_iter=1000)
    log.fit(x_train, pd.DataFrame(y_train).values.ravel())
    y_pre = log.predict(x_test)
    accuracy = accuracy_score(y_test, y_pre)
    return accuracy


# 通过knn计算 accuracy
def accuracy_model_KNN(x_train, x_test, y_train, y_test, n_neibor):
    knn = KNeighborsClassifier(n_neighbors=n_neibor)
    knn.fit(x_train, pd.DataFrame(y_train).values.ravel())
    y_pre = knn.predict(x_test)
    accuracy = accuracy_score(y_test, y_pre)
    return accuracy


# 通过 朴素贝叶斯 计算 accuracy
def accuracy_model_GaussianNB(x_train, x_test, y_train, y_test):
    gaussian = GaussianNB()
    gaussian.fit(x_train, pd.DataFrame(y_train).values.ravel())
    y_pre = gaussian.predict(x_test)
    accuracy = accuracy_score(y_test, y_pre)
    return accuracy


# 通过 随机森林 计算 accuracy
def accuracy_model_RandomForest(x_train, x_test, y_train, y_test):
    rf = RandomForestClassifier()
    rf.fit(x_train, pd.DataFrame(y_train).values.ravel())
    y_pre = rf.predict(x_test)
    accuracy = accuracy_score(y_test, y_pre)
    return accuracy

sess = 1
total = 54
folderpath_train = f'D:\\rzkx_lab\\data\\train_sess{sess}'
folderpath_test = f'D:\\rzkx_lab\\data\\test_sess{sess}'
x_train = []
y_train = []

# 获取所有的 x_train 数据和 y_train 数据
for i in range(1, total + 1):
    i = f'0{i}' if i < 10 else str(i)
    tr_data_x = f'{folderpath_train}\sess0{sess}_subj{i}_FT_tr_x.mat'
    te_data_x = f'{folderpath_test}\sess0{sess}_subj{i}_FT_te_x.mat'
    tr_data_y = f'{folderpath_train}\sess0{sess}_subj{i}_FT_tr_y.mat'
    te_data_y = f'{folderpath_test}\sess0{sess}_subj{i}_FT_te_y.mat'
    x_train.append(one_mat(tr_data_x, te_data_x))
    y_train.append(one_mat(tr_data_y, te_data_y))

train_test_x_data = left_one(x_train)
train_test_y_data = left_one(y_train)

# 获取所有的 accuracy
accuracy_svm = []
accuracy_logic = []
accuracy_knn = []
accuracy_GaussianNB = []
accuracy_randomforest = []
for i in range(len(train_test_x_data)):
    train_x_data, test_x_data = train_test_x_data[i]
    train_y_data, test_y_data = train_test_y_data[i]
    accuracy_svm.append(accuracy_model_svm(train_x_data, test_x_data, train_y_data, test_y_data))
    accuracy_logic.append(accuracy_model_logic(train_x_data, test_x_data, train_y_data, test_y_data))
    accuracy_knn.append(accuracy_model_KNN(train_x_data, test_x_data, train_y_data, test_y_data, 6))
    accuracy_GaussianNB.append(accuracy_model_GaussianNB(train_x_data, test_x_data, train_y_data, test_y_data))
    accuracy_randomforest.append(accuracy_model_RandomForest(train_x_data, test_x_data, train_y_data, test_y_data))

print(f'总被试人数: {total}')
print('svm: ', statistics.mean(accuracy_svm), statistics.stdev(accuracy_svm))
print('逻辑回归: ', statistics.mean(accuracy_logic), statistics.stdev(accuracy_logic))
print('knn: ', statistics.mean(accuracy_knn), statistics.stdev(accuracy_knn))
print('朴素贝叶斯: ', statistics.mean(accuracy_GaussianNB), statistics.stdev(accuracy_GaussianNB))
print('随机森林: ', statistics.mean(accuracy_randomforest), statistics.stdev(accuracy_randomforest))


x_labels = []
for i in range(total):
    x_labels.append(f'sub-{i + 1}')

line = (
    Bar(init_opts=opts.InitOpts(width="1600px"))
    .add_xaxis(xaxis_data=x_labels)
    .add_yaxis(series_name='svm', y_axis=accuracy_svm)
    .add_yaxis(series_name='逻辑回归', y_axis=accuracy_logic)
    .add_yaxis(series_name='knn', y_axis=accuracy_knn)
    .add_yaxis(series_name='朴素贝叶斯', y_axis=accuracy_GaussianNB)
    .add_yaxis(series_name='随机森林', y_axis=accuracy_randomforest)
    .set_global_opts(title_opts=opts.TitleOpts(title=f'跨被试结果---session{sess}'), xaxis_opts=opts.AxisOpts(name='被试者编号'), yaxis_opts=opts.AxisOpts(name='accuracy(%)'))
    .set_series_opts(label_opts=opts.LabelOpts(is_show=False))
)

line.render()

5、模型结果评估指标：

准确率：准确率是指模型正确预测的数量所占全部预测的数量的比例，其公式为：

其中TP是正确做出positive的数量，TN是正确做出Negative的数量，FP是指错误地将Negative判定为Positive的数量，FN是指错误地将Positive判定为Negative的数量。

6、结果（准确度+标准差）做表格：

这个的话你们就按照自己的结果和审美来做。注意加减标准差就行。

举例：

完整代码：

链接：https://pan.baidu.com/s/17bDPJYiSWP15-Prx_FTcZA?pwd=2268
提取码：2268

【Cheng】

关注

23
点赞
踩
21

收藏

觉得还不错? 一键收藏
0
评论
认知科学实践：openBMI运动想象(二分类)

项目邀请54名健康受试者(年龄24-35岁;25名女性)参加了实验。其中38名被试者为初始脑机接口使用者，其他都有脑机接口实验的经验。所有被试都没有神经学、精神病学和其他相关疾病的病史。受试者坐在60(±5)cm的椅子上面对21英寸液晶显示器(刷新率: 60 Hz；分辨率: 1600 × 1200)。水平和垂直的近似视角分别为 37.7度和28.1度。
复制链接

扫一扫