MFCC原理与Matlab实现:语音特征提取全解析
发布时间:2026/9/9 20:13:17 作者:尧图编辑部 阅读量:1,286

简介一份面向语音识别、情感分析、语音合成等场景的MATLAB代码完整实现了梅尔倒谱系数(MFCC)的标准提取流程涵盖预加重、分帧、加窗、快速傅里叶变换、梅尔滤波器组、对数运算、离散余弦变换以及Delta与Delta-Delta动态特征计算。代码注释细致、结构清晰适合语音处理初学者系统学习也可作为研究者和开发者的高效特征提取工具。压缩包共23个文件以MATLAB脚本(.m)为主包含主函数、对比脚本与工具函数另配有音频样本(.wav)、MFCC特征文件(.mfc)、PDF/TXT说明文档、配置文件及演示截图包体仅291KB便于快速查阅和二次开发。已有3014人学习附带的HTK对比测试脚本可帮助验证MFCC提取结果的一致性用户还可根据实际需求灵活调整帧长、帧移、滤波器数量等参数为后续语音建模提供可靠特征。 做语音识别或者说话人识别绕不开的一个东西就是梅尔倒谱系数也就是常说的MFCC。刚入行那会儿我看论文里一堆公式差点被劝退后来自己用Matlab老老实实跑了一遍流程才发现它的设计思路其实特别直白就是模拟人耳在不同频率下听觉敏感度不一样这件事。这篇文章不整虚的直接从原理到Matlab代码一步步拆开讲清楚每一行在干嘛、为什么这么干。无论你是正在做毕设、搞声学项目还是只是好奇语音特征怎么来的照着抄就能用。1. 整体思路拆解MFCC到底在解决什么问题1.1 为什么选择倒谱域特征语音信号里头真正对识别有用的信息并不是原始波形本身而是声道形状变化带来的共振峰特性。人类发声时声带振动产生声源激励经过口鼻腔的滤波作用才形成语音。这个滤波作用就是声道传递函数MFCC要做的就是把声源激励的影响和声道滤波的影响分离开。怎么分离这里用的是信号与系统里的老方法倒谱分析。把语音信号做傅里叶变换得到频谱对频谱取对数再做一次傅里叶变换就进入了倒谱域。在倒谱域里声道滤波的慢变成分和声源激励的快变成分会被拉开到不同的倒谱区间这样截取合适的倒谱系数就能只保留声道参数丢掉基频和激励细节。对我们做识别的人来说剩下的就是一组稳定又紧凑的特征向量很适合送进分类器。1.2 梅尔尺度让机器学会“听”单纯的倒谱系数还不够因为人耳在低频段分辨能力强高频段分辨能力差这种特性用线性频率无法表达。梅尔尺度就是一张把线性频率映射成符合人耳感知的刻度表公式很经典mel(f) 2595 * log10(1 f / 700)可以看到这是一个对数式变换把1000Hz以下的频段拉伸把高频压缩。换句话说我们的机器如果直接用线性谱做识别等于浪费了大量在高频段本就不敏感的感知信息还会让低频段的关键特征被遮掩。在MFCC流程里加梅尔滤波器组就是用一组在梅尔轴上等间距、在频率轴上随着中心频率升高而加宽的三角滤波器把频谱进行加权聚合模拟人耳的临界频带效应。选Matlab来实现这套流程优势很明显。信号处理工具箱自带了设计滤波器、分帧等函数而且数组操作效率高对音频读入、矩阵化计算支持很全面。我自己实测同样一套代码在Python要装librosa这套第三方库虽然也行但Matlab对矩阵逻辑的暴露更直观对把原理搞清楚非常有帮助。2. MFCC提取全流程拆解2.1 从波形到特征的七个关键步骤完整的一套MFCC提取流程业界基本已经标准化了每一步都有明确的物理含义和参数设置经验值。这里我按照常用的HTK和Kaldi开源工具的思路来讲个人用的话可以微调。第一步是预加重。语音的声门激励有约-6dB/倍频程的滚降从频谱上看高频部分能量衰减得厉害。预加重就是用一个一阶高通滤波器y(n) x(n) - αx(n-1)其中α一般取0.95到0.97之间。它不是为了改变声音听感而是让频谱变得平坦保证后续做FFT时高频细节不会被数值精度吃掉。第二步是分帧。语音信号是典型的非平稳信号但在极短的时间段内可以看作平稳的这个极短时间一般是20到40毫秒。我们把信号切成一段段短帧帧长常用25ms帧移10ms。为什么要重叠因为一帧边缘加窗之后能量会被削弱相邻帧重叠可以补偿信息丢失保证特征序列平滑过渡。第三步是加窗。对每一帧信号乘以窗函数最常用的是汉明窗w(n) 0.54 - 0.46 * cos(2πn / (N-1))。窗函数的目的是抑制帧边缘的不连续跳变避免FFT之后出现频谱泄漏。很多新手会跳过这步结果发现提取出来的特征噪声很大原因就是没加窗或者窗类型不对频谱被旁瓣污染了。第四步是FFT。对加窗后的信号做快速傅里叶变换取幅度谱。注意是幅度谱不是功率谱有些实现会用幅度的平方差别不大但幅度谱的数值范围更稳一点。第五步是梅尔滤波器组滤波。设计一组三角滤波器在梅尔轴上均匀分布每个滤波器和频谱做加权求和得到M个频带能量。这个M是滤波器数量常用26个或40个。第六步是取对数。人耳对声音响度的感知是指数型的取log之后能量特征更贴近人的主观感受同时对输入信号的幅度缩放不敏感也就是具备了一定鲁棒性。第七步是DCT离散余弦变换。对M个对数能量做DCT得到L个倒谱系数L一般取13或者20。DCT的作用是把滤波器组造成的相关性去除让特征各维度尽量独立后面接GMM或HMM做概率建模时效果好很多。注意不用IDCT那种含虚部的做法DCT-II系数的实数域操作更常用。2.2 滤波器组的构造细节梅尔滤波器组是整个流程的核心也是很多人最容易写错的地方。设计思路是确定采样率fs和FFT点数NFFT把频率范围0到fs/2转换成梅尔刻度在这段梅尔轴上均匀取M2个点再反变换回频率轴上对应的物理频率每个三角滤波器的起点、峰值点、终点就是这M2个点中相邻的三个。有个细节值得专门拎出来讲滤波器组的起始端点和末端端点尤其是底部的第一个滤波器起点在0Hz附近这在音频中会引入直流分量。实际数据中直流倒谱系数C0常常反映的是信号能量有的应用会保留有的会丢弃。做说话人识别时C0保留有价值做语音识别时丢掉也不影响太多。我通常会把C0单独拿出来作为能量特征和其他系数一起使用效果更灵活。3. Matlab核心代码实现3.1 一个完整可运行的MFCC函数先给出一份我一直在用的核心提取函数。我没有依赖任何第三方工具箱只用Matlab自带功能基于2022b版本开发验证老版本也都兼容。代码思路参考了HTK的实现但做了简化方便理解。function [mfcc, deltaMfcc, ddeltaMfcc] my_mfcc(audio, fs, params) % 输入 audio: 一维列向量fs: 采样率params: 结构体参数 % 输出 mfcc: 每帧特征向量行数为帧数 % ---------- 默认参数 ---------- if nargin 3 || isempty(params) params.frameLen 25; % 帧长 ms params.frameShift 10; % 帧移 ms params.numFilters 26; % 梅尔滤波器数量 params.numCoeffs 13; % 保留的低维倒谱系数个数不含C0 params.freqRange [0, fs/2]; % 频率范围 params.filterAlpha 0.97; % 预加重系数 params.dctNorm ortho; % DCT正交归一化 end % ---------- 预加重 ---------- audio filter([1, -params.filterAlpha], 1, audio); % ---------- 分帧 ---------- frameLenSmp round(params.frameLen / 1000 * fs); frameShiftSmp round(params.frameShift / 1000 * fs); [frames, numFrames] frame_signal(audio, frameLenSmp, frameShiftSmp); % ---------- 加窗 ---------- win hamming(frameLenSmp, periodic); frames_win frames .* repmat(win, 1, numFrames); % ---------- FFT和幅度谱 ---------- nfft 2^nextpow2(frameLenSmp); magSpec abs(fft(frames_win, nfft, 1)); magSpec magSpec(1:nfft/21, :); % 只取单边谱 % ---------- 梅尔滤波器组 ---------- melBank mel_filterbank(nfft, fs, params.numFilters, params.freqRange); melSpec melBank * magSpec; % ---------- 对数能量 ---------- logMelSpec log(melSpec eps); % ---------- DCT ---------- dctMat dctmtx(params.numFilters); % 生成DCT矩阵 mfccRaw dctMat(1:params.numCoeffs1, :) * logMelSpec; mfcc mfccRaw; % 转置成每行一帧 % ---------- 一阶和二阶差分 ---------- if nargout 1 deltaMfcc delta_coeff(mfcc, 2); end if nargout 2 ddeltaMfcc delta_coeff(deltaMfcc, 2); end end3.2 辅助函数分帧和滤波器组设计接下来是配套的分帧函数。有人喜欢用内建的buffer函数但那个函数对数据尾部不足一帧的处理太粗暴我这里先做了一个向后兼容性好的版本逻辑简单性能也不错。function [frames, numFrames] frame_signal(x, frameLen, frameShift) lenx length(x); numFrames floor((lenx - frameLen) / frameShift) 1; indices (0:frameLen-1); frameIdx 1 (0:numFrames-1) * frameShift; frames zeros(frameLen, numFrames); for i 1:numFrames frames(:, i) x(indices frameIdx(i), 1); end end然后用代码实现梅尔滤波器组的生成。这里采用纯手写方式方便检查和修改。关键是把频率映射到梅尔轴再均匀分成M2段画出一组重叠的三角滤波器。function melBank mel_filterbank(nfft, fs, numFilters, freqRange) % 定义梅尔映射 mel (f) 2595 * log10(1 f / 700); invMel (m) 700 * (10.^(m / 2595) - 1); fmin freqRange(1); fmax freqRange(2); melMin mel(fmin); melMax mel(fmax); melPoints linspace(melMin, melMax, numFilters 2); hzPoints invMel(melPoints); % 映射到FFT bin索引 binPoints floor((nfft 1) * hzPoints / fs); melBank zeros(numFilters, nfft/2 1); for m 2:numFilters1 leftBin binPoints(m-1); centerBin binPoints(m); rightBin binPoints(m1); for k max(1, leftBin):centerBin if leftBin centerBin melBank(m-1, k) 1; else melBank(m-1, k) (k - leftBin) / (centerBin - leftBin); end end for k centerBin1:rightBin melBank(m-1, k) (rightBin - k) / (rightBin - centerBin); end end end最后是一阶差分的实现这个在语音识别里非常重要单纯静态系数对音调变化不敏感加上差分就包含了动态特征。function delta delta_coeff(mfcc, deltaWindow) [numFrames, dim] size(mfcc); delta zeros(size(mfcc)); for i 1:numFrames denom 0; for t -deltaWindow:deltaWindow denom denom t^2; if i t 1 i t numFrames delta(i, :) delta(i, :) t * mfcc(i t, :); end end if denom 0 delta(i, :) delta(i, :) / denom; end end end3.3 主程序调用演示有了上面的函数主程序就很清爽了。我常用两种读音频的方式一种是处理小的wav文件一种是处理从麦克风采集的实时数据。这里给出文件路径的版本% 主脚本 [audio, fs] audioread(speech.wav); % 若多声道则取单声道 if size(audio, 2) 1 audio mean(audio, 2); end % 提取前先归一化到[-1, 1]audioread本身已经是归一化数据但有时来自别处的是int类型 audio double(audio) / (max(abs(audio)) eps); params.frameLen 25; params.frameShift 10; params.numFilters 26; params.numCoeffs 12; % 加上C0一共13维 [mfcc, delta, ddelta] my_mfcc(audio, fs, params); % 可视化特征热图 figure; imagesc(mfcc); xlabel(帧号); ylabel(MFCC维度); colorbar; title(MFCC特征热图); % 保存为文本方便喂给分类器 save(mfcc_feature.mat, mfcc, delta, ddelta, fs);我给这个函数测试过一段4秒的录音采样率16kHz帧长25ms帧移10ms。算下来一帧400个采样点总帧数约390帧在普通笔记本上整个提取过程不到0.3秒。把这13维静态、13维一阶差分、13维二阶差分拼接起来特征维度就是39维这就是很多经典语音识别系统所谓的标准39维特征。4. 实际测试效果和参数怎么调4.1 不同采样率和音频长度的表现MFCC对采样率不算特别敏感但从16kHz降到8kHz时高频信息会被截断影响识别的上限。实际项目中电话语音一般是8kHz这种情况下滤波器组只覆盖0到4kHz识别电话ARS系统完全够用。宽带语音会议系统用16kHz更稳妥。音频长度上一句话2秒到3秒的分帧数量足够提取几十帧可靠特征。但如果你做的是一个唤醒词检测语音片段可能只有0.8秒帧数严重不足这时候建议把帧移从10ms缩短到5ms或者把滤波器数量从26降到20保证每帧特征的稳定性。我在一个车载唤醒词项目里做过对比把帧移缩短后识别准确率提升了大概2.5个百分点。4.2 参数选择的经验和原则帧长和帧移是最先需要敲定的参数。25ms/10ms是语音识别界的黄金组合基本覆盖了绝大多数场景。遇到说话速度特别快的语音比如儿童朗读或者语速极快的方言可以试试20ms/5ms。滤波器的数量决定了频带分辨率一般语音识别用26到40个。如果目的是提取特征去喂给深度神经网络做声学模型用40个滤波器会带来更好的频域细节因为DNN建模能力强能吃下高维特征。滤波器输出取对数以后建议加一个eps非线性保护不然会出现log(0)的NaN。这是我见过新手最多的问题之一在纯静音段尤其容易出现。4.3 对噪声的鲁棒性调节还要提一嘴CMVN也就是倒谱均值方差归一化。MFCC的系数序列在时间轴上做一个简单的高通滤波或减去均值能大幅减轻通道效应和加性噪声的影响。实际使用中我习惯在原始MFCC后面紧跟一步CMVN再做差分很多公开语音数据集上这一步能把识别率提升好几个点代价几乎为零。Matlab里实现就是一行的均值减除mfccCmn (mfcc - mean(mfcc)) ./ (std(mfcc) eps);5. 常见问题排查与避坑记录5.1 高频问题速查表我整理一下这些年实际踩过、也见别人踩过的坑按出现频率排了个序。很多问题表面看起来代码没问题但细查下来都是这些老套路。现象可能原因排查方法所有帧的MFCC几乎全零分帧时索引没有加1第一帧起点错位检查frameIdx是否从1开始特征图有条纹状干扰没有加窗或用了矩形窗频谱泄漏改用汉明窗或汉宁窗某一帧出现Inf或NaN信号中有NaN或log(0)加eps保护输入前先检查滤波器组面积不相等滤波器边界用round取整时导致的能量偏差对每个滤波器做归一化除以面积差分系数跳变剧烈差分窗口参数设置过大或数据有脉冲噪声用deltaWindow2并做平滑高频特征能量异常偏大预加重系数太小高频被放大过头检查filterAlpha是否误设成1.2这类调用dctmtx报错老版本Matlab没有该函数手动构造DCT矩阵或改用dct函数5.2 我重点想说的几个Bug细节第一个坑是分帧时音频不够一帧。分帧代码里numFrames的计算默认了音频足够长如果音频时长小于一帧frameLenSmp比lenx还大indices frameIdx(i)直接越界报错。实际处理短语音片段时一定要先补零if lenx frameLenSmp audio [audio; zeros(frameLenSmp - lenx, 1)]; end第二个坑是滤波器组的频率范围要和采样率匹配。假如你的fs是44100Hz但freqRange写死成0到8000那么8000Hz以上的能量全被丢掉了这不是在高频做抑制而是直接丢失信息。正确做法是在调用时用freqRange [0, fs/2]或者按实际有效带宽动态计算。第三个坑是Matlab自带的dct函数和dctmtx生成矩阵的差异。dct(x)是按信号处理的DCT-II定义计算输出长度和输入长度一致dctmtx(n)是生成一个n×n的正交矩阵使用时通常取前L行做截断。两者在系数缩放上有细微差异用dctmtx时如果选了ortho归一化得到的数值和标准HTK中mfcc的幅度会差一个倍数。我自己习惯统一用dctmtx加上ortho归一化这样各维度能量比较均衡输入分类器更友好。还有一个很隐蔽的问题是FFT点数和帧长不一致。通常nfft 2^nextpow2(frameLenSmp)也就是比帧长略大一点的2的幂次这没有问题。但如果帧长本身就是2的幂次nextpow2返回的还是同一数量级nfft等于frameLenSmp此时单边谱长度是frameLenSmp/21滤波器组矩阵的列数必须动态适配。我在设计melBank时直接用了nfft/21所以整个链路是自洽的。你要是自己改帧长就要同步检查这些维度。要是你用的是旧版Matlab还没有Signal Processing Toolboxhamming函数会报错。替代方案是手动生成汉明窗0.54 - 0.46cos(2pi*(0:N-1)/(N-1))我早期做嵌入式验证时就是这么干的性能也没有差别。5.3 一个实用的调试技巧我调MFCC代码时的经验是不要只盯着最终特征图看要逐步检查中间结果。比如做完预加重之后画出前1000个采样点的波形高频部分应该明显抬升做完分帧加窗后随机挑一帧做FFT频谱包络应该没有明显的刺状旁瓣滤波器组输出那一层把melSpectrum画出来应该是一条条平滑的包络线。中间任何一步不对劲都能定位出来比最后看到一份乱七八糟的特征图再去倒推要省时间得多。结语我做MFCC这几年最大的感受就是它看起来公式多实际步骤拆开全是信号处理基本功的组合。只要肯花一个下午手写一遍预加重、分帧、滤波器组你对到底该保留哪些系数、放弃哪些系数会形成直觉。最后再分享一个小技巧训练分类器之前把每段音频的MFCC按帧做零均值归一化就是减均值除方差这个方法比调那些复杂的数据增强参数更管用实测在简单命令词识别任务里差错率能降一截。希望这份Matlab代码能让你少走点弯路。本文还有配套的精品资源点击获取