简介：面向运动生物力学、康复工程与表面肌电信号处理学习者的Matlab实践包，聚焦由sEMG数据估算肌肉激活度这一常见分析任务，适合已具备基础Matlab操作、希望快速复现激活度曲线的本科生、研究生及科研人员。压缩包共5个文件，约135KB，含2个m脚本分别承担主流程与激活度计算函数，2个mat数据文件提供运动过程及最大自主收缩（MVC）条件下的sEMG样本，另有1个txt说明文档辅助理解数据组织与运行逻辑。已有3509人学习下载，说明该小体量案例在入门与课堂演示中较受关注。读者可直接基于现成数据跑通预处理、归一化与激活度求解流程，对照输出曲线检查算法效果，并以此为基础替换自己的实验数据、调整参数或扩展至不同肌肉与动作场景，减少从零搭建脚本和整理MVC基准数据的时间。

1. 从一段贴片电极信号到能画图的激活度曲线

做过肌电采集的开发者大概都经历过这个场景：贴片电极位置调了半天，采集卡驱动装好了，录了一段 sEMG 原始信号，打开一看全是密密麻麻的高频毛刺，根本看不出某块肌肉到底什么时候发力、发了多大力。肌肉激活度 Matlab 代码及数据这个资源，干的就是把这段原始电压信号翻译成 0 到 1 之间激活度曲线的事。它面向的是运动生物力学、康复评估、人机交互这群需要量化肌肉发力的人，也适合刚接触表面肌电信号处理、想照着现成流程跑一遍再动手改参数的学生和工程开发者。整个链路不长：读数据、滤波去噪、整流平滑、归一化映射，Matlab 脚本加一份配套数据就能复现，改造成自己的处理管线也不难。

2. 处理链路拆解：从原始 sEMG 到激活度曲线的四步流水线

这套代码的价值不在某个玄学算法，而在把 sEMG 处理里那几个必须做、又容易做错的环节串成了一条可复现的线。拿到资源先别急着点运行，把数据流捋清楚，后面调参数才知道每一步在动什么。

2.1 目录里都有什么

这类资源包的典型结构，一般会拆成数据目录和脚本目录两块：

真实包里的文件名可能略有出入，但功能划分基本就是这几块。我一般会先只跑main.m，确认能出图，再逐个文件进去看。这个顺序很重要，一上来就钻滤波器系数，很容易在数据读取的格式问题上卡住，血泪经验。

2.2 数据读取与采样率对齐

sEMG 的一切后续处理都建立在采样率正确的前提上。采样率写错，滤波截止频率全部失真，出来的激活度曲线时间轴也是错的。常见做法是在脚本开头把采样率显式声明成常量，而不是从文件里隐式读取：

% 采样率必须与采集设备设置一致，写错后面全废 fs = 2000; % 单位 Hz，sEMG 常用 1000~2000 raw = load('data/emg_raw.mat'); sig = raw.emg; % 单通道列向量，单位 V % 若数据是多通道矩阵，先取目标通道 if size(sig, 2) > size(sig, 1) sig = sig'; % 统一成列向量，避免后续维度翻车 end sig = sig - mean(sig); % 去直流偏置，电极极化常带来零点漂移 t = (0:length(sig)-1) / fs; % 时间轴，秒

fs是整条链路的基准参数，后面每个滤波器的归一化频率都由它算出来。sig - mean(sig)这一步别省，采集设备直流偏置会让信号整体抬高，不处理的话整流之后包络基线偏移，归一化出来的激活度永远到不了 0。多通道数据统一转成列向量，是为了让后面所有函数只处理一种形状，少一类维度报错。

2.3 带通滤波与工频陷波

原始 sEMG 的有效能量集中在 20~500 Hz，低于 20 Hz 是运动伪迹和电极漂移，高于 500 Hz 基本是噪声。50 Hz 工频干扰在国内环境几乎是必现的，所以陷波也得跟上：

% 20~500 Hz 带通，4 阶零相位滤波避免相位延迟 [b_bp, a_bp] = butter(4, [20 500]/(fs/2), 'bandpass'); sig_f = filtfilt(b_bp, a_bp, sig); % 50 Hz 工频陷波，Q 值控制陷波宽度 wo = 50/(fs/2); bw = wo/35; [b_notch, a_notch] = iirnotch(wo, bw); sig_f = filtfilt(b_notch, a_notch, sig_f);

butter的第一个参数 4 是阶数，阶数越高过渡带越陡，但相位非线性也更明显；用filtfilt而不是filter就是为了零相位，代价是首尾各有约三倍阶数的瞬态段需要丢弃。[20 500]/(fs/2)是归一化到奈奎斯特频率，fs 一旦变化这里自动跟着变。iirnotch的bw是带宽参数，Q 值调大陷波越窄、对邻近频段伤害越小，但工频漂移时就压不干净，这个权衡得根据实际采集环境试。

2.4 整流、平滑与激活度归一化

滤波后的信号是零均值振荡，要变成“发力大小”必须整流再求包络。低通截止频率是这里最敏感的参数：

% 全波整流 sig_rect = abs(sig_f); % 低通求包络，截止频率决定平滑程度 [b_lp, a_lp] = butter(4, 6/(fs/2), 'lowpass'); envelope = filtfilt(b_lp, a_lp, sig_rect); % 归一化到 0~1 激活度 ref = max(envelope); % 或用 MVC 试验的最大值 activation = envelope / ref; activation = min(max(activation, 0), 1); % 硬截断防越界

6 Hz 这个截止频率是常见取值，对应约 150 ms 的包络响应时间，能跟得上大多数肢体动作又不会抖成锯齿。用max(envelope)自归一化适合只看相对发力模式的场景；做跨被试、跨次对比时，标准做法是让受试者先做最大随意收缩（MVC），用那次的最大包络值当分母，否则不同人之间的激活度没有可比性。最后那行截断是后悔药，滤波瞬态偶尔会让归一化值轻微越过 1，画图时看着很难受。

3. 跑通与参数标定：一次完整复现的实操顺序

代码能跑和结果能信是两件事。这章按我自己复现别人处理管线的习惯，把顺序和验证点摆出来。

3.1 环境准备与首次运行

Matlab 版本不用太新，带 Signal Processing Toolbox 就能跑，butter、filtfilt、iirnotch都在这个工具箱里。第一次运行按这个顺序：

1. 把整个包解压到一个不含中文和空格的路径下，Matlab 对中文路径的处理时好时坏，这是老坑。
2. 在 Matlab 里cd到解压目录，把data/、脚本目录都addpath进去。
3. 打开main.m，先确认开头的fs和你的数据匹配。
4. 逐段运行（Ctrl+Enter），不要直接 F5 一把梭，中间每步都看一眼变量。

第三步的采样率核对是重点。资源包里的数据如果是别人采集的，采样率往往只在文档里写一句，脚本里可能是个魔数。我一般会在读数据后立刻length(sig)/fs算一下总时长，如果算出几百秒而实际只录了几秒，采样率多半写错了。

3.2 用一段已知信号验证链路

在信任这套代码之前，建议先造一段假信号把链路验一遍，这是我最推荐的排错手段：

% 构造已知信号：5 Hz 正弦 + 50 Hz 干扰 + 白噪声 t = (0:1/fs:4-1/fs)'; clean = 0.5 * sin(2*pi*5*t); % 模拟 5 Hz 发力节律 hum = 0.3 * sin(2*pi*50*t); % 工频干扰 noise = 0.05 * randn(size(t)); test_sig = clean + hum + noise; % 走一遍同一条处理链 test_out = process_emg(test_sig, fs); % 检查：激活度应在 5 Hz 节律上起伏，且被压在 0~1 assert(all(test_out >= 0 & test_out <= 1), '激活度越界');

process_emg是把前面滤波、整流、平滑、归一化封成一个函数的写法，方便复用。这段自检的关键在于：输入信号的节律是已知的，如果输出的包络峰值位置和 5 Hz 正弦对不上，说明滤波器的相位或参数出了问题。assert那行是防御性的，跑批处理时能第一时间把越界数据揪出来，不至于画完图才发现。

3.3 关键参数怎么调

跑通之后大概率要按自己的数据调参，这几个参数按敏感度从高到低排：

我的习惯是先固定带通，只调包络低通，因为后者对最终曲线的观感影响最直接。调到能看清发力起止点、又不丢动作细节为止，这个平衡点因任务而异：做静态等长收缩可以压到 3 Hz，做快速动态动作可能得放到 10 Hz。

4. 常见排查：激活度曲线不对时先看这几处

跑不通和跑出来不对是两类问题，后者更折磨人。下面这几条是我在 sEMG 处理里翻车最多的场景，按现象、原因、解决的顺序列出来。

4.1 现象：曲线整体贴在 0 附近不动

原因通常是归一化分母取错了。如果ref用了整段数据的max，而这段数据恰好没包含真正的最大发力时刻，激活度会被整体压低。另一种可能是保留下限截止频率时误把有效频段滤掉了。

解决：先画原始信号确认有没有明显发力段；把ref改成用 MVC 参考段的最大值；检查带通下限是否被误设成 50 Hz 以上。

4.2 现象：曲线上有规律的锯齿抖动

锯齿的周期如果是整秒或半秒，多半是工频残留；如果是无规律的密集抖动，是包络低通截止太高。

解决：先pwelch看一眼频谱，确认 50 Hz 或 100 Hz 处有没有尖峰；有就检查陷波器是否真的生效（filtfilt顺序写反过的人不少）；没有就把包络低通从 6 Hz 往 3 Hz 降。

4.3 现象：首尾段出现巨大尖峰

这是filtfilt的瞬态效应。零相位滤波本质上正向反向各滤一遍，首尾各有一段不稳定区，长度约为滤波器阶数的数倍。

解决：裁剪掉首尾各3*order个采样点，或者用filtfilt的边界延拓选项。我在批处理时会统一裁掉首尾 0.2 秒，简单粗暴但有效。

4.4 现象：不同被试的激活度没法比

症状是 A 的曲线整体比 B 高一截，但两人实际发力差不多。原因是每个被试的电极阻抗、皮下脂肪厚度不同，绝对幅值本来就不可比。

解决：放弃自归一化，改用 MVC 归一化；同时确认每个被试的电极位置记录一致，位置差一点幅值就差一截。

4.5 现象：Matlab 报维度不匹配

十有八九是信号从行向量和列向量之间来回倒。butter返回的系数、filtfilt的输入、abs的输出，只要有一处形状不一致就会在某个中间步骤炸掉。

解决：在读数据后立刻统一成列向量，所有函数只接受列向量；调试时在每步后面加size()打印，比看报错堆栈快得多。

5. 二次开发：把单机脚本改成可批处理的处理管线

资源包给的是单文件脚本，真要用起来基本都得改造成能批量跑多被试、多通道的管线。我的改造思路分三层：参数外置、函数封装、结果落盘。

第一层，把所有魔数抽成配置结构体。fs、带通上下限、包络截止、归一化方式全部塞进一个cfg，主脚本只读配置不写数字。这样换一套数据只改配置，代码一行不动：

cfg.fs = 2000; cfg.bandpass = [20 500]; cfg.notch_freq = 50; cfg.env_lp = 6; cfg.order = 4; cfg.norm_mode = 'mvc'; % 'mvc' | 'self' % 批处理：遍历目录下所有数据文件 files = dir('data/*.mat'); for i = 1:numel(files) S = load(fullfile('data', files(i).name)); act = process_emg(S.emg, cfg); % cfg 整体传进去 save(fullfile('out', files(i).name), 'act'); end

cfg.norm_mode用字符串做开关，process_emg内部用switch分支处理，比传一堆布尔标志清晰。批处理循环里每次结果单独存文件，好处是中途出错不会毁掉前面已算完的结果，坏处是文件会很多，得配个汇总脚本。

第二层，函数签名统一成(signal, cfg)两参数。我见过有人把 fs 单独传、滤波参数打包传、归一化参数又另传，调用点乱成一团。统一成 config 结构体后，加参数只需改 config 默认值，调用点零改动。

第三层，加一个结果校验环节。批处理最容易出的问题不是算错，是算出一堆 NaN 或者全零文件还没人发现：

% 每个文件算完做一次合法性检查 bad = isnan(act) | act < 0 | act > 1; if any(bad) warning('%s: %d 个点越界', files(i).name, sum(bad)); end if all(act == 0) warning('%s: 激活度全零，检查输入信号', files(i).name); end

warning而不是error，是因为批处理里个别文件有问题不该中断整批，但要留下痕迹。跑完扫一遍 warning 列表，比事后逐个打开文件看图省事得多。

改造完还有一点值得说：这套代码的处理链是单向的，信号进、激活度出，中间没有回传。如果要接入实时应用（实时反馈、外骨骼控制），需要把filtfilt换成支持流式的filter加状态缓存，代价是失去零相位特性，得有延迟补偿。这个转换不是调个参数的事，是架构层面的改动，动手前想清楚。

从那以后我每次改这类信号处理代码，都会先在合成信号上跑一遍再碰真实数据，因为合成的信号我知道正确答案长什么样，真实数据出问题时你根本分不清是算法错还是数据脏。希望这套拆解能帮到你，少踩几个我已经踩过的坑。

本文还有配套的精品资源，点击获取