ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

表面肌电信号归一化:MVC与min-max方法及MATLAB实现

表面肌电信号归一化:MVC与min-max方法及MATLAB实现 简介面向需要处理表面肌电信号sEMG的科研与工程人员这份MATLAB源码针对肌电信号幅值差异大、难以直接比较的问题实现了归一化处理并将归一化前后的波形以图形方式清晰显示方便用户观察信号整体变化。适合新手及有一定经验的开发人员快速掌握肌电信号的预处理与绘图方法也可作为生物医学信号处理课程设计或毕业设计的实用参考。资源包大小仅770B共1个文件核心为一个MATLAB脚本代码经过测试校正百分百可运行涵盖了从表面肌电信号读取、归一化参数计算、数值映射到最终对比图绘制的完整流程步骤清晰即使初学者也能轻松理解关键逻辑。用户还可在此基础上灵活调整归一化方式如极值归一化、均值方差标准化或添加滤波、特征提取等处理便于二次开发与算法比较资源由达摩老生出品质量有保障。目前已有1721人学习下载适用于康复评估、人机交互、运动分析等场景是一份小而精的肌电信号处理入门资源。1. 表面肌电归一化不是可选项为什么幅值必须折算到同一尺度表面肌电信号sEMG的原始幅值只要电极位置偏移 1 cm、贴合松紧不同或者皮下脂肪厚度有差异同一个受试者在相同力量水平下测到的电压就可能相差一倍以上。直接把两条 sEMG 曲线叠在一起谈“谁的激活更强”是没有意义的除非先折算到同一个参考尺度上。归一化处理解决的就是这个问题把实验或握力任务中采集到的原始表面肌电信号通过基准收缩或信号自身幅值区间转换为可比较的相对值。这个资源里的emg_jizhangli.m就是围绕“先消除噪声与漂移、再做 min-max 归一化、再用统一坐标显示”这条线完成处理的。适合做实验数据预处理、康复评估以及需要批处理肌电数据的 MATLAB 使用者不管现在是在校做课题还是在写仪器采集程序这套处理思路都能直接移植到自己的信号链路上。2. MVC 归一化与 min-max 归一化方法选择和信号预处理2.1 不同归一化方法的实际含义表面肌电归一化在生理信号处理里常用两种参考最大自主收缩maximum voluntary contraction, MVC和信号自身范围。MVC 归一化需要让受试者先做几次持续 3 秒左右的最大力量收缩取稳定段的 RMS 或平均整流值作为 100% 基准之后所有信号都表示成“MVC%”。它的优点是结果有明确的生理含义适合比较不同肌肉、不同受试者的激活程度。缺点是需要额外采集 MVC 数据并且受试者如果不会正确发力MVC 值容易偏低后续所有归一化值都会超过 100%图形直接失效。min-max 归一化则直接以当前这一段信号的幅值最小值作为 0、最大值作为 1计算简单适合同一段信号内部相对变化的分析。emg_jizhangli.m中比较稳妥的做法是先用带通滤波和陷波把噪声清掉再做 min-max 归一化这样可以不需要额外测量 MVC又能把握力过程中从静息到最大强度的变化范围压缩到 0-1 之间。要注意的是min-max 结果会被偶发的尖峰伪迹拉偏所以滤波和平滑步骤的质量直接决定归一化曲线是否可信。2.2 滤波参数怎么取带通 20-500 Hz 与 50 Hz 陷波表面肌电的有效能量主要分布在 20-500 Hz。低于 20 Hz 的部分主要是运动伪迹和基线漂移高于 500 Hz 的部分大多是测量电路的高频噪声。常见做法是先用零相位带通滤波器处理一次再用窄带陷波器把 50 Hz 工频干扰去掉。MATLAB 中filtfilt比filter更合适因为它是零相位滤波不会造成波形整体延迟多通道叠加显示时也不会出现通道之间的时间错位。fs 2000; % 采样率根据采集设备修改 f_low 20; % 高通截止频率单位 Hz f_high 500; % 低通截止频率单位 Hz [b, a] butter(4, [f_low/(fs/2), f_high/(fs/2)], bandpass); emg_filtered filtfilt(b, a, raw_emg); % 零相位带通滤波执行完这段代码后emg_filtered中仍可能残留 50 Hz 工频成分尤其是台式采集设备和未屏蔽线缆。再用iirnotch把 50 Hz 及附近 3 Hz 带宽内的成分衰减掉wo 50 / (fs/2); bw 3 / (fs/2); [b_notch, a_notch] iirnotch(wo, bw); emg_clean filtfilt(b_notch, a_notch, emg_filtered);参数上有几点经验四阶巴特沃斯在大多数 sEMG 场景下够用阶数太高会在通带边缘引入振铃如果采样率只有 1000 Hzf_high应降到 400 Hz因为 500 Hz 已经接近奈奎斯特频率的一半滤波器过渡带会变得不可控。陷波带宽 3 Hz 是折中值太窄收敛慢太宽会把 48-52 Hz 附近的真实肌电一起削掉。做完滤波后再执行emg_clean - mean(emg_clean)去掉直流偏置后续归一化才不会被静息基线抬高。2.3 归一化之前是否要整流和平滑严格来说min-max 归一化既可以直接作用于滤波后的双极性信号也可以作用于整流后的信号。直接作用于双极性信号的问题在于负半周的幅值会让最小值接近负的峰值归一化结果会呈现以 0.5 为中心的剧烈振荡看不出肌肉收缩的包络。因此更合理的流程是先把滤波后的信号取绝对值再做平滑得到反映肌肉激活程度的包络线。平滑窗口通常选 20-100 ms。窗口太短包络仍然抖动窗口太长会丢掉动作电位发放的细节。对于握力这类力量变化较慢的动作我一般用 50 ms 窗口对应代码是movmean(rect_signal, round(0.05*fs))。如果最终目的是做频谱分析就不要整流和包络如果目的是做图形显示和力量等级对比就一定要做整流和平滑。这里不要把平滑和低通滤波混为一谈movmean是滑动平均能用很直观的方式把包络拉出来。步骤常用参数作用注意事项带通滤波20-500 Hz4 阶去除运动伪迹与高频噪声采样率低时降低 f_high陷波滤波50 Hz3 Hz 带宽去除工频干扰确认当地电网频率为 50 Hz去直流mean 或 detrend消除基线偏移静息段较短时优先 detrend整流abs把负相位翻正用于包络显示平滑50 ms 移动平均或 RMS得到激活包络不要用于频谱分析3. 实现 emg_jizhangli.m从原始 EMG 到归一化曲线3.1 数据读取与变量命名资源中的emg_jizhangli.m是一个可直接运行的 MATLAB 脚本它的核心数据流可以抽象为加载、滤波、归一化、绘图四个环节。假设数据文件和脚本放在同一目录常见做法是通过load命令把信号读入工作区raw_data load(grip_emg.mat); raw_emg raw_data.emg; % 假设数据文件中变量名为 emg fs 2000;load返回一个结构体访问字段名要和grip_emg.mat中保存的变量名一致。如果数据是 CSV 格式改用readmatrix(grip_emg.csv)第一列是时间第二列是肌电。读者如果拿到的是没带数据文件的纯脚本也可以先用一段模拟信号验证归一化流程t 0:1/fs:5; % 生成 5 秒时间轴 raw_emg 0.45*sin(2*pi*3*t) 0.2*randn(size(t));这里的目的是把关注点放在归一化本身而不是卡在数据格式上。3.2 滤波和归一化的完整代码下面是一段与脚本思路一致的参考实现把上一章的预处理步骤串联起来。实际采集时放大器的直流偏置可能比较高运动伪迹也会带来缓慢漂移所以先取静息段的前 0.5 秒平均值并减去保证 min 值反映的是真实基底噪声raw_emg raw_emg - mean(raw_emg(1:round(0.5*fs))); % 去直流 [b, a] butter(4, [20/(fs/2), 500/(fs/2)], bandpass); emg_f filtfilt(b, a, raw_emg); [bn, an] iirnotch(50/(fs/2), 3/(fs/2)); emg_f filtfilt(bn, an, emg_f); rect_emg abs(emg_f); smooth_emg movmean(rect_emg, round(0.05*fs)); min_val min(smooth_emg); max_val max(smooth_emg); normalized_emg (smooth_emg - min_val) / (max_val - min_val eps);这段代码的关键是分母中的eps。当一段信号完全静息、包络几乎为零时max_val - min_val可能接近 0加入eps可以避免出现Inf或NaN。处理完成后normalized_emg被限制在 0 到 1 之间0 表示当前段的最低激活水平1 表示本次任务中的最大激活水平。3.3 为什么不能用原始信号直接除最大值有些教程会把归一化写成norm raw_emg / max(abs(raw_emg))。这种写法在双极性信号上会出现两个问题一是正峰是 1负峰是 -1图形看上去是幅度调制而不是激活程度二是基线漂移会让最大值偏高整条曲线被压低力量差异看不出来。即使先取abs再除以最大值也只是做了等比例缩放并没有把基线放到 0静息阶段的噪声仍然会在图上占据较大视觉比例。因此在emg_jizhangli.m中先整流、平滑再通过 min 和 max 做线性映射才能得到一条从 0 附近开始、峰值接近 1 的激活曲线。不同受试者之间的差异此时会更直观地反映在曲线的上升斜率、平台宽度和回落陡度上而不是绝对幅值。平滑窗口适用场景说明20 ms快速挥动、力量突变包络抖动大峰值偏高50 ms握力、等长收缩包络稳定适合演示100 ms慢速力量变化趋势清楚但细节丢失4. 绘图显示归一化结果坐标轴、对比与伪影识别4.1 基础绘图与坐标轴设置归一化完成后的图形显示不是简单一句plot就结束。表面肌电归一化曲线的横轴是时间纵轴刻度在 0 到 1 之间纵轴留白过大会让曲线显得很平留白过窄又会让曲线顶部被截断。我会用下面这种方式一次画出两行子图上一行是带通和陷波后的原始肌电下一行是同一时段的归一化包络方便核对信号对齐情况t (0:length(normalized_emg)-1) / fs; figure(Color, w); subplot(2,1,1); plot(t, emg_f, Color, [0.30, 0.30, 0.30]); ylabel(原始肌电 (mV)); xlim([0, t(end)]); grid on; subplot(2,1,2); plot(t, normalized_emg, LineWidth, 1.2, Color, [0.85, 0.33, 0.10]); ylim([-0.05, 1.05]); ylabel(归一化幅值); xlabel(时间 (s)); grid on;ylim([-0.05, 1.05])留出 5% 边距既能完整看到 0 和 1又不会让曲线贴在边框上。上一行的原始肌电不要用强饱和色波形密集时会变成一片黑色[0.30, 0.30, 0.30]这样的深灰已经足够。4.2 多通道和多段收缩的对比实际握力实验经常同时采集多个通道比如指浅屈肌、桡侧腕屈肌、肱桡肌。把每个通道的normalized_emg放到同一张图里配合图例可以快速看出发力时各肌肉的启动顺序channels {ch1_normalized, ch2_normalized, ch3_normalized}; labels {指浅屈肌, 桡侧腕屈肌, 肱桡肌}; colors lines(3); figure(Color, w); hold on; for k 1:length(channels) plot(t, channels{k}, LineWidth, 1, Color, colors(k, :)); end hold off; legend(labels, Location, northeast); ylim([-0.05, 1.05]); ylabel(归一化幅值); xlabel(时间 (s));这里有一个容易忽视的细节不同通道的 min 和 max 是各自独立计算的因此两个通道即使原始幅度差 3 倍归一化后也都能到 1不能据此判断哪块肌肉更“强”。要比较肌肉之间的强度差异需要所有通道共用一个参考值比如后面提到的 MVC。如果只是观察启动顺序各通道各自归一化是合理的但描述“谁发力更大”时就会产生误导。所以emg_jizhangli.m的纵轴只标“归一化幅值”避免把相对值直接当成绝对值解释。4.3 从图形上判断归一化是否成功图形是质量检查的第一道关卡。如果归一化后的曲线在静息段不断出现窄尖峰说明陷波带宽太窄或运动伪迹没有去干净需要回头检查带通滤波器的f_low。如果曲线在最大发力处出现平顶也就是有一段持续等于 1说明该段信号已经发生过载削波min-max 归一化的最大值不是真实肌电峰值而是放大器饱和电压必须回退到原始采集数据重新检查。反过来如果发现静息段整段都是 0不是“无激活”的理想表现通常是平滑窗口太短静息噪声没有被纳入包络min_val取到了接近零的极小值。这时把平滑窗口加大到 80 ms曲线通常会变得更连续。这些判断做完之后再保存图形才有意义。绘图参数推荐值作用ylim[-0.05, 1.05]给 0 和 1 留出边距LineWidth1 至 1.5防止多通道重叠时看不清网格on辅助比较时间点和峰值位置图例位置northeast标注通道或实验状态5. 验证 MVC 区间自动定位最大自主收缩的小技巧min-max 归一化虽然能把曲线压缩到 0-1但当你想把结果解释为“接近最大自主收缩的百分比”时只取整段信号的最大值并不可靠。如果最大收缩只是短暂尖峰中等力量收缩会被显示成接近 100%如果受试者没有真正发力最大值又没有达到生理上限归一化曲线整体会被抬高。更实用的做法是识别出收缩的稳定平台段用平台段的平均值作为 100% 基准。先从包络信号计算移动 RMS找到最高平台的中心位置rms_win sqrt(movmean(emg_clean.^2, fs)); % 1 秒窗口 RMS [~, idx] max(rms_win); % 峰值位置 half round(fs / 2); % 前后各取 0.5 秒 mvc_idx max(1, idx-half) : min(length(rms_win), idxhalf); mvc_level mean(rms_win(mvc_idx)); % 作为 100% 基准得到mvc_level后把所有时刻的短窗 RMS 转换成百分比rms_short sqrt(movmean(emg_clean.^2, round(0.05*fs))); pct_mvc rms_short / mvc_level * 100; plot(t, pct_mvc); ylabel(MVC %); ylim([0, 120]);这里的参数选择和普通 min-max 有本质区别MVC 用的是峰值附近的邻域平均值而不是瞬时最大值。瞬时最大值叠加了噪声和电极运动干扰峰值中心 0.5 秒的窗口才能代表“能维持住”的最大能力。两个移动窗口长度要分开检测 MVC 用 1 秒窗口能看到稳定力量平台绘制百分比曲线用 50 ms 窗口保留动作细节。实际使用中容易遇到一个坑当数据里有多次重复握力时max(rms_win)会落在第一次或最后一次发力处而第一次发力往往偏小最后一次可能因疲劳而撑不住。这时只搜索整个时长的中间 80% 区域start_idx round(0.1 * length(rms_win)); end_idx round(0.9 * length(rms_win)); [~, idx] max(rms_win(start_idx:end_idx)); idx idx start_idx - 1;这样识别出的 MVC 区间更稳定再把它单独存成mvc_segment.mat下一次对同一批数据进行归一化时直接调用避免手工截取带来的批次差异。本文还有配套的精品资源点击获取
返回列表