语音信号处理之语音特征提取（1）机器学习的语音处理

本文首先是将Speech Processing for Machine Learning: Filter banks, Mel-Frequency Cepstral Coefficients (MFCCs) and What's In-Between这篇2016年4月21号的文章翻译了一下，而后讲出一点本身的理解和注释。

机器学习的语音处理：滤波器组，梅尔频率倒谱系数（MFCC）以及二者之间的区别

语音处理在任何语音系统中都扮演着重要的角色，不管它是自动语音识别（ASR）仍是说话者识别等等。 长期以来，梅尔频率倒谱系数（MFCC）是很是受欢迎的功能。 可是最近，滤波器组正变得愈来愈流行。 在这篇文章中，我将讨论滤波器组和MFCC，以及为何滤波器组变得愈来愈流行。

计算滤波器组和MFCC涉及相同的过程，在这两种状况下，都将计算滤波器组，并经过一些额外的步骤能够得到MFCC。 简而言之，信号会经过预加剧滤波器。 而后将其切成（重叠的）帧，并将窗口函数应用于每一个帧； 以后，咱们在每一个帧上进行傅立叶变换（或更具体地说是短时傅立叶变换），并计算功率谱； 而后计算滤波器组。 为了得到MFCC，可将离散余弦变换（DCT）应用于滤波器组，以保留多个所得系数，而其他系数则被丢弃。 两种状况的最后一步是平均归一化。

 

创建

       在本文中，我今后处使用了一个16位PCM wav文件，称为“ OSR_us_000_0010_8k.wav”，其采样频率为8000 Hz。 wav文件是干净的语音信号，包括单个语音，在其中说出一些句子，中间有一些停顿。 为简单起见，我使用了信号的前3.5秒，大体至关于wav文件中的第一句话。
       我将使用Python 2.7.x，NumPy和SciPy。 这篇文章中使用的一些代码基于该存储库中可用的代码。

import numpy
import scipy.io.wavfile
from scipy.fftpack import dct

sample_rate, signal = scipy.io.wavfile.read('OSR_us_000_0010_8k.wav') # File assumed to be in the same directory
signal = signal[0:int(3.5 * sample_rate)] #Keep the first 3.5 seconds

'''
显示这个图片
plt.figure(figsize=(20,5)) # 修改坐标显示
plt.xlabel("Time(s)")
plt.ylabel("Amplitude")
plt.plot(x , signal)
plt.show()
'''

原始信号在时域中具备如下形式：

时域信号

强化前

第一步是对信号应用预加剧滤波器，以放大高频。 预加剧滤波器在几种方面颇有用：（1）平衡频谱，由于与低频相比，高频一般具备较小的幅度；（2）避免在傅立叶变换操做期间出现数值问题；（3）还可改善信号 噪声比（SNR）。

可使用如下等式中的一阶滤波器将预加剧滤波器应用于信号xx：

y（t）= x（t）-αx（t-1）y（t）= x（t）-αx（t-1）

可使用如下代码轻松实现，其中滤波器系数（αα）的典型值为0.95或0.97，pre_emphasis = 0.97：

emphasized_signal = numpy.append(signal[0], signal[1:] - pre_emphasis * signal[:-1])

预加剧在现代系统中的影响不大，主要是由于除避免了不该该成为问题的傅立叶变换数值问题外，大多数预加剧滤波器的动机均可以使用均值归一化（在本文后面讨论）来实现。 在现代FFT实现中。
预加剧后的信号在时域中具备如下形式：

预加剧后的时域信号

 

构图

       通过预加剧后，咱们须要将信号分红短时帧。 此步骤的基本原理是信号中的频率会随时间变化，所以在大多数状况下，对整个信号进行傅立叶变换是没有意义的，由于咱们会随时间丢失信号的频率轮廓。 为避免这种状况，咱们能够安全地假设信号中的频率在很短的时间内是固定的。 所以，经过在此短时帧上进行傅立叶变换，咱们能够经过串联相邻帧来得到信号频率轮廓的良好近似值。
       语音处理中的典型帧大小为20毫秒至40毫秒，连续帧之间有50％（+/- 10％）重叠。 流行的设置是帧大小为25毫秒，frame_size = 0.025和10毫秒跨度（重叠15毫秒），frame_stride = 0.01。

frame_length, frame_step = frame_size * sample_rate, frame_stride * sample_rate # Convert from seconds to samples
signal_length = len(emphasized_signal)
frame_length = int(round(frame_length))
frame_step = int(round(frame_step))
num_frames = int(numpy.ceil(float(numpy.abs(signal_length - frame_length)) / frame_step)) # Make sure that we have at least 1 frame

pad_signal_length = num_frames * frame_step + frame_length
z = numpy.zeros((pad_signal_length - signal_length))
pad_signal = numpy.append(emphasized_signal, z) # Pad Signal to make sure that all frames have equal number of samples without truncating any samples from the original signal

indices = numpy.tile(numpy.arange(0, frame_length), (num_frames, 1)) + numpy.tile(numpy.arange(0, num_frames * frame_step, frame_step), (frame_length, 1)).T
frames = pad_signal[indices.astype(numpy.int32, copy=False)]

 

窗口

将信号切成帧后，咱们对每一个帧应用诸如汉明窗之类的窗口函数。 Hamming窗口具备如下形式：

w [n] = 0.54-0.46cos（2πnN-1）w [n] = 0.54-0.46cos（2πnN-1）

其中，0≤n≤N-10≤n≤N-1，NN为窗口长度。 绘制先前的方程式将得出如下曲线：

Hamming Window

咱们有几个缘由须要将窗函数应用于这些帧，特别是要抵消FFT所作的数据无限的假设并减小频谱泄漏。

frames *= numpy.hamming(frame_length)
# frames *= 0.54 - 0.46 * numpy.cos((2 * numpy.pi * n) / (frame_length - 1)) # Explicit Implementation **

傅立叶变换和功率谱

如今，咱们能够在每一个帧上执行NN点FFT来计算频谱，这也称为短时傅立叶变换（STFT），其中NN一般为256或512，NFFT = 512； 而后使用如下公式计算功率谱（周期图）：

P = | FFT（xi）| 2NP = | FFT（xi）| 2N

其中，xixi是信号xx的第i个帧。 这能够用如下几行实现：

mag_frames = numpy.absolute(numpy.fft.rfft(frames, NFFT)) # Magnitude of the FFT
pow_frames = ((1.0 / NFFT) * ((mag_frames) ** 2)) # Power Spectrum

滤波器

计算滤波器组的最后一步是将三角滤波器（一般为40个滤波器，在Mel等级上为nfilt = 40）应用于功率谱以提取频带。 梅尔音阶的目的是模仿低频的人耳对声音的感知，方法是在较低频率下更具判别力，而在较高频率下则具备较少判别力。 咱们可使用如下公式在赫兹（ff）和梅尔（mm）之间转换：

m = 2595log10（1 + f700）m =2595log10⁡（1 + f700）

f = 700（10m / 2595-1）f = 700（10m / 2595-1）

滤波器组中的每一个滤波器都是三角形的，在中心频率处的响应为1，并朝着0线性减少，直到它到达两个相邻滤波器的中心频率处，响应为0，以下图所示：

梅尔秤上的滤波器组

能够经过如下方程式建模（今后处获取）：

Hm（k）=⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪0k-f（ m-1）f（m）-f（m-1）1f（m + 1）-kf（m + 1）-f（m）0k <f（m-1）f（m-1）≤k< f（m）k = f（m）f（m）<k≤f（m + 1）k> f（m + 1）Hm（k）= {0k <f（m-1）k-f（m -1）f（m）-f（m-1）f（m-1）≤k<f（m）1k = f（m）f（m + 1）-kf（m + 1）-f（m ）f（m）<k≤f（m + 1）0k> f（m + 1）

low_freq_mel = 0
high_freq_mel = (2595 * numpy.log10(1 + (sample_rate / 2) / 700)) # Convert Hz to Mel
mel_points = numpy.linspace(low_freq_mel, high_freq_mel, nfilt + 2) # Equally spaced in Mel scale
hz_points = (700 * (10**(mel_points / 2595) - 1)) # Convert Mel to Hz
bin = numpy.floor((NFFT + 1) * hz_points / sample_rate)

fbank = numpy.zeros((nfilt, int(numpy.floor(NFFT / 2 + 1))))
for m in range(1, nfilt + 1):
f_m_minus = int(bin[m - 1]) # left
f_m = int(bin[m]) # center
f_m_plus = int(bin[m + 1]) # right

for k in range(f_m_minus, f_m):
fbank[m - 1, k] = (k - bin[m - 1]) / (bin[m] - bin[m - 1])
for k in range(f_m, f_m_plus):
fbank[m - 1, k] = (bin[m + 1] - k) / (bin[m + 1] - bin[m])
filter_banks = numpy.dot(pow_frames, fbank.T)
filter_banks = numpy.where(filter_banks == 0, numpy.finfo(float).eps, filter_banks) # Numerical Stability
filter_banks = 20 * numpy.log10(filter_banks) # dB

将滤波器组应用于信号的功率谱（周期图）后，咱们获得如下频谱图：

信号的频谱图

若是梅尔级滤波器组是所需的功能，那么咱们能够跳过以归一化。

梅尔频率倒谱系数（MFCC）

       事实证实，在上一步中计算出的滤波器组系数是高度相关的，这在某些机器学习算法中可能会出现问题。 所以，咱们能够应用离散余弦变换（DCT）去相关滤波器组系数，并生成滤波器组的压缩表示。 一般，对于自动语音识别（ASR），结果倒谱系数2-13将保留，其他的将被丢弃； num_ceps =12。丢弃其余系数的缘由是，它们表明滤波器组系数的快速变化，而这些细微的细节对自动语音识别（ASR）毫无帮助。

mfcc = dct(filter_banks, type=2, axis=1, norm='ortho')[:, 1 : (num_ceps + 1)] # Keep 2-13

      能够将正弦提高器1应用于MFCC，以下降对较高MFCC的强调，这被认为能够改善噪声信号中的语音识别。

(nframes, ncoeff) = mfcc.shape
n = numpy.arange(ncoeff)
lift = 1 + (cep_lifter / 2) * numpy.sin(numpy.pi * n / cep_lifter)
mfcc *= lift #*

生成的MFCC：

MFCCs

平均归一化

如前所述，为了平衡频谱并改善信噪比（SNR），咱们能够简单地从全部帧中减去每一个系数的平均值。

filter_banks -= (numpy.mean(filter_banks, axis=0) + 1e-8)

均值归一化滤波器组：

归一化滤波器组

对于MFCC一样如此：

mfcc -= (numpy.mean(mfcc, axis=0) + 1e-8)

均值归一化MFCC：

标准化MFCC

滤波器组与MFCC

到目前为止，已根据其动机和实现方式讨论了计算滤波器组和MFCC的步骤。有趣的是，计算滤波器组所需的全部步骤都是由语音信号的性质和人类对此类信号的感知所激发的。相反，某些机器学习算法的局限性促使了计算MFCC所需的额外步骤。须要离散余弦变换（DCT）来使滤波器组系数去相关，该过程也称为白化。特别是，当高斯混合模型-隐马尔可夫模型（GMMs-HMM）很是流行而且MFCC和GMMs-HMM一块儿发展成为自动语音识别（ASR）2的标准方式时，MFCC便很受欢迎。随着语音系统中深度学习的出现，人们可能会质疑MFCC是否仍然是正确的选择，由于深度神经网络不太容易受到高度相关的输入的影响，所以离散余弦变换（DCT）再也不是必要的步骤。值得注意的是，离散余弦变换（DCT）是线性变换，所以是不但愿的，由于它会丢弃语音信号中的某些信息，这些信息是高度非线性的。
怀疑是否须要进行傅立叶变换是明智的。鉴于傅立叶变换自己也是线性运算，所以忽略它并尝试在时域中直接从信号中学习可能会有所帮助。确实，最近的一些工做已经尝试了这一点，并报告了积极的结果。可是，傅立叶变换操做是一个很难学习的操做，而且能够说会增长得到相同性能所需的数据量和模型复杂性。此外，在进行短时傅立叶变换（STFT）时，咱们假设信号在此短期内保持稳定，所以傅立叶变换的线性不会形成严重问题。

结论

在本文中，咱们探讨了计算梅尔级滤波器组和梅尔频率倒谱系数（MFCC）的过程。 讨论了该过程当中每一个步骤的动机和实现。 咱们也提出了与MFCC相比，滤波器组愈来愈受欢迎的缘由。
tl; dr：若是机器学习算法不易受到高度相关的输入的影响，请使用梅尔级滤波器组。 若是机器学习算法易受相关输入的影响，请使用MFCC。

1.提高在倒频谱域中进行过滤。 请注意，频谱和倒谱中的符号分别使用了滤波和提高。 ↩
2.本文对此主题进行了精彩的讨论。 ↩

重要名词理解：

滤波器组：

人耳对声音频谱的响应是非线性的，经验代表：若是咱们可以设计一种前端处理算法，以相似于人耳的方式对音频进行处理，能够提升语音识别的性能。FilterBank就是这样的一种算法。FBank特征提取要在预处理以后进行，这时语音已经分帧，咱们须要逐帧提取FBank特征。

傅里叶变换：

人类听觉系统其实就像是在作傅里叶变换，经过耳蜗（能够视为一组滤波器）将时域声压转换为独立的不一样频率的份量，可是，准确地说，更接近于短时傅里叶变换（STFT）。

傅里叶变换处理非平稳信号（频率随时间变化的信号）有天生缺陷。它只能获取一段信号整体上包含哪些频率的成分，可是对各成分出现的时刻并没有所知。所以时域相差很大的两个信号，可能频谱图同样。
对于这样的非平稳信号，只知道包含哪些频率成分是不够的，咱们还想知道各个成分出现的时间。知道信号频率随时间变化的状况，各个时刻的瞬时频率及其幅值——这也就是时频分析。
一个简单可行的方法就是——加窗。我又要套用方沁园同窗的描述了，“把整个时域过程分解成无数个等长的小过程，每一个小过程近似平稳，再傅里叶变换，就知道在哪一个时间点上出现了什么频率了。”这就是短时傅里叶变换。

 梅尔频率倒谱系数（MFCC）：

参考：https://zhuanlan.zhihu.com/p/27416870

均值归一化：

参考：https://www.jianshu.com/p/ef3534ddda15

第一部分构图用的完整代码：

import numpy
import scipy.io.wavfile
from scipy.fftpack import dct
import matplotlib as mpl
import matplotlib.pyplot as plt

sample_rate, signal = scipy.io.wavfile.read('OSR_us_000_0010_8k.wav') # File assumed to be in the same directory
signal = signal[0:int(3.5 * sample_rate)] #Keep the first 3.5 seconds


plt.figure(figsize=(20,5)) # 修改坐标显示
plt.xlabel("Time(s)")
plt.ylabel("Amplitude")
plt.plot(signal)
plt.show()

参考文章：https://www.jianshu.com/p/b416d5617b0c