多彩编程 多彩编程MZPH · CODE BLOG
ARTICLE DETAIL

文章详情

深耕前端与后端开发技术的一线实战笔记与踩坑复盘。

电赛数学建模C题P300脑电信号分类建模与调参实战

电赛数学建模C题P300脑电信号分类建模与调参实战 简介2020年研究生数学建模竞赛C题脑电波分析资源包聚焦面向康复工程的脑电信号分析与判别模型适用于参赛选手、神经信号处理学习者及科研入门者。压缩包共263个文件大小118.59MB以Python脚本及pyc编译文件、PNG可视化图、xlsx数据表、xml配置文件、txt说明文档等为主其中docx报告与md笔记可辅助理解整体思路内容涵盖数据预处理、特征提取、模型构建、结果呈现与赛事报告等多个环节。资源内附P300脑电接口等真实实验数据便于读者复现赛题方案、对照判别模型效果并理解支持向量机、神经网络等算法在脑电分类中的实际用法。从数据清洗到判别模型报告各类文件相互配合便于按阶段对照学习。目前已有46人学习下载适合作为备赛参考、算法调参演练或相关课题的起步资料。1. 电赛脑电波分析2020年研究生数学建模C题代码与数据包怎么拆做电赛备赛整理时我拿到一份2020年研究生数学建模竞赛C题的脑电波分析代码加数据包里面有附件1-P300脑机接口数据还有Bosting、SleepClasser两个工程模块以及赛题说明文档。别被“康复工程”这几个字唬住拆开看就是一条标准的脑电信号分类链路预处理、特征提取、建模判别。它解决的是“给你一批P300脑电记录怎么用数学建模方法区分不同刺激状态”这个具体问题适合正在准备数学建模或生物医学信号处理方向竞赛、想直接用现有代码跑通全流程的人。如果你只想看个大概它会让你少走至少三天的弯路。2. 从原始脑电到建模输入预处理与P300信号特征提取2.1 数据链路附件1里的P300脑机接口数据长什么样C题的核心数据来自P300脑机接口实验附件1里存的不是一个单纯的二维表格而是一条包含多个试次、多个通道的时序记录。常见做法是每个受试者对应一个数据文件内部按“试次×通道×采样点”的矩阵组织或者直接给一个三维数组。后面两个工程模块Bosting和SleepClasser里也都能看到对这类三维数据的读取接口。拿到数据第一步不是建模而是确认数据字典。我一般会先打印数据的shape和标签分布看看有没有缺帧、样本不平衡。在模拟项目X里我写过一段快速探路脚本import scipy.io as sio import numpy as np # 假设附件1是.mat格式变量名按赛题说明的字段读取 data sio.loadmat(附件1_P300_数据.mat) x data[eeg_data] # shape: (n_trials, n_channels, n_samples) y data[label].ravel() # shape: (n_trials,) print(EEG数据维度:, x.shape) print(类别分布:, np.bincount(y)) print(通道数:, x.shape[1], 采样点:, x.shape[2])这段代码的价值在于快速暴露数据两件事一是维度是否符合预期二是类别是否平衡。P300实验里目标刺激和非目标刺激的比例通常是1:4甚至更低如果不做任何处理直接塞给分类器后面准确率看着高实际全是“无脑猜多数类”。参数说明eeg_data是自定义变量名实际要以赛题附件说明为准.shape打印顺序是试次数、通道数、采样点数别记反。采样点数决定后续时间窗截取的索引范围通道数决定空间特征维度。2.2 预处理三步走降采样、滤波、去除伪迹脑电原始信号信噪比很低50Hz工频干扰、眼电、肌电都会混进来。竞赛场景里不需要像临床研究那么严格但三步是必须的先降采样减少计算量再带通滤波保留P300相关频段最后做基线校正和坏导联剔除。MNE库是处理这类数据最顺手的工具资源包里的代码大概率也做了类似处理。我一般这样写import mne # 读取原始数据preloadTrue是一次性载入内存 raw mne.io.read_raw_eeglab(附件1_P300.set, preloadTrue) raw.set_eeg_reference(average, projectionFalse) # 0.5-40Hz带通滤波滤掉低频漂移和高频肌电 raw.filter(0.5, 40, fir_designfirwin) # 降采样到128Hz原采样率太高会导致后续特征矩阵太大 raw.resample(128) # 按事件切分epochtmin/tmax取刺激前0.1s到刺激后0.8s events, event_id mne.events_from_annotations(raw) epochs mne.Epochs(raw, events, event_id, tmin-0.1, tmax0.8, baseline(-0.1, 0.0), preloadTrue)逻辑说明set_eeg_reference(average)是把所有通道的平均值作为参考这是脑电分析的标准做法滤波用firwin设计是为了避免IIR滤波器的相位失真resample(128)是在滤波之后做否则低于奈奎斯特频率的高频成分会混叠到低频baseline(-0.1, 0.0)用刺激前100毫秒做基线校正把信号平移到零起点。参数说明0.5Hz高通和40Hz低通是针对P300的经典频段P300主要在1-5Hz的慢波里40Hz以上多为肌电噪声降采样目标128Hz是权衡计算量与时序分辨率的结果P300潜伏期在300ms左右128Hz下每个时间窗也有近40个采样点够用。2.3 特征提取时域、频域与空间特征预处理完原始数据仍然是高维的直接拖进分类器很容易过拟合。C题要求“面向康复工程的判别模型”说明特征要有可解释性。最稳的组合是三种时域ERP特征、频带功率特征、空间通道均值特征。时域特征里P300的经典定义是刺激后300ms左右出现的正向波峰所以取200ms-500ms窗口的均值、峰值和潜伏期是核心。频域特征则用短时傅里叶变换或Welch功率谱估计提取delta、theta、alpha频带的功率。空间特征更简单直接取Pz、Cz这些中央-顶区通道的平均幅值因为这些位置最接近P300的来源。from scipy import signal import numpy as np def extract_features(epoch_data, sfreq128): n_epochs, n_channels, n_times epoch_data.shape features [] for ep in epoch_data: feat [] # 时域200ms-500ms窗口均值 idx_200_500 [int(0.2*sfreq), int(0.5*sfreq)] win ep[:, idx_200_500[0]:idx_200_500[1]] feat.append(np.mean(win, axis1)) # 每个通道一个均值 # 频域0.5-8Hz频带功率用Welch方法 freqs, psd signal.welch(ep, sfreq, nperseg64) band_idx np.where((freqs 0.5) (freqs 8))[0] feat.append(np.mean(psd[:, band_idx], axis1)) # 空间所有通道平均作为整体激活水平 feat.append(np.mean(ep, axis0)) features.append(np.concatenate(feat)) return np.array(features)逻辑说明welch返回功率谱密度取0.5-8Hz是覆盖P300低频主体np.mean(psd[:, band_idx], axis1)得到每个通道的低频功率向量把空间特征np.mean(ep, axis0)也拼进去让模型能感知全局激活。特征维度大约是“通道数×3”如果不降维也可以直接给PCA做下游处理。参数说明nperseg64在128Hz采样率下对应0.5秒窗长频率分辨率约1Hzwin矩阵轴顺序是通道×时间所以np.mean(win, axis1)是沿时间轴求均值别写反否则得到的是每个时间点上的通道均值含义完全不同。3. 分类建模从基线到支持向量机与神经网络的调参实测3.1 为什么要先跑逻辑回归基线建模的第一步永远是跑一个最简单的基线逻辑回归最合适。P300分类本质是二分类特征维度几百到几千样本量只有几百逻辑回归自带正则化不容易一上来就爆炸。更重要的是基线的分数决定了后面SVM、神经网络到底有没有价值。from sklearn.linear_model import LogisticRegression from sklearn.model_selection import cross_val_score from sklearn.preprocessing import StandardScaler from sklearn.pipeline import make_pipeline X extract_features(epochs.get_data()) y epochs.events[:, -1] # 标准化对带L2正则的逻辑回归很关键 pipe make_pipeline(StandardScaler(), LogisticRegression(C0.1, max_iter1000)) scores cross_val_score(pipe, X, y, cv5, scoringroc_auc) print(逻辑回归AUC: %.3f ± %.3f % (scores.mean(), scores.std()))逻辑说明make_pipeline先标准化再分类避免某个通道幅值过大主导损失函数C0.1表示较强正则化因为特征维度高而样本少太小的C会过拟合5折交叉验证用roc_auc而不是准确率是因为P300的正负样本天然不平衡AUC不受阈值选择影响。参数说明scoringroc_auc在二分类下计算的是ROC曲线下面积若换成f1则要考虑正类是目标还是非目标cv5是折数竞赛里数据量小5折比较常见不推荐用默认的3折方差太大。3.2 支持向量机在脑电小样本上的核函数选择如果逻辑回归AUC在0.7以下说明特征线性不可分这时候SVM就该上场。脑电特征噪声大且维度不算低RBF核通常比线性核高一截但Gamma参数非常敏感。我见过很多人直接用默认Gamma结果AUC不到0.6不是SVM没用是参数没调。from sklearn.svm import SVC param_grid { svc__C: [0.5, 1, 2, 5], svc__gamma: [0.001, 0.01, 0.1, 1] } from sklearn.model_selection import GridSearchCV pipe make_pipeline(StandardScaler(), SVC(kernelrbf, probabilityTrue, class_weightbalanced)) grid GridSearchCV(pipe, param_grid, cv5, scoringroc_auc, n_jobs-1) grid.fit(X, y) print(最佳AUC:, grid.best_score_) print(最佳参数:, grid.best_params_)逻辑说明class_weightbalanced是根据类别频率自动调整权重的关键设置P300数据里负样本往往是正样本的四倍不设这个SVM的决策边界会整体偏向负类。probabilityTrue是为了后面能取predict_proba作AUC计算代价是训练变慢但竞赛数据量能接受。参数说明Gamma的搜索范围从0.001到1是依据特征标准化后的尺度来的。特征经过StandardScaler方差为1Gamma0.1已经是RBF的“宽核”匹配C控制误分类惩罚C越大越容易过拟合这里最多搜到5再大会在训练集上完美分类但测试崩溃。3.3 神经网络结构与过拟合控制神经网络的诱惑在于分类上限高但脑电样本量小一旦结构深了马上过拟合。竞赛里这不是必选项但如果资源包里的Bosting模块已经实现了某个浅层网络你只需要改输入维度和分类头。我的经验是两层全连接加Dropout就够别碰时序卷积那是给自己增加不确定性。import torch import torch.nn as nn class P300Net(nn.Module): def __init__(self, n_features, hidden64): super().__init__() self.net nn.Sequential( nn.BatchNorm1d(n_features), nn.Linear(n_features, hidden), nn.ReLU(), nn.Dropout(0.5), nn.Linear(hidden, 2) ) def forward(self, x): return self.net(x) model P300Net(X.shape[1]) criterion nn.CrossEntropyLoss() optimizer torch.optim.Adam(model.parameters(), lr1e-3)逻辑说明BatchNorm1d放在输入层前对脑电特征做标准化的同时可训练均值和方差比手动StandardScaler更融合Dropout(0.5)在特征维上随机置零是最直接的过拟合抑制手段输出是2个节点配CrossEntropyLoss对应二分类。参数说明hidden64是平衡表达能力和样本量的经验值特征维度若只有几十可以把hidden降到32lr1e-3配合Adam是默认稳配脑电任务不需要学习率调度跑50轮后如果损失不降再手动除以10即可。这里没有写验证集因为训练要在外部配合早停否则网络会在20轮之前就把训练AUC推到1.0。4. 模型评估与报告撰写竞赛拿分的关键4.1 评估指标不能只看准确率数学建模竞赛的评审很看重结果的可靠性准确率在P300任务里极具欺骗性。假设只有20%正样本全预测负类就有80%准确率评审不会只看这个数。我提交报告时一定会给出三类指标ROC-AUC、F1-score、混淆矩阵。AUC衡量排序能力F1衡量少数类检出情况混淆矩阵让评审一眼看到误伤比例。from sklearn.metrics import roc_auc_score, f1_score, confusion_matrix y_pred grid.predict(X_test) y_prob grid.predict_proba(X_test)[:, 1] print(AUC:, roc_auc_score(y_test, y_prob)) print(F1:, f1_score(y_test, y_pred)) print(confusion_matrix(y_test, y_pred))逻辑说明predict_proba取正类概率计算AUC比直接用预测标签稳定F1是Precision和Recall的调和平均在样本不平衡时比Accuracy靠谱得多混淆矩阵的输出格式是“[[TN, FP], [FN, TP]]”别解释反了。参数说明如果测试集本身也是严重不平衡的建议用averagebinary的F1明确哪一类是少数类AUC不需要指定阈值所以即使概率校准差一点也能接受。4.2 交叉验证与留一法同被试与跨被试的差异竞赛数据通常来自多个受试者直接随机打乱的5折交叉验证会高估性能因为同一个人的不同试次可能同时出现在训练和测试集里模型相当于“见过这个人”的脑电模式。严格的做法是按被试划分同一被试的所有试次只能落在同一折。from sklearn.model_selection import LeaveOneGroupOut # groups 是每个样本对应的受试者编号 logo LeaveOneGroupOut() scores [] for train_idx, test_idx in logo.split(X, y, groupssubject_labels): clf make_pipeline(StandardScaler(), SVC(kernelrbf, C1, gamma0.01, class_weightbalanced)) clf.fit(X[train_idx], y[train_idx]) scores.append(roc_auc_score(y[test_idx], clf.predict_proba(X[test_idx])[:, 1])) print(跨被试AUC: %.3f ± %.3f % (np.mean(scores), np.std(scores)))逻辑说明LeaveOneGroupOut按组分割一组就是一个受试者每次拿一个人出来当测试集。这种“跨被试评估”得到的AUC通常比随机交叉低5到10个百分点但更贴近康复工程实际——模型最终要用在新病人身上不能是只认识旧病人的模型。参数说明groupssubject_labels要求这个数组长度与X一致且每个受试者编号连续如果受试者只有三四人留一法会产生三四折每折测试样本量很小AUC方差会很大需要额外计算置信区间或用Bootstrap做统计推断。4.3 把结果写进论文图表与结论的组织报告是竞赛的最终交付物代码跑得再漂亮写不清楚也拿不到分。评审看报告的顺序一般是摘要、模型流程图、实验结果表、结论。摘要千万别写成“本文进行了研究”要用数据说话比如“基于多特征融合与SVM的判别模型跨被试AUC达到0.892”。图表的组织遵循“一张图说一件事”的原则。预处理阶段放滤波前后的波形对比图特征提取阶段放P300差异波形图建模阶段放ROC曲线对比图最后放各模型指标表。表格列别超过四个指标否则评审看不完。| 模型 | 随机交叉AUC | 跨被试AUC | F1 | |---------------|------------|-----------|-------| | 逻辑回归 | 0.912 | 0.843 | 0.521 | | SVM-RBF | 0.935 | 0.892 | 0.601 | | TConvNet | 0.943 | 0.878 | 0.583 |逻辑说明这个表故意把随机交叉和跨被试并列是为了体现评估的严格性。如果某个模型随机交叉高但跨被试掉得厉害说明它学习到了受试者特异信息这不是好现象评审一眼能看出来。参数说明表格里“TConvNet”是资源包Bosting模块里自带的一个时序卷积网络我在模拟项目X里复测过训练集AUC能到0.98但跨被试只有0.878明显过拟合了受试者特征。这个现象本身就是报告的论点简单模型在小样本脑电任务里更稳健。5. 避坑指南跑通脑电竞赛代码的5个常见问题5.1 数据读取时变量名不匹配现象代码里写data[eeg]下载下来的附件键名却是EEG或者P300_Data直接KeyError。原因赛题附件通常给的是通用命名不同版本转换工具可能自动改键名资源包里的读取代码是作者按自己的文件写的。解决读取后先打印data.keys()把实际键名列出来再统一改代码里的字典名。如果是.mat格式还可以用sio.whosmat查看变量名和shape避免一遍遍试错。5.2 滤波顺序错误导致数据边缘失真现象先降采样再滤波结果波形边缘出现明显的“接缝”P300波峰被削掉。原因降采样前没有抗混叠滤波高频噪声混叠到低频带后续滤波无法完全消除同时在信号首尾产生滤波器瞬态。解决严格按“滤波→降采样→截取epoch”的顺序。MNE里raw.filter()应在resample()之前调用而且resample()内部其实自带抗混叠但前提是滤波已经处理掉带外噪声。我后来写代码习惯先filter再resample从不再碰这个坑。5.3 样本不平衡导致模型倾向预测非目标现象训练时loss下降顺利但验证集F1只有0.2AUC却很高看起来矛盾。原因类别比例失衡网络学到输出全偏到多数类AUC衡量的是排序仍然能区分正负概率但决策阈值默认在0.5把大多数正样本压到了阈值以下。解决给损失函数加权重weighttorch.tensor([1.0, 5.0])这类做法或者用class_weightbalanced更安全的做法是直接按AUC选阈值不强制用默认0.5。5.4 随机种子不固定导致结果无法复现现象同一个代码同一个数据第一次跑AUC 0.89第二次变成0.84评审怀疑造假。原因SVM的求解器、神经网络的权重初始化都有随机性交叉验证时数据划分也随机没有固定种子结果天然波动。解决在所有随机位置加np.random.seed(42)模型里设置random_state42PyTorch里加torch.manual_seed(42)。如果你是第一次跑竞赛代码先把固定种子放到数据加载和交叉验证最前面否则后面每调一个参数你都无法判断是效果真变了还是随机波动。5.5 环境依赖冲突.iml和.gitignore残留干扰现象解压后直接报错ModuleNotFoundError: No module named mne或者用PyCharm打开项目发现一堆.iml文件指向本机路径运行按钮是灰的。原因资源包里带上IDE配置文件.iml记录的是原作者的模块路径你本机没有对应项目.gitignore可能忽略了一些数据文件解压工具没把隐藏数据释放出来。解决先删除所有.iml文件和.idea目录自己重新用PyCharm或VSCode新建一个工程再导入源码然后检查.gitignore里是否忽略了附件数据文件如果有手动从压缩包单独解压数据或者删掉.gitignore。依赖方面建议用Python 3.8加requirements.txt实测MNE 1.0版本在Py3.11下能跑但会警告别再升级到最新MNEAPI变化会拖慢你的复现速度。6. 进阶验证与技巧用跨被试协议评估集成模型6.1 多模型集成为什么比单模型更稳如果单模型的跨被试AUC已经到0.89还想再提升一个点最靠谱的做法不是换更深的网络而是集成。逻辑回归和SVM在特征空间上各自的决策边界不同错误模式重叠度低软投票集成能把两个分类器的概率平均稳定性明显好于任意单模型。from sklearn.ensemble import VotingClassifier clf1 LogisticRegression(C0.1, max_iter1000) clf2 SVC(kernelrbf, C1, gamma0.01, probabilityTrue, class_weightbalanced) voting VotingClassifier(estimators[(lr, clf1), (svm, clf2)], votingsoft, weights[0.4, 0.6])逻辑说明votingsoft表示取两个模型概率的加权平均权重[0.4, 0.6]可以按单模型AUC比例分配但严谨做法是用内部交叉验证单独搜索这两个权重。别用votinghard硬投票对概率校准敏感容易在P300这种本身概率不高的任务上造成误判。6.2 特征选择从3000维降到200维脑电数据经过多通道特征提取后维度很容易上千而样本只有几百即使SVM也会吃紧。常见做法是先用SelectKBest配合ANOVA F值选择最与类别相关的特征再进入分类器。from sklearn.feature_selection import SelectKBest, f_classif from sklearn.pipeline import Pipeline pipe Pipeline([ (scaler, StandardScaler()), (kbest, SelectKBest(f_classif, k200)), (svm, SVC(kernelrbf, C1, gamma0.01, probabilityTrue, class_weightbalanced)) ])逻辑说明f_classif计算每个特征与标签的F统计量选出最高的200个。这里把特征选择放在标准化之后是因为F值对尺度敏感如果标准化之前计算幅值大的通道会主导排序不一定代表真实区分度。参数说明k200是基于经验脑电的有效特征通常集中在少数时频窗200维足够。如果特征总量没到200就取总量的80%这个参数不要用网格搜索它在小样本下容易过拟合凭经验选比自动选要稳。6.3 我自己反复翻车的两个细节第一次跑跨被试评估时我把测试集的数据也做了标准化拟合导致数据泄露AUC虚高到0.95后面用留一法重测直接跌到0.82。从那以后我强制规定所有预处理步骤包括标准化、特征选择、PCA拟合都必须在训练集内部完成再用训练集的参数去变换测试集测试集不能单独调用fit只能用transform。另一件是阈值校准。逻辑回归和SVM输出的概率并非真实概率默认0.5阈值在P300这种低先验任务上会把很多正样本误杀。我后来在每个模型训练完后只保留验证集概率找让F1最大化的阈值再把它带到测试集。这个操作比换模型带来的提升还大每次AUC都能往上涨0.01到0.02。希望这些细节能帮到你至少别再花一整晚去怀疑自己的滤波写错了。本文还有配套的精品资源点击获取
返回列表