【SLP·Python】基于 DTW GMM HMM 三种方法实现的语音分类识别系统-程序员宅基地

技术标签: python  机器学习  hmm  语音识别  

前言


Github 项目地址:https://github.com/Sherry-XLL/Digital-Recognition-DTW_HMM_GMM

在孤立词语音识别(Isolated Word Speech Recognition) 中,DTWGMMHMM 是三种典型的方法:

  • 动态时间规整(DTW, Dyanmic Time Warping)
  • 高斯混合模型(GMM, Gaussian Mixed Model)
  • 隐马尔可夫模型(HMM, Hidden Markov Model)

本文并不介绍这三种方法的基本原理,而是侧重于 Python 版代码的实现,针对一个具体的语音识别任务——10 digits recognition system,分别使用 DTW、GMM 和 HMM 建立对 0~9 十个数字的孤立词语音分类识别模型。

一、音频数据处理


1. 数据集

本文中的语音识别任务是对 0 ~ 9 这十个建立孤立词语言识别模型,所使用的数据集类似于语音版的 MNIST 手写数字数据库。原始数据集是 20 个人录制的数字 0 ~ 9 的音频文件,共 200 条 .wav 格式的录音,同样上传到了 Github 项目 中,分支情况如下:

records
digit_0
digit_1
digit_2
digit_3
digit_4
digit_5
digit_6
digit_7
digit_8
digit_9
1_0.wav
......
20_0.wav

records 文件夹下包括 digit_0 ~ digit_9 十个子文件夹。每个子文件夹代表一个数字,里面有 20 个人对该数字的录音音频,1_0.wav 就代表第 1 个人录数字 0 的音频文件。

为了比较不同方法的性能,我们将数据集按照 4:1 的比例分为训练集和测试集。我们更想了解的是模型在未知数据上的表现,于是前16个人的数据都被划分为了训练集,共 160 条音频文件;17-20 这四个人的音频为测试集,共 40 条 .wav 格式的录音。

2. 音频预处理

preprocess.py referencing
https://github.com/rocketeerli/Computer-VisionandAudio-Lab/tree/master/lab1

当录数字的音频时,录音前后会有微小时间的空隙,并且每段音频的空白片段长度不一。如果忽略这个差异,直接对原音频进行建模运算,结果难免会受到影响,因此本文选择基于双阈值的语音端点检测对音频进行预处理。

对于每段 .wav 音频,首先用 wave 读取,并通过 np.frombufferreadframes 得到二进制数组。然后,编写计算能量 (energy) 的函数 calEnergy 和计算过零率 (zero-crossing rate,ZCR) 的函数 calZeroCrossingRate,帧长 framesize 设置为常用的 256:

def calEnergy(wave_data):
    """
    :param wave_data: binary data of audio file
    :return: energy
    """
    energy = []
    sum = 0
    for i in range(len(wave_data)):
        sum = sum + (int(wave_data[i]) * int(wave_data[i]))
        if (i + 1) % 256 == 0:
            energy.append(sum)
            sum = 0
        elif i == len(wave_data) - 1:
            energy.append(sum)
    return energy


def calZeroCrossingRate(wave_data):
    """
    :param wave_data: binary data of audio file
    :return: ZeroCrossingRate
    """
    zeroCrossingRate = []
    sum = 0
    for i in range(len(wave_data)):
        sum = sum + np.abs(int(wave_data[i] >= 0) - int(wave_data[i - 1] >= 0))
        if (i + 1) % 256 == 0:
            zeroCrossingRate.append(float(sum) / 255)
            sum = 0
        elif i == len(wave_data) - 1:
            zeroCrossingRate.append(float(sum) / 255)
    return zeroCrossingRate

代入音频的二进制数据得到能量和过零率数组之后,通过语音端点检测得到处理之后的音频数据,端点检测函数 endPointDetect 主要为三个步骤:

  • 先根据 energy 和 zeroCrossingRate 计算高能阈值 (high energy threshold, MH)
  • 第一步是计算较高的短时能量作为高能阈值 (high energy threshold, MH),可用于识别语音的浊音部分,得到初步检测数据 A。
  • 第二步是计算较低的能量阈值 (low energy threshold, ML),使用该阈值搜索两端,第二次检测从 A 得到 B。
  • 第三步是计算过零率阈值 (zero crossing rate threshold, Zs),在考虑过零率的基础上进一步处理第二步的结果 B,返回处理好的数据 C。

最后,将编号 1-16 被试者的音频文件存入 processed_train_records ,编号 17-20 存入 processed_test_records

3. MFCC 特征提取

音频文件无法直接进行语音识别,需要对其进行特征提取。在语音识别(Speech Recognition)领域,最常用到的语音特征就是梅尔倒谱系数(Mel-scale Frequency Cepstral Coefficients,MFCC)。

在 Python 中,可以用封装好的工具包 librosapython_speech_features 实现对 MFCC 特征的提取。本文使用 librosa 提取给定音频的 MFCC 特征,每帧提取 39 维 MFCC 特征:

import librosa
def mfcc(wav_path, delta = 2):
    """
    Read .wav files and calculate MFCC
    :param wav_path: path of audio file
    :param delta: derivative order, default order is 2
    :return: mfcc
    """
    y, sr = librosa.load(wav_path)
    # MEL frequency cepstrum coefficient
    mfcc_feat = librosa.feature.mfcc(y=y, sr=sr, n_mfcc=13)
    ans = [mfcc_feat]
    
    # Calculate the 1st derivative
    if delta >= 1:
        mfcc_delta1 = librosa.feature.delta(mfcc_feat, order=1, mode ='nearest')
        ans.append(mfcc_delta1)
        
    # Calculate the 2nd derivative
    if delta >= 2:
        mfcc_delta2 = librosa.feature.delta(mfcc_feat, order=2, mode ='nearest')
        ans.append(mfcc_delta2)
        
    return np.transpose(np.concatenate(ans, axis=0), [1,0])

可以参考官网 librosapython-speech-features 了解更多音频特征提取工具。

二、Isolated Speech Recognition - DTW


1. DTW 算法步骤

动态时间规整(DTW, Dyanmic Time Warping) 本质上是一种简单的动态规划算法,它可以分为两个步骤,一个是计算两种模式的帧之间的距离,即得到帧匹配距离矩阵;另一个是在帧匹配距离矩阵中找到最优路径。算法步骤如下:

  • 首先,使用欧几里德距离 (Euclidean distance) 初始化成本矩阵 (cost matrix)。
  • 其次,使用动态规划计算成本矩阵,将每个数据的向量与模板向量进行比较,动态转移方程为
    c o s t [ i , j ] + = m i n ( [ c o s t [ i , j ] , c o s t [ i + 1 , j ] , c o s t [ i , j + 1 ] ] ) . cost[i, j] += min([cost[i, j], cost[i+1, j], cost[i, j+1]]). cost[i,j]+=min([cost[i,j],cost[i+1,j],cost[i,j+1]]).
  • 第三,计算规整路径 (warp path),最小距离为标准化距离。将第一个训练序列设置为标准 (template),合并所有训练序列并求其平均值,返回 0-9 位数字的模型。
  • 测试过程中,使用 DTW 衡量每条音频到 0-9 这十个模板的”距离“,并选择模板距离最小的数字作为音频的预测值。

2. DTW 函数编写

编写 dtw 函数的关键在于两段 MFCC 序列的动态规划矩阵和规整路径的计算。令输入的两段 MFCC 序列分别为 M1 和 M2,M1_len 和 M2_len 表示各自的长度,则 cost matrix 的大小为 M1_len * M2_len,先用欧式距离对其进行初始化,再根据转移式计算 cost matrix:

def dtw(M1, M2):
    """
    Compute Dynamic Time Warping(DTW) of two mfcc sequences.
    :param M1, M2: two mfcc sequences
    :return: the minimum distance and the wrap path
    """
    # length of two sequences
    M1_len = len(M1)
    M2_len = len(M2)
    cost_0 = np.zeros((M1_len + 1, M2_len + 1))
    cost_0[0, 1:] = np.inf
    cost_0[1:, 0] = np.inf
    
    # Initialize the array size to M1_len * M2_len
    cost = cost_0[1:, 1:]
    for i in range(M1_len):
        for j in range(M2_len):
            cost[i, j] = calEuclidDist(M1[i], M2[j])
            
    # dynamic programming to calculate cost matrix
    for i in range(M1_len):
        for j in range(M2_len):
            cost[i, j] += min([cost_0[i, j], \
                               cost_0[min(i + 1, M1_len - 1), j], \
                               cost_0[i, min(j + 1, M2_len - 1)]])

欧式距离的计算函数 calEuclidDist 可以用一行代码完成:

def calEuclidDist(A, B):
    """
    :param A, B: two vectors
    :return: the Euclidean distance of A and B
    """
    return sqrt(sum([(a - b) ** 2 for (a, b) in zip(A, B)]))

cost matrix 计算完成后还需计算 warp path,MFCC 序列长度为 1 时的路径可以单独分情况讨论,最小代价的 path_1 和 path_2 即为所求,最后返回数组 path:

# calculate the warp path
   if len(M1) == 1:
       path = np.zeros(len(M2)), range(len(M2))
   elif len(M2) == 1:
       path = range(len(M1)), np.zeros(len(M1))
   else:
       i, j = np.array(cost_0.shape) - 2
       path_1, path_2 = [i], [j]
       # path_1, path_2 with the minimum cost is what we want
       while (i > 0) or (j > 0):
           arg_min = np.argmin((cost_0[i, j], cost_0[i, j + 1], cost_0[i + 1, j]))
           if arg_min == 0:
               i -= 1
               j -= 1
           elif arg_min == 1:
               i -= 1
           else:
               j -= 1
           path_1.insert(0, i)
           path_2.insert(0, j)
       # convert to array
       path = np.array(path_1), np.array(path_2)

3. 模型训练与预测

在模型训练阶段,对于每个数字的 16 条音频文件,现计算 MFCC 序列到 mfcc_list 列表中,将第一个音频文件的 MFCC 序列设置为标准,计算标准模板和每条模板之间的动态时间规整路径,再对其求平均进行归一化,将结果 appendmodel 列表中。返回的 model 包含 0-9 这 10 个数字的 DTW 模型:

 # set the first sequence as standard, merge all training sequences
 mfcc_count = np.zeros(len(mfcc_list[0]))
 mfcc_all = np.zeros(mfcc_list[0].shape)
 for i in range(len(mfcc_list)):
     # calculate the wrap path between standard and each template
     _, path = dtw(mfcc_list[0], mfcc_list[i])
     for j in range(len(path[0])):
         mfcc_count[int(path[0][j])] += 1
         mfcc_all[int(path[0][j])] += mfcc_list[i][path[1][j]]

 # Generalization by averaging the templates
 model_digit = np.zeros(mfcc_all.shape)
 for i in range(len(mfcc_count)):
     for j in range(len(mfcc_all[i])):
         model_digit[i][j] = mfcc_all[i][j] / mfcc_count[i]
         
 # return model with models of 0-9 digits
 model.append(model_digit)

测试阶段比较简单,将训练得到的 model 和需要预测音频的路径输入到函数 predict_dtw 中,分别计算 model 中每位数字到音频 MFCC 序列的最短距离 min_distmin_dist 所对应的数字即为该音频的预测标签。

def predict_dtw(model, file_path):
    """
    :param model: trained model
    :param file_path: path of .wav file
    :return: digit
    """
    # Iterate, find the digit with the minimum distance
    digit = 0
    min_dist, _ = dtw(model[0], mfcc_feat)
    for i in range(1, len(model)):
        dist, _ = dtw(model[i], mfcc_feat)
        if dist < min_dist:
            digit = i
            min_dist = dist
    return str(digit)

三、Isolated Speech Recognition - GMM


1. GMM 模型训练流程

参考 sklearn 高斯混合模型中文文档:https://www.scikitlearn.com.cn/0.21.3/20/#211

高斯混合模型(GMM, Gaussian Mixed Model) 是一个假设所有的数据点都是生成于有限个带有未知参数的高斯分布所混合的概率模型,包含了数据的协方差结构以及隐高斯模型的中心信息,同样可以用于解决本文的数字识别任务。

scikit-learn 是基于 Python 语言的机器学习工具,sklearn.mixture 是其中应用高斯混合模型进行非监督学习的包,GaussianMixture 对象实现了用来拟合高斯混合模型的期望最大化 (EM, Expectation-Maximum) 算法。

如果想要直接调用封装好的库,只需在代码中加入 from sklearn.mixture import GaussianMixtureGaussianMixture.fit 便可以从训练数据中拟合出一个高斯混合模型,并用 predict 进行预测,详情参见官方文档

GMM 模型的训练和预测过程与 DTW 类似,此处不再赘述,流程图如下所示:
GMM 训练流程图

2. EM 算法步骤

参考文档:https://blog.csdn.net/nsh119/article/details/79584629?spm=1001.2014.3001.5501

GMM 模型的核心在于期望最大化 (EM, Expectation-Maximum) 算法,模型参数的训练步骤主要为 Expectation 的 E 步和 Maximum 的 M 步:

(1) 取参数的初始值开始迭代。
(2) E 步:依据当前模型参数,计算分模型 k k k 对观测数据 y j y_j yj 的相应度
γ ^ j k = α j ϕ ( y j ∣ θ k ) Σ k = 1 K α k ϕ ( y j ∣ θ k ) , j = 1 , 2 , . . . , N ; k = 1 , 2 , . . . , K \hat{\gamma}_{jk}=\dfrac{\alpha_j\phi(y_j|\theta_k)}{\Sigma^K_{k=1}\alpha_k\phi(y_j|\theta_k)},j=1,2,...,N; k=1,2,...,K γ^jk=Σk=1Kαkϕ(yjθk)αjϕ(yjθk)j=1,2,...,N;k=1,2,...,K
(3) M 步:计算新一轮迭代的模型参数
μ ^ k = Σ j = 1 N γ ^ j k y j Σ j = 1 N γ ^ j k , k = 1 , 2 , . . . , K σ ^ k 2 = Σ j = 1 N γ ^ j k ( y j − μ k ) 2 Σ j = 1 N γ ^ j k , k = 1 , 2 , . . . , K α ^ k = Σ j = 1 N γ ^ j k N , k = 1 , 2 , . . . , K \begin{aligned} \hat{\mu}_k&=\dfrac{\Sigma^N_{j=1}\hat{\gamma}_{jk}y_j}{\Sigma^N_{j=1}\hat{\gamma}_{jk}},k=1,2,...,K \\ \hat{\sigma}^2_k&=\dfrac{\Sigma^N_{j=1}\hat{\gamma}_{jk}(y_j-\mu_k)^2}{\Sigma^N_{j=1}\hat{\gamma}_{jk}},k=1,2,...,K \\ \hat{\alpha}_k&=\dfrac{\Sigma^N_{j=1}\hat{\gamma}_{jk}}{N},k=1,2,...,K \end{aligned} μ^kσ^k2α^k=Σj=1Nγ^jkΣj=1Nγ^jkyjk=1,2,...,K=Σj=1Nγ^jkΣj=1Nγ^jk(yjμk)2k=1,2,...,K=NΣj=1Nγ^jkk=1,2,...,K

3. GMM 实现代码

sklearn.mixture/_gaussian_mixture.py:https://github.com/scikit-learn/scikit-learn/blob/main/sklearn/mixture/_gaussian_mixture.py

GMM 模型实现可以参考 GaussianMixture源码,本文对 sklearn 的实现方式进行了简化,关键在于 GMM 这个类的函数实现。

我们实现的 GMM 类包含三个参数 (parameter) 和三个属性 (attribute):n_components 表示模型中状态的数目;mfcc_data 为音频提取的 MFCC 数据;random_state 是控制随机数的参数;means 表示在每个状态中 mixture component 的平均值;var 为方差;weights 是每个 component 的参数矩阵。

class GMM:
    """Gaussian Mixture Model.
    Parameters
    ----------
    n_components : int
        Number of states in the model.
    mfcc_data : array, shape (, 39)
        Mfcc data of training wavs.
    random_state: int
        RandomState instance or None, optional (default=0)
    Attributes
    ----------
    means : array, shape (n_components, 1)
        Mean parameters for each mixture component in each state.
    var : array, shape (n_components, 1)
        Variance parameters for each mixture component in each state.
    weights : array, shape (1, n_components)
        Weights matrix of each component.
    """

根据均值和方差可以计算对数高斯概率:

def log_gaussian_prob(x, means, var):
    """
    Estimate the log Gaussian probability
    :param x: input array
    :param means: The mean of each mixture component.
    :param var: The variance of each mixture component.
    :return: the log Gaussian probability
    """
    return (-0.5 * np.log(var) - np.divide(np.square(x - means), 2 * var) - 0.5 * np.log(2 * np.pi)).sum()

GMM 类中包含五个部分:__init__ 用于初始化;e_step 是 Expectation 步骤;m_step 是 Maximization 步骤;train 用于调用 E 步和 M 步进行参数训练;log_prob 计算每个高斯模型的对数高斯概率。

def __init__(self, mfcc_data, n_components, random_state=0):
    # Initialization
    self.mfcc_data = mfcc_data
    self.means = np.tile(np.mean(self.mfcc_data, axis=0), (n_components, 1))
    # randomization
    np.random.seed(random_state)
    for k in range(n_components):
        randn_k = np.random.randn()
        self.means[k] += 0.01 * randn_k * np.sqrt(np.var(self.mfcc_data, axis=0))
    self.var = np.tile(np.var(self.mfcc_data, axis=0), (n_components, 1))
    self.weights = np.ones(n_components) / n_components
    self.n_components = n_components

def e_step(self, x):
    """
    Expectation-step.
    :param x: input array-like data
    :return: Logarithm of the posterior probabilities (or responsibilities) of
             the point of each sample in x.
    """
    log_resp = np.zeros((x.shape[0], self.n_components))
    for i in range(x.shape[0]):
        log_resp_i = np.log(self.weights)
        for j in range(self.n_components):
            log_resp_i[j] += log_gaussian_prob(x[i], self.means[j], self.var[j])
        y = np.exp(log_resp_i - log_resp_i.max())
        log_resp[i] = y / y.sum()
    return log_resp

def m_step(self, x, log_resp):
    """
    Maximization step.
    :param x: input array-like data
    :param log_resp: Logarithm of the posterior probabilities (or responsibilities) of
                  the point of each sample in data
    """
    self.weights = np.sum(log_resp, axis=0) / np.sum(log_resp)
    denominator = np.sum(log_resp, axis=0, keepdims=True).T
    means_num = np.zeros_like(self.means)
    for k in range(self.n_components):
        means_num[k] = np.sum(np.multiply(x, np.expand_dims(log_resp[:, k], axis=1)), axis=0)
    self.means = np.divide(means_num, denominator)
    var_num = np.zeros_like(self.var)
    for k in range(self.n_components):
        var_num[k] = np.sum(np.multiply(np.square(np.subtract(x, self.means[k])),
                                        np.expand_dims(log_resp[:, k], axis=1)), axis=0)
    self.var = np.divide(var_num, denominator)

def train(self, x):
    """
    :param x: input array-like data
    """
    log_resp = self.e_step(x)
    self.m_step(x, log_resp)

def log_prob(self, x):
    """
    :param x: input array-like data
    :return: calculated log gaussian probability of single modal Gaussian
    """
    sum_prob = 0
    for i in range(x.shape[0]):
        prob_i = np.array([np.log(self.weights[j]) + log_gaussian_prob(x[i], self.means[j], self.var[j])
                           for j in range(self.n_components)])
        sum_prob += np.max(prob_i)
    return sum_prob

四、Isolated Speech Recognition - HMM


1. 隐马尔可夫模型

参考 hmmlearn 文档:https://hmmlearn.readthedocs.io/en/stable

隐马尔可夫模型(HMM, Hidden Markov Model) 是一种生成概率模型,其中可观测的 X X X 变量序列由内部隐状态序列 Z Z Z 生成。隐状态 (hidden states) 之间的转换假定为(一阶)马尔可夫链 (Markov chain) 的形式,可以由起始概率向量 π π π 和转移概率矩阵 A A A 指定。可观测变量的发射概率 (emission probability) 可以是以当前隐藏状态为条件的任何分布,参数为 θ \theta θ 。HMM 完全由 π π π A A A θ \theta θ 决定。

HMM 有三个基本问题:

  • 给定模型参数和观测数据,估计最优隐状态序列 (estimate the optimal sequence of hidden states)
  • 给定模型参数和观测数据,计算模型似然 (calculate the model likelihood)
  • 仅给出观测数据,估计模型参数 (estimate the model parameters)

第一个和第二个问题可以通过动态规划中的维特比算法 (Viterbi algorithm) 和前向-后向算法 (Forward-Backward algorithm) 来解决,最后一个问题可以通过迭代期望最大化 (EM, Expectation-Maximization) 算法求解,也称为 Baum-Welch 算法。

如果想要直接调用封装好的库,只需在代码中加入 from hmmlearn import hmm,详情参见官方文档源码

2. 代码实现

对于本文的数字语音识别任务,HMM 模型的训练预测过程与 DTW 和 GMM 类似,不再重复说明。HMM 模型实现可以参考 hmmlearn源码,本文的实现关键在于 HMM 这个类中的函数。与 GMM 的 EM 算法不同的是,在 HMM 的实现中我们使用的是维特比算法 (Viterbi algorithm) 和前向-后向算法 (Forward-Backward algorithm)。

HMM 类包含两个参数 (parameter) 和四个属性 (attribute):n_components 表示模型中状态的数目;mfcc_data 为音频提取的 MFCC 数据;means 表示在每个状态中 mixture component 的平均值;var 为方差;trans_mat 是状态之间的转移矩阵。

class HMM():
    """Hidden Markov Model.
    Parameters
    ----------
    n_components : int
        Number of states in the model.
    mfcc_data : array, shape (, 39)
        Mfcc data of training wavs.
    Attributes
    ----------
    init_prob : array, shape (n_components, )
        Initial probability over states.
    means : array, shape (n_components, 1)
        Mean parameters for each mixture component in each state.
    var : array, shape (n_components, 1)
        Variance parameters for each mixture component in each state.
    trans_mat : array, shape (n_components, n_components)
        Matrix of transition probabilities between states.

HMM 类中包含五个部分:

  • __init__ 对参数进行初始化 (Initialization);
  • viterbi Viterbi 算法用于解决 HMM 的解码问题,即找到观测序列中最可能的标记序列,本质上是一种动态规划算法 (dynamic programming);
  • forward 计算所有时间步中所有状态的对数概率 (log probabilities);
  • forward_backward 是参数估计的前向-后向算法 (Forward-Backward algorithm);
  • train 用于调用 forward_backwardviterbi 进行参数训练;
  • log_prob 计算每个模型的对数概率 (log likelihood)。
def viterbi(self, data):
    # Viterbi algorithm is used to solve decoding problem of HMM
    # That is to find the most possible labeled sequence of the observation sequence
    for index, single_wav in enumerate(data):
        n = single_wav.shape[0]
        labeled_seq = np.zeros(n, dtype=int)
        log_delta = np.zeros((n, self.n_components))
        psi = np.zeros((n, self.n_components))
        log_delta[0] = np.log(self.init_prob)
        for i in range(self.n_components):
            log_delta[0, i] += log_gaussian_prob(single_wav[0], self.means[i], self.var[i])
        log_trans_mat = np.log(self.trans_mat)
        for t in range(1, n):
            for j in range(self.n_components):
                temp = np.zeros(self.n_components)
                for i in range(self.n_components):
                    temp[i] = log_delta[t - 1, i] + log_trans_mat[i, j] + log_gaussian_prob(single_wav[t], self.means[j], self.var[j])
                log_delta[t, j] = np.max(temp)
                psi[t, j] = np.argmax(log_delta[t - 1] + log_trans_mat[:, j])
        labeled_seq[n - 1] = np.argmax(log_delta[n - 1])
        for i in reversed(range(n - 1)):
            labeled_seq[i] = psi[i + 1, labeled_seq[i + 1]]
        self.states[index] = labeled_seq

def forward(self, data):
    # Computes forward log-probabilities of all states at all time steps.
    n = data.shape[0]
    m = self.means.shape[0]
    alpha = np.zeros((n, m))
    alpha[0] = np.log(self.init_prob) + np.array([log_gaussian_prob(data[0], self.means[j], self.var[j]) for j in range(m)])
    for i in range(1, n):
        for j in range(m):
            alpha[i, j] = log_gaussian_prob(data[i], self.means[j], self.var[j]) + np.max(
                np.log(self.trans_mat[:, j].T) + alpha[i - 1])
    return alpha

def forward_backward(self, data):
    # forward backward algorithm to estimate parameters
    gamma_0 = np.zeros(self.n_components)
    gamma_1 = np.zeros((self.n_components, data[0].shape[1]))
    gamma_2 = np.zeros((self.n_components, data[0].shape[1]))
    for index, single_wav in enumerate(data):
        n = single_wav.shape[0]
        labeled_seq = self.states[index]
        gamma = np.zeros((n, self.n_components))
        for t, j in enumerate(labeled_seq[:-1]):
            self.trans_mat[j, labeled_seq[t + 1]] += 1
            gamma[t, j] = 1
        gamma[n - 1, self.n_components - 1] = 1
        gamma_0 += np.sum(gamma, axis=0)
        for t in range(n):
            gamma_1[labeled_seq[t]] += single_wav[t]
            gamma_2[labeled_seq[t]] += np.square(single_wav[t])
    gamma_0 = np.expand_dims(gamma_0, axis=1)
    self.means = gamma_1 / gamma_0
    self.var = (gamma_2 - np.multiply(gamma_0, self.means ** 2)) / gamma_0
    for j in range(self.n_components):
        self.trans_mat[j] /= np.sum(self.trans_mat[j])

五、实验结果


Github 项目地址:https://github.com/Sherry-XLL/Digital-Recognition-DTW_HMM_GMM

1. 实验环境

  macOS Catalina Version 10.15.6, Python 3.8, PyCharm 2020.2 (Professional Edition).

2. 文件简介

  dtw.py: Implementation of Dynamic Time Warping (DTW)
  gmm.py: Implementation of Gaussian Mixture Model (GMM)
  hmm.py: Implementation of Hidden Markov Model (HMM)

  gmm_from_sklearn.py: Train gmm model with GaussianMixture from sklearn
  hmm_from_hmmlearn.py: Train hmm model with hmm from hmmlearn

  preprocess.py: preprocess audios and split data
  processed_test_records: records with test audios
  processed_train_records: records with train audios
  records: original audios
  utils.py: utils function

若想运行代码,只需在命令行输入如下指令,以 DTW 为例:

  eg:
  python preprocess.py (mkdir processed records)
  python dtw.py 

3. 结果对比

各个方法的数据集和预处理部分完全相同,下面是运行不同文件的结果:

python dtw.py
----------Dynamic Time Warping (DTW)----------
Train num: 160, Test num: 40, Predict true num: 31
Accuracy: 0.78
python gmm_from_sklearn.py:
---------- GMM (GaussianMixture) ----------
Train num: 160, Test num: 40, Predict true num: 34
Accuracy: 0.85
python gmm.py
---------- Gaussian Mixture Model (GMM) ----------
confusion_matrix: 
 [[4 0 0 0 0 0 0 0 0 0]
 [0 4 0 0 0 0 0 0 0 0]
 [0 0 4 0 0 0 0 0 0 0]
 [0 0 0 4 0 0 0 0 0 0]
 [0 0 0 0 4 0 0 0 0 0]
 [0 0 0 0 0 4 0 0 0 0]
 [0 0 0 0 0 0 4 0 0 0]
 [0 0 0 0 0 0 0 4 0 0]
 [0 0 0 0 0 0 0 0 4 0]
 [0 0 0 0 0 0 0 0 0 4]]
Train num: 160, Test num: 40, Predict true num: 40
Accuracy: 1.00
python hmm_from_hmmlearn.py
---------- HMM (GaussianHMM) ----------
Train num: 160, Test num: 40, Predict true num: 36
Accuracy: 0.90
python hmm.py
---------- HMM (Hidden Markov Model) ----------
confusion_matrix: 
 [[4 0 0 0 0 0 0 0 0 0]
 [0 4 0 0 0 0 0 0 0 0]
 [0 0 3 0 1 0 0 0 0 0]
 [0 0 0 4 0 0 0 0 0 0]
 [0 0 0 0 4 0 0 0 0 0]
 [0 0 0 0 1 3 0 0 0 0]
 [0 0 0 0 0 0 3 0 1 0]
 [0 0 0 0 0 0 0 4 0 0]
 [0 0 0 1 0 0 1 0 2 0]
 [0 0 0 0 0 0 0 0 0 4]]
Train num: 160, Test num: 40, Predict true num: 35
Accuracy: 0.875
Method DTW GMM from sklearn Our GMM HMM from hmmlearn Our HMM
Accuracy 0.78 0.85 1.00 0.90 0.875

值得注意的是,上面所得正确率仅供参考。我们阅读源码就会发现,不同文件中 n_components 的数目并不相同,最大迭代次数 max_iter 也会影响结果。设置不同的超参数 (hyper parameter) 可以得到不同的正确率,上表并不是各种方法的客观对比。事实上 scikit-learn 的实现更完整详细,效果也更好,文中三种方法的实现仅为基础版本。

完整代码详见笔者的 Github 项目,如有问题欢迎在评论区留言指出。

六、参考链接


  1. https://github.com/Sherry-XLL/Digital-Recognition-DTW_HMM_GMM
  2. https://github.com/rocketeerli/Computer-VisionandAudio-Lab/tree/master/lab1
  3. http://librosa.github.io/librosa/core.html
  4. https://python-speech-features.readthedocs.io/
  5. http://yann.lecun.com/exdb/mnist/
  6. https://www.scikitlearn.com.cn/0.21.3/20/#211
  7. https://scikit-learn.org/stable/developers/advanced_installation.html#install-bleeding-edge
  8. https://blog.csdn.net/nsh119/article/details/79584629?spm=1001.2014.3001.5501
  9. https://github.com/scikit-learn/scikit-learn/blob/main/sklearn/mixture/_gaussian_mixture.py
  10. https://github.com/hmmlearn/hmmlearn
  11. https://hmmlearn.readthedocs.io/en/stable
  12. https://hmmlearn.readthedocs.io/en/stable/api.html#hmmlearn-hmm
版权声明:本文为博主原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接和本声明。
本文链接:https://blog.csdn.net/Sherry_ling/article/details/118713802

智能推荐

校验class文件的魔数,前四个字节_判断class魔数-程序员宅基地

文章浏览阅读1.4k次。在class文件开头的四个字节, 存放着class文件的魔数, 这个魔数是class文件的标志,他是一个固定的值: 0XCAFEBABE 。 也就是说他是判断一个文件是不是class格式的文件的标准, 如果开头四个字节不是0XCAFEBABE, 那么就说明它不是class文件, 不能被JVM识别。校验一下啦:package java2.os;import com.google.comm_判断class魔数

Linux命令之添加新用户useradd_useradd创建用户-程序员宅基地

文章浏览阅读1.7w次,点赞10次,收藏61次。`useradd` 命令可以用于添加新用户。_useradd创建用户

H3C设备系列问题-程序员宅基地

文章浏览阅读611次。一、h3c交换器和交换机的Telnet或SSH登录用户名和密码忘记了,怎么办?处理步骤:  1.使用Console线连接交换机或路由器的Console口,确保笔记本已连上设备,在设备启动过程中根据提示按crtl+B进入bootware界面,选择<8>跳过 Console 密码认证。    2.#选择<0> Reboot重启后进入命令行。  3.将Con..._failed to login because the idle timer expired

android-获取Bitmap的方法_android bitmap资源地址-程序员宅基地

文章浏览阅读569次。https://blog.csdn.net/taily_duan/article/details/52219607从资源存放路径分:(1)图片放在sdcard中Bitmap imageBitmap = BitmapFactory.decodeFile(path);// (path 是图片的路径,跟目录是/sdcard)..._android bitmap资源地址

Flowable简单了解_cmmn 容易理解-程序员宅基地

文章浏览阅读968次。Flowable简单了解基本介绍Flowable定义Flowable分类Flowable BPMN 业务流程引擎Flowable DMN 决策引擎Flowable CMMN 案例模型引擎Flowable Form 表单引擎基本使用流程图xml文件依赖设置基本代码基本介绍Flowable定义Flowable 是一个使用 Java 编写的轻量级业务流程引擎,使用 Apache V2 license 协议开源。2016 年 10 月,Activiti 工作流引擎的主要开发者离开 Alfresco 公司并在 _cmmn 容易理解

【公众号系列】12306购票送温暖-程序员宅基地

文章浏览阅读139次。公众号:SAP Technical本文作者:matinal原文出处:http://www.cnblogs.com/SAPmatinal/ 原文链接:【公众号系列】12306购票送温暖写在前面据说12306有了新功能,叫候补购票,有砖家表示,这种购票方式要比市面上所谓的抢票软件效率高、速度快。砖家的说法一定有道理,拿个小本本记下来。自从12306网上购票上...

随便推点

Without SSH/JSP/Servlet,不走寻常路,Java可以更酷-程序员宅基地

文章浏览阅读77次。标题的构思来源于Rod Johnson的那本"Without EJB"以及CCTV5中一句耳熟能详的广告词,不过此文并不是用来批判SSH(Struts、Spring、Hibernate)/JSP/Servlet的,也不是为某品牌做广告,而是用来分享这将近一年来的研究心得。去年圣诞节时曾在JavaEye发过一两篇文章,不过现在找不到了,文章内容提到要在3个月左右的时间内设计出一个..._without servlet

c语言数据结构图书馆系统_c语言数据结构实现图书管理系统-程序员宅基地

文章浏览阅读2.9k次,点赞7次,收藏37次。#include&amp;amp;amp;lt;stdio.h&amp;amp;amp;gt;#include&amp;amp;amp;lt;stdlib.h&amp;amp;amp;gt;#include&amp;amp;amp;lt;Windows.h&amp;amp;amp;gt;#include &amp;amp;amp;lt;string.h&amp;amp;amp;gt;#include &amp;amp_c语言数据结构实现图书管理系统

堆的创建和基本操作-程序员宅基地

文章浏览阅读1.6k次。1、堆的概念如果有一个关键码的集合 K={K0,K1,K2,…,Kn-1},把所有完全二叉树的顺序存储方式存储在一个一维数组中,并满足:Ki <= K2i +1 且 Ki <= K2i + 2 (Ki >= K2i+1 且 Ki >= K2i+2) i = 0,1,2…,则称为 小堆(或大堆)。将根节点最大的堆叫做最大堆或大根堆,根节点最小的堆叫做最小堆或小根堆。堆的性..._堆的创建

wps软件打不开共享超链接_wps版本word超链接无法打开指定文件-程序员宅基地

文章浏览阅读1.9w次。wps版本word超链接无法打开指定文件卡饭网本站整理2018-09-09Word文档中有时候经常使用超链接,这样有助阅读时快速跳转到指定的位置。然而某天发现word文档中插入的指向excel的超链接打不开了,怎么办呢方法一:链接文件位置发生变化,导致无法打开建议一般超链接的文件以及本文件放在一个文件夹中还可以通过设置超链接基础来解决这个问题。超链接基础在文件属性的摘要中,对于不同版本的offic..._wps发给别人超链接打不开

基于Flot可放缩的折线图-程序员宅基地

文章浏览阅读57次。Flot初步Flot是一个免费开源的图标插件,可以用它开发出功能强大的图表系统。下面着重讲解在Asp.net中如何使用这个插件做出功能强大的图表应用。关于Flot,可以在这里查看现有的例子(或者是这里的例子),可以在这里查看现有的API。Flot的官方页面在这里,在里面我们可以下载需要使用的插件。下面一段话是摘自Flot官网,藉此对Flot有个初步的映像:Flot is a pure ..._flot折线图

对ionic build android中所遇到的问题的总结_ionic exception in thread "main" java.lang.runtime-程序员宅基地

文章浏览阅读7.1k次。因为实习岗位要求,我不得不加入前端大军。公司开发的项目采用ionic。好在据说这样的Hybrid框架简单易学,只好半路出家,扛起枪就上了。但是在最初环境搭建和项目build的过程中遇到了不少问题。前前后后花了一周的时间才处理好。现在说说我遇到的问题,希望能够给大家一些参考,少走一些弯路。Exception in thread "main" java.lang.RuntimeExcept_ionic exception in thread "main" java.lang.runtimeexception: gradle distribu

推荐文章

热门文章

相关标签