背景：

报了深蓝学院的语音识别课程，这里做学习记录

第二课：语音信号处理-->特征提取

1.数字信号处理

• 1）傅立叶变换

2.常用特征提取

• 1）Fbank特征一般用于DNN训练
• 2）MFCC特征一般用于对角GMM训练，各维度之间相关性小

Step1.预加重（pre-emphasis）

• 为什么需要预加重？

提高信号高频部分的能量，高频信号在传递过程中，衰减较快，
但是高频部分又蕴含很多对语音识别有利的特征，
因此，在特征提取部分，需要提高高频部分能量

预加重滤波器是一个一阶高通滤波器，给定时域输入信号𝑥[𝑛]，预加重之后的信号为:
𝑦 [𝑛] = 𝑥[𝑛] − 𝛼𝑥[𝑛 − 1]

import librosa
import numpy as np
from scipy.fftpack import dct

#preemphasis config 
alpha = 0.97

# Enframe config
frame_len = 400      # 25ms, fs=16kHz
frame_shift = 160    # 10ms, fs=15kHz
fft_len = 512

# Mel filter config
num_filter = 23
num_mfcc = 12

# Read wav file
wav, fs = librosa.load('./test.wav', sr=None)

def preemphasis(signal, coeff=alpha):
    """perform preemphasis on the input signal.

        :param signal: The signal to filter.
        :param coeff: The preemphasis coefficient. 0 is no filter, default is 0.97.
        :returns: the filtered signal.
    """
    return np.append(signal[0], signal[1:] - coeff * signal[:-1])

Step2.加窗（windowing）分帧

• 为什么需要分帧？

语音信号为非平稳信号，其统计属性是随着时间变化的，以汉语为例，一句话中包含很多声母和韵母，不同的拼音，发音的特点很明显是不一样的；

• 但是！语音信号又具有短时平稳的属性，比如汉语里一个声母或者韵母，往往只会持续几十到几百毫秒，在这一个发音单元里，语音信号表现出明显的稳定性，规律性（可以自己使用Audition观察一段语音）
  •在进行语音识别的时候，对于一句话，识别的过程也是以较小的发音单元（音素、字、字节）为单位进行识别，因此用滑动窗来提取短时片段，

• 帧长、帧移、窗函数的概念，对于采样率为16kHz的信号，帧长、帧移一般为25ms、10ms，即400和160个采样点

def enframe(signal, frame_len=frame_len, frame_shift=frame_shift, win=np.hamming(frame_len)):
    """Enframe with Hamming widow function.

        :param signal: The signal be enframed
        :param win: window function, default Hamming
        :returns: the enframed signal, num_frames by frame_len array
    """
    
    num_samples = signal.size
    num_frames = np.floor((num_samples - frame_len) / frame_shift)+1
    frames = np.zeros((int(num_frames),frame_len))
    for i in range(int(num_frames)):
        frames[i,:] = signal[i*frame_shift:i*frame_shift + frame_len] 
        frames[i,:] = frames[i,:] * win

    return frames

Step3.傅里叶变换

• 将上一步分帧之后的语音帧，由时域变换到频域，取DFT系数的模，得到谱特征

def get_spectrum(frames, fft_len=fft_len):
    """Get spectrum using fft
        :param frames: the enframed signal, num_frames by frame_len array
        :param fft_len: FFT length, default 512
        :returns: spectrum, a num_frames by fft_len/2+1 array (real)
    """
    cFFT = np.fft.fft(frames, n=fft_len)
    valid_len = int(fft_len / 2 ) + 1
    spectrum = np.abs(cFFT[:,0:valid_len])
    return spectrum

Step4.梅尔滤波器组和对数操作

• DFT得到了每个频带上信号的能量，但是人耳对频率的感知不是等间隔的，近似于对数函数
• 将线性频率转换为梅尔频率，梅尔频率和线性频率转换关系

mel(f) = 2595 * log(1 + f /700)

• 梅尔三角滤波器组：根据起始频率、中间频率和截止频率，确定各滤波器系数

三角滤波图示（采样率为16k，对应的fbank的滤波外形）

filter_bank_m(k) = [ k - f(m-1) ] / [ f(m) - f(m-1)] , f(m-1) < k < f(m)
filter_bank_m(k) = [ f(m+1) - k] / [ f(m) - f(m-1)] , f(m) < k < f(m+1)

def fbank(spectrum, num_filter = num_filter):
    """Get mel filter bank feature from spectrum
        :param spectrum: a num_frames by fft_len/2+1 array(real)
        :param num_filter: mel filters number, default 23
        :returns: fbank feature, a num_frames by num_filter array 
        DON'T FORGET LOG OPRETION AFTER MEL FILTER!
    """

    feats=np.zeros((spectrum.shape[0], num_filter))
    freq_h = fs // 2 # 最大截止频率
    mel_h = f_mel(freq_h) # 最大梅尔刻度
    mel_lst = np.linspace(0, mel_h, num_filter+2) # mel刻度等间隔
    freq_lst = mel_f(mel_lst) # mel对应的频率
    fft_mel = np.floor((fft_len + 1) * freq_lst / fs) # mel对应fft的位置
    
    bank = np.zeros((num_filter, int(fft_len//2)+1))
    for i in range(1, num_filter+1):
        left = int(fft_mel[i-1]) # 三角波左侧
        center = int(fft_mel[i]) # 三角波顶点
        right = int(fft_mel[i+1]) # 三角波右侧
        for j in range(left, center):
            bank[i-1, j] = (j - left) / (center -left)
        for j in range(center, right):
            bank[i-1, j] = (right - j) / (right - center)

    feats = np.dot(spectrum, bank.T) # 频谱与滤波器点乘
    feats = np.where(feats == 0, np.finfo(float).eps, feats) # 避免log（0）
    return np.log(feats)

Step5.动态特征计算

• 一阶差分（Delta,Δ），类比速度，最简单的一阶差分计算方法
  Δ t = [c(t + 1) − c(t − 1)]/2

def mfcc(fbank, num_mfcc = num_mfcc):
    """Get mfcc feature from fbank feature
        :param fbank: a num_frames by  num_filter array(real)
        :param num_mfcc: mfcc number, default 12
        :returns: mfcc feature, a num_frames by num_mfcc array 
    """

    feats = dct(fbank, type=2, axis=1, norm="ortho")[:, 1 : (num_mfcc + 1)] 

    return feats

• 二阶差分（Delta delta, ΔΔ），类比加速度，简单计算方法
  ΔΔ t = [Δ(t + 1) − Δ(t − 1)]/2

3. 代码实践

主要代码在上面有详细实现，以下是一些实现的感想：

• 1）主要是写fbank的滤波器，mfcc就是在fbank的基础上做了一个离散余弦变换，可以直接调用scipy；

• 2）对于fbank，核心是三角滤波器, 直接按照公式，循环计算就可以了