CN102800322B - 一种噪声功率谱估计与语音活动性检测方法 - Google Patents

Abstract

本发明涉及一种噪声功率谱估计与语音活动性检测方法，该方法基于一阶回归的序贯隐马尔可夫模型SHMM描述语音在每一频率分量上的时序相关性，最终推演出语音的在这个频率子带上的出现概率和噪声的功率谱信息，1)对于语音信号在每一个频率分量上提取对数幅度谱包络，并建立一个对应的二元隐马尔可夫模型，每一个状态均由高斯分布表示；2)对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，采用极大似然估计算法建立一个初始化的模型；3)在得到初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新每一频带的HMM模型，依次递推得到噪声值和语音信号的出现概率。

Description

一种噪声功率谱估计与语音活动性检测方法

技术领域

本发明涉及语音信号处理的技术领域，具体的说，本发明涉及一种基于序贯隐马尔可夫模型的噪声谱估计与语音活动性检测方法。其中，语音活动性检测是在时间维度上判断语音出现与否的算法，它既能以“是”或“否”的形式回答存在性，也可以用语音出现概率来描述语音的存在性。

背景技术

语音活动性检测和噪声功率谱估计是降噪算法必不可少的组成部分，它们的性能直接影响降噪算法的性能，特别是在恶劣的噪声环境下，它们间接影响了语音处理系统(如语音识别、话者识别以及语言识别器)的性能。

大多数的语音应用系统不得不面对环境噪声干扰。前人提出了很多方法去除噪声对语音系统的干扰，几乎所有的方法都依赖于语音活动性检测和噪声功率谱估计。这两个模块存在着紧密的联系，它们的准确性直接影响系统的整体抗噪性能。传统的估计方法尽管性能良好，但仍然有两个地方值得提高：

1、充分利用连续语音/非语音信号在某一频率分量上的时序相关性，现有的算法对于时间相关性的利用不够充分，它们往往采用比较简单的一阶回归平滑器对幅度谱包络进行平滑，而且平滑器的平滑因子是固定的。而语音信号本身就是一个分段平稳的信号，它的统计特征包括时序相关性，都是随着时间在不断变化的，一个固定的模型无法反映这种时变特性。如果我们能够采用自适应的模型对时序相关性进行建模，那么算法的性能无疑将会得到提升。这一方法在以往的文献资料中没有提及。

2、传统的序贯HMM的参数自适应采用高阶的回归平均方式，当前的HMM参数集依赖于上一时刻的模型、当前的观察值以及过去多个时刻的观察值，这种参数回归的方式计算量巨大。如果能够在精度损失不大的情况下把这种高阶回归简化为一阶回归，那么，算法的计算效率将大为提高。基于一阶回归的序贯HMM算法，在以往的文献资料中也没有提及。

另外，传统的解决方案是基于半监督学习的方式。在起始阶段，一股的系统需要做出“噪声起始”的假设，即假设句子的开头总是存在一段非语音信号。这段非语音信号可以理解为人工标注的背景噪声样本，从这些标注样本中建立噪声的初始化模型，这是一种监督学习方法。其缺陷是：在一些应用中这个假设难以得到满足，比如当句子以语音信号开头，那么将导致噪声模型的初始化失败，继而使得语音检测和噪声功率谱估计都不准确。这种初始化的方法在中国申请号201010178166.4的专利中已经公开。

发明内容

本发明目的的目的在于，为提供一种基于序贯隐马尔可夫模型的噪声谱估计与语音活动性检测方法，该方法利用隐马尔可夫模型对语音信号在某频率分量上存在的时序相关性进行建模，在某个频率分量上的对数功率谱包络可以看作一个马尔可夫链，这个链条在语音“出现”和“不出现”两种状态间跳转，对于每个状态，采用一个高斯分布描述其功率谱的分布，再根据HMM的前向因子，可以导出语音在某个时频点的出现概率。

为实现上述发明目的，本发明提供了一种噪声功率谱估计与语音活动性检测方法，该方法基于一阶回归的序贯隐马尔可夫模型SHMM描述语音在每一频率分量上的时序相关性，并采用增量学习的方式来逐步更新SHMM，最终，推演出语音的在这个频率子带上的出现概率和噪声的功率谱信息，以准确反映语音的时序统计特征，该方法包括下列步骤：

1)对于语音信号在每一个频率分量上提取对数幅度谱包络，并建立一个对应的二元隐马尔可夫模型，其中，一个分量表示语音能量的分布，另一个分量是噪声能量的分布，每一个状态均由高斯分布表示；

2)对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，采用极大似然估计算法建立一个初始化的模型；

3)在得到初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新每一频带的HMM模型，依次递推得到噪声值和语音信号的出现概率。

该方法的具体步骤包括：

1)对于语音信号在每一个频率分量上提取对数幅度谱包络，对于一个频率分量上的对数幅度谱时间序列x_l＝{x₁，x₂，...，x_l}，建立一个隐马尔可夫模型s_l＝{s₁，s₂，...，s_l}，s_t∈{0，1}是其对应的状态序列，1表示语音出现状态，0表示噪声出现状态，λ_l表示从序列x_l中获取的模型参数估值，那么，对于给定一个参数集λ_l，对应的观察值序列x_l的概率密度函数可以表示为：

$p\left(x_l\right|{\mathrm{\λ}}_l)=\underset{s_l}{\mathrm{\Σ}}p\left(s_l\right|{\mathrm{\λ}}_l)p\left(x_l\right|{\mathrm{\λ}}_l,s_l);$

其中，p(s_l|λ_l)表示状态序列s_l出现的先验概率，高斯分量表示为：

$p\left(s_l\right|{\mathrm{\λ}}_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}{a}_{s_{t-1},s_t};$

这里的

表示状态转移概率，

表示初始状态概率，p(x_l|λ_l，s_l)表示给定状态s_l和参数集λ_l的情况下观察值序列x_l的似然度：

$p\left(x_l\right|{\mathrm{\λ}}_l,s_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}b\left(x_t\right|s_t,{\mathrm{\λ}}_l);$

其中，

$b\left(x_t\right|s_t,{\mathrm{\λ}}_l)=\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\κ}}_{s_t,l}}}\exp \{-\frac{1}{2}{(x_t-{\mathrm{\μ}}_{s_t,l})}^2/{\mathrm{\κ}}_{s_t,l}\};$

变化；

这个模型中的μ_0，l就是我们想要估计的噪声，同时，我们可以推导出语音信号在第l帧的某个频点上出现的概率为

2)对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，代入步骤1)的HMM模型对每一个频点上初始化一个隐马尔可夫模型λ_M，下标M表示初始化的时间窗口长度，l≥M；

3)在得到初始化的模型λ_M之后，从第M+1帧开始，HMM模型采用增量学习的方法，逐帧更新SHMM模型，依次递推得到λ_l；并得出噪声值μ_0，l和语音信号在第l帧的某个频点上的出现概率。

作为上述技术方案的一种改进，所述的步骤1)中提取一帧幅度谱的步骤包括：

首先，对这一帧的数字化声音信号做预处理，设每帧长度为F点，先补零到N点，N≥F，N＝2^j，j为整数，且j≥8，进行N点离散傅里叶变换，得到离散谱

其中，y_l，n表示缓存中第l帧的第n个采样点，Y_l，k表示缓存中第i帧的第k个傅里叶变换值(k＝0，1，...，N-1)；那么，它的幅度值可以计算为式中，b(r)为加窗函数。所述的预处理包括加窗、或/和预加重；所述的加窗函数采用汉宁窗或哈宁窗。

作为上述技术方案的一种改进，所述的步骤2)中HMM的初始化，在某个频点上具体初始化的步骤包括：

步骤201)：通过聚类的方法将M个样本分成两类：和

其中，M₀+M₁＝M，均值较大的一类用上标(1)表示，另一类用上标(0)表示；所述的步骤201)中的聚类的方法采用LBG非监督聚类或者模糊聚类方法；

两类的均值为

能量较小的一类的均值为其中，

两类的方差分别为： ${\stackrel{\‾}{\mathrm{\κ}}}_{0,M}=\frac{1}{M_0}{\mathrm{\Σ}}_{j=1}^{M_0}{(x_{i_j}-{\stackrel{\‾}{\mathrm{\μ}}}_{0,M})}^2,{\stackrel{\‾}{\mathrm{\κ}}}_{1,M}=\frac{1}{M_1}{\mathrm{\Σ}}_{j=1}^{M_1}{(x_{i_j}-{\stackrel{\‾}{\mathrm{\μ}}}_{1,M})}^2;$

两类的初始化权重系数为： ${\stackrel{\‾}{a}}_{00,M}={\stackrel{\‾}{a}}_{01,M}={\stackrel{\‾}{a}}_{11,M}={\stackrel{\‾}{a}}_{10,M}=0.5;$

计算新模型的似然度，

并开始迭代运算；在以下的迭代过程中，老的模型参数集合表示为λ′_M，新的模型参数为：开始迭代以前，

L′设置为一个很大的负数，初始化前向因子，

初始化后向因子，

步骤202)：计算前向因子： ${\stackrel{\‾}{F}}_l\left(z\right)=\underset{y}{\mathrm{\Σ}}{\stackrel{\‾}{F}}_{l-1}\left(z\right){\stackrel{\‾}{a}}_{y,z,M}b\left(x_l\right|y,{\stackrel{\‾}{\mathrm{\λ}}}_M),z,y\∈\left\{\mathrm{0,1}\right\};$

步骤203)：计算后向因子： ${\stackrel{\‾}{B}}_l\left(z\right)=\underset{y}{\mathrm{\Σ}}{\stackrel{\‾}{B}}_{l+1}\left(y\right){\stackrel{\‾}{a}}_{z,y,M}b\left(x_{l+1}\right|y,{\stackrel{\‾}{\mathrm{\λ}}}_M),z,y\∈\left\{\mathrm{0,1}\right\};$

步骤204)：计算噪声和语音出现概率： $p\left(z\right|x_l,{\stackrel{\‾}{\mathrm{\λ}}}_M)=\frac{{\stackrel{\‾}{F}}_l\left(z\right){\stackrel{\‾}{B}}_l\left(z\right)}{{\mathrm{\Σ}}_z{\stackrel{\‾}{F}}_l\left(z\right){\stackrel{\‾}{B}}_l\left(z\right)},z\∈\left\{\mathrm{0,1}\right\};$

步骤205)：如果

停止跌代，其中ζ是接近于零但大于零的小数；

步骤206)：计算转移概率：

$p(s_{l-1}=y,s_l=z|x_l,{\stackrel{\‾}{\mathrm{\λ}}}_M)=\frac{{\stackrel{\‾}{F}}_{l-1}\left(y\right){\stackrel{\‾}{B}}_l\left(z\right){\stackrel{\‾}{a}}_{\mathrm{yz},M}b\left(x_l\right|z,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{{\mathrm{\Σ}}_z{\stackrel{\‾}{F}}_{l-1}\left(y\right){\stackrel{\‾}{B}}_l\left(z\right){\stackrel{\‾}{a}}_{\mathrm{yz},M}b\left(x_l\right|z,{\stackrel{\‾}{\mathrm{\λ}}}_M)};$

步骤207)：计算新的初始化概率 ${\mathrm{\π}}_z^{\′}=p(s_1=z|x_1,{\stackrel{\‾}{\mathrm{\λ}}}_M);$

步骤208)：计算新的均值 ${\mathrm{\μ}}_{z,M}^{\′}={\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)x_t}{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}}_{{\stackrel{\‾}{\mathrm{\λ}}}_M{\stackrel{\‾}{\mathrm{\λ}}}_M{\stackrel{\‾}{\mathrm{\λ}}}_M};$

步骤209)：对新的均值进行约束：μ′_1，M＝max{μ′_0，M，μ′_0，M+δ}，其中，δ是一个常数，取值范围在0至100之间；

步骤210)：计算新方差 ${\mathrm{\κ}}_{z,M}^{\′}=\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M){(x_t-{\stackrel{\‾}{\mathrm{\μ}}}_{z,M})}^2}{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)};$

步骤211)：对新方差进行约束，κ′_1，M＝max{κ′_0，M，κ′_1，M}；

步骤212)：计算新的转移概率， $a_{\mathrm{yz},M}^{\′}=\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{{\mathrm{\Σ}}_{t=1}^M{\mathrm{\Σ}}_{z}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)};$

步骤213)：计算新模型的似然度 $\stackrel{\‾}{L}=\log \left(p\left(x_M\right|{\mathrm{\λ}}_M^{\′})\right);$

步骤214)：如果满足条件

终止迭代，其中，ε是一个很小的数字，如果

迭代跳转至步骤202)。

在上述的HMM参数的建模过程中，分别对HMM的均值、权重、方差以及转移概率进行约束。需要指出的是：在初始化过程中，这里的转移概率所发挥的作用与专利201010178166.4中的权重系数相当，由于201010178166.4中的权重系数在初始化过程中用作分母项，所以它必须在初始化过程中进行约束。而本专利中的转移概率不存在此问题。

作为上述技术方案的一种改进，所述的步骤3)中的HMM的序贯更新是在建立初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新HMM模型，其迭代过程可以表述为：在每一个频点上，已知λ_l和当前的观察值x_l，推断λ_l+1；对于第l+1帧进行傅立叶变换，得到Y_l+1，k，其中，0≤k＜N；在每一个频点上，计算幅度值对于每一个频点，在l+1帧的参数更新步骤如下：

步骤301)：计算前向因子， $F_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\underset{y}{\mathrm{\Σ}}{F}_{l|{\mathrm{\λ}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l),$ z∈{0，1}；

步骤302)：计算语音和噪声出现概率， ${\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(z\right)}{{\mathrm{\Σ}}_z{F}_{l+1|{\mathrm{\λ}}_l}\left(z\right)},z\∈\left\{\mathrm{0,1}\right\};$

步骤303)：计算条件转移概率，

${\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}{{\mathrm{\Σ}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)};$

步骤304)：计算平均噪声语音出现概率， ${\widetilde{\mathrm{\γ}}}_{l+1}\left(z\right)=\mathrm{\α}{\widetilde{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right);$

步骤305)：计算时间依赖平滑因子， ${\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)=\frac{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)}{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)};$

步骤306)：计算状态均值， ${\mathrm{\μ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\μ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]x_{l+1};$

步骤307)：对新的状态均值进行约束：μ_1，l+1＝max{μ_1，l+1，μ_0，l+1+δ}，l≥M；

步骤308)：计算新的状态方差， ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2;$

步骤309)：对新状态方差进行约束，κ_1，l+1＝max{κ_0，l+1，κ_1，l+1}，l≥M；

步骤310)：计算平均转移概率， ${\widetilde{\mathrm{\ξ}}}_{l+1}(y,z)=\mathrm{\α}{\widetilde{\mathrm{\ξ}}}_l(y,z)+(1-\mathrm{\α}){\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z);$

步骤311)：计算状态概率， $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,1-z)};$

步骤312)：对新的转移概率进行约束，a_01，l＝max{a_01，l，η}，a_00，l＝1-a_01，l，a_10，l＝max{a_10，l，η}，a_11，l＝1-a_10，l，l≥M；

从以上子步骤中，获得了λ_l+1中所有的参数，从而得到了相应的语音出现概率γ_l+1|λl(1)和噪声信号的功率谱估值μ_0，l+1。

所述的步骤3)中的HMM模型采用的增量学习方法包括：递推权重系数、递推均值和递椎方差；

其中，所述的递推均值：式中，

是一个依赖于语音出现概率的平滑因子，小于1但接近于1；

所述的递推方差： ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2;$

所述的递推转移概率：

或者a_yz，l+1＝βa_yz，l+(1-β)ξ_l+1|λl(y，z)；式中，β是一个小于1但接近于1的平滑因子，例如β＝0.99。

所述的基于一阶回归的序贯隐马尔可夫模型的参数递推方法为：

计算HMM的前向因子： $F_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\underset{y}{\mathrm{\Σ}}{F}_{l|{\mathrm{\λ}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l),z\∈\left\{\mathrm{0,1}\right\};$

计算语音和噪声出现概率， ${\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(z\right)}{{\mathrm{\Σ}}_z{F}_{l+1|{\mathrm{\λ}}_l}\left(z\right)},z\∈\left\{\mathrm{0,1}\right\};$

计算条件转移概率， ${\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}{{\mathrm{\Σ}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)};$

计算平均噪声语音出现概率， ${\widetilde{\mathrm{\γ}}}_{l+1}\left(z\right)=\mathrm{\α}{\widetilde{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right);$

计算时间依赖平滑因子， ${\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)=\frac{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)}{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)};$

计算均值， ${\mathrm{\μ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\μ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]x_{l+1};$

计算新的方差， ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2;$

计算平均转移概率， ${\widetilde{\mathrm{\ξ}}}_{l+1}(y,z)=\mathrm{\α}{\widetilde{\mathrm{\ξ}}}_l(y,z)+(1-\mathrm{\α}){\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z);$

计算转移概率， $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l(y,z)}}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,1-z)}.$

上述技术方案中，采用的保障二元统计模型稳健运行的约束机制包括：

1)在初始阶段，当语音的平均出现概率小于某个固定门限ζ时，

估计算法停止跌代。该约束在步骤205中实施。

2)为防止隐马尔可夫模型的状态跳转停止，对模型的转移状态进行约束。a_01，l＝max{a_01，l，η}，a_00，l＝1-a_01，l，a_10，l＝max{a_10，l，η}，a_11，l＝1-a_10，l，l≥M。该约束在步骤312中实施。

3)在追踪过程中，对均值的约束，μ_1，l+1＝max{μ_1，l+1，μ_0，l+1+δ}，l≥M。该约束在步骤307中实施。

4)对方差的约束，κ_1，l+1＝max{κ_0，l+1，κ_1，l+1}，l≥M。该约束在步骤309中实施。

本发明涉及一种基于序贯隐马尔可夫模型(Sequential Hidden Markov Model，SHMM)的噪声功率谱估计与语音活动性检测方法，包括下列步骤：1)对于语音信号在每一个频率分量上的对数幅度特征，建立一个SHMM模型，2)对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，代入步骤1)的SHMM模型进行初始化，得到初始化的模型λ_M；3)在得到初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新SHMM模型。模型中的噪声状态平均值就是当前的噪声估计值，估计过程中的语音出现概率表示语音在时频域的活动性。递推的方法是：根据当前的观察值x_l和上一时刻的模型参数集λ_l-1，估计当前时刻的模型参数集λ_l。由此，依次得到每一时刻某个频率分量上的噪声功率谱和语音出现的概率。本发明是谱估计与语音活动性检测的紧耦合解决方案，能够增强语音应用系统对噪声环境的适应性；本发明不依赖于“噪声起始”假设；并且，本发明还能够提供语音活动性在时频二维空间上的描述。本专利是在已申请的专利号201010178166.4的专利基础之上发展而来的，由于采用了更加模型精确的建模方法，本专利的性能优于201010178166.4，但计算复杂度高于201010178166.4。

与现有技术相比，本发明具有下列技术效果：

基于语音信号在某频率分量上存在时序相关性，本发明利用隐马尔可夫模型对这种时序相关性进行建模，在某个频率分量上的对数功率谱包络可以看作一个马尔可夫链，这个链条在语音“出现”和“不出现”两种状态间跳转，对于每个状态，采用一个高斯分布描述其功率谱的分布。为了简化计算，本发明还提出了一阶回归的序贯HMM追踪方法，它的参数随着输入信号不断变化。其中HMM的语音“不出现”状态的均值就是噪声功率谱的估计值，根据HMM的前向因子，可以导出语音在某个时频点的出现概率。

本发明是一种语音活动性检测和噪声功率谱估计紧耦合的方案，能够增强语音应用系统对噪声环境的适应性；并且，本发明还能够提供语音活动性在时频二维空间上的描述，有利于对噪声进行进一步的精细化处理。

附图说明

图1本发明的噪声谱估计与语音活动性检测方法的的流程图；

图2为实例比较了本发明的SHMM噪声估计算法、经典的最小统计算法(MS)、最小控制回归平均算法(MCRA)、以及它的提高版本IMCRA效果比对图。

具体实施方式

本发明提出了一种基于序贯隐马尔可夫模型的噪声功率谱估计与语音活动性检测方法。

如图1所示，包括下列步骤：

1)对于语音信号在每一个频点上的对数幅度特征，建立一个HMM模型，数学表达式如下：

$p\left(x_l\right|{\mathrm{\λ}}_l)=\underset{s_l}{\mathrm{\Σ}}\underset{t=1}{\overset{l}{\mathrm{\Π}}}{a}_{s_{t-1},s_t}\underset{t=1}{\overset{l}{\mathrm{\Π}}}b\left(x_t\right|s_t,{\mathrm{\λ}}_l)$

这里的

表示状态转移概率，

表示初始状态概率，其中高斯分量表示为：

$b\left(x_t\right|s_t,{\mathrm{\λ}}_l)=\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\κ}}_{z,l}}}\exp \{-\frac{1}{2}{(x_t-{\mathrm{\μ}}_{z,l})}^2/{\mathrm{\κ}}_{z,l}\}$

其中，x_l表示第l帧的某个频点上的对数幅度谱，z＝0表示语音不出现状态，z＝1表示出现状态。μ_z，k和κ_z，k分别表示均值和方差，参数集λ_l＝{μ_0，l，μ_1，l，κ_1，l，κ_0，l，a_01，l，a_10，l，a_00，l，a_11，l，π₀，π₁}。

2)对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，代入步骤1)的GMM模型进行初始化，得到初始化的模型λ_0，k；初始化过程采用约束EM算法；M表示初始化窗口的长度。

3)在得到初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新HMM模型，依次递推得到λ_l。并得出噪声值μ_0，l和语音信号在第l帧的某个频点上的出现概率。

其中i＝1，2，3，……。

其中，所述GMM的增量学习方法包括递推权重系数、递推均值和递推方差；

其中前向因子递推方法为： $F_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\underset{y}{\mathrm{\Σ}}{F}_{l|{\mathrm{\λ}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l),$ z∈{0，1}。

语音和噪声出现概率递推方法为： ${\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(z\right)}{{\mathrm{\Σ}}_z{F}_{l+1|{\mathrm{\λ}}_l}\left(z\right)},z\∈\left\{\mathrm{0,1}\right\}$

条件转移概率递推方法为： ${\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}{{\mathrm{\Σ}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}$

平均噪声语音出现概率递推方法为： ${\widetilde{\mathrm{\γ}}}_{l+1}\left(z\right)=\mathrm{\α}{\widetilde{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)$

时间依赖平滑因子递推方法为： ${\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)=\frac{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)}{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)}$

状态均值递推方法为： ${\mathrm{\μ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\μ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]x_{l+1}$

状态方差递推方法为： ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2$

平均转移概率递推方法为： ${\widetilde{\mathrm{\ξ}}}_{l+1}(y,z)=\mathrm{\α}{\widetilde{\mathrm{\ξ}}}_l(y,z)+(1-\mathrm{\α}){\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)$

状态概率递推方法为： $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,1-z)}$

序贯隐马尔可夫模型的最大特点是可以在线追踪某频率分量上语音出现的时序相关性，它将某频率分量上的功率普包络看作一个在语音与非语音状态之间切换的马尔可夫链。它采用非监督的方式构建初始化模型。具体说来，它具有以下几个特点：

●由于采用HMM，可以采用Viterbi译码的方式，给出一个时间序列上，语音是否出现的最优估计。

●在初始化阶段，不依赖于噪声起始假设，所以该发明的应用范围比一股的解决方案应用更加宽泛。

●语音活动性是一个“时间——频率”的二维信息，而其它的语音活动性检测算法仅仅描述了语音在时间维度上的存在性。

在一个实施例中，非监督学习框架的载体是二元隐马尔可夫模型(Hidden MarkovModel，缩写为HMM)。其中一个分量表示语音能量的分布，另一个分量是噪声能量的分布。在每一个频率分量上提取对数幅度谱包络，并建立一个对应的HMM。首先采用EM算法初始化HMM，然后采用增量学习的方式来逐步更新HMM。根据HMM模型，分别推演出语音的在这个子带上的出现概率和噪声的功率谱信息。在HMM参数的建模过程中，分别对HMM的均值、权重、方差以及转移概率进行约束。其中，对于HMM参数的序贯估计方法，具体包括递推权重系数、递推均值和递推方差以及递推的计算。

1)递推均值：

其中

是一个依赖于语

音出现概率的平滑因子，小于1但接近于1。

2)递推方差， ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2;$

3)递推转移概率，

或者

其中β是一个小于1但接近于1的平滑因子，例如β＝0.99。

下面结合一个优选实施例对本发明做进一步地描述。

本发明的原理如下：

噪声估计过程在每一个频率分量上并行运行，所以，以下的描述中省略掉频率分量索引k。对于语音信号在每一个频率分量上的对数幅度谱时间序列x_l＝{x₁，x₂，...，x_l}，建立一个隐马尔可夫模型，s_l＝{s₁，s₂，...，s_l}，s_t∈{0，1}是其对应的状态序列，1表示语音出现状态，0表示噪声出现状态，λ_l表示从序列x_l中获取的模型参数估值，那么对于给定一个参数集λ_l，对应的观察值序列x_l的概率密度函数可以表示为：

$p\left(x_l\right|{\mathrm{\λ}}_l)=\underset{s_l}{\mathrm{\Σ}}p\left(s_l\right|{\mathrm{\λ}}_l)p\left(x_l\right|{\mathrm{\λ}}_l,s_l)$

其中，p(s_l|λ_l)表示状态序列s_l出现的先验概率高斯分量表示为：

$p\left(s_l\right|{\mathrm{\λ}}_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}{a}_{s_{t-1},s_t}$

这里的

表示状态转移

表示初始状态概率，p(x_l|λ_l，s_l)表示给定状态s_l和参数集λ_l的情况下观察值序列x_l的似然度：

$p\left(x_l\right|{\mathrm{\λ}}_l,s_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}b\left(x_t\right|s_t,{\mathrm{\λ}}_l)$

其中

$b\left(x_t\right|s_t,{\mathrm{\λ}}_l)=\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\κ}}_{z,l}}}\exp \{-\frac{1}{2}{(x_t-{\mathrm{\μ}}_{z,l})}^2/{\mathrm{\κ}}_{z,l}\}$

这里κ_z，l表示高斯分布方差，μ_z，l表示均值，s_l＝z，λ_l＝{μ_0，l，μ_1，l，κ_1，l，κ_0，l，a_01，l，a_10，l，a_00，l，a_11，l，π₀，π₁}，参数集中的初始概率π_z不随着时间变化。

这个模型中的μ_0，l就是我们想要估计的噪声。同时，我们可以推导出语音信号在第l帧的某个频点上出现的概率为γ_t|λl(z)＝p(s_t＝z|x_t，λ_l)。

基于上述原理，根据本发明的一个实施例，所述噪声功率谱估计和语音活动性检测方法包括下列步骤：

步骤100：设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的幅度谱。提取一帧幅度谱的方法如下：

首先对这一帧的数字化声音信号做预处理(根据系统实际情况，可以包括加窗、预加重等)，设每帧长度为F点，先补零到N点(其中N≥F，N＝2^j，j为整数且j≥8)，进行N点离散傅里叶变换，得到离散谱

其中y_l，n表示缓存中第l帧的第n个采样点，Y_l，k表示缓存中第i帧的第k个傅里叶变换值(k＝0，1，...，N-1)。那么，它的幅度值可以计算为

b(r)为加窗函数(如汉宁窗，哈宁窗等等)，注意这里的k在下面的描述中被省略。

步骤200：HMM的初始化。在每一个频点上初始化一个隐马尔可夫模型λ_M，其中下标M表示初始化的时间窗口长度，初始化过程采用约束EM算法，在某个频点上，具体初始化步骤如下：

步骤201：通过聚类的方法(例如LBG非监督聚类，或者模糊聚类等等)将M个样本分成两类：

和

其中M₀+M₁＝M，均值较大的一类用上标(1)表示，另一类用上标(0)表示。两类的均值为

能量较小的一类的均值为其中

两类的方差分别为：两类的初始化权重系数：

计算新模型的似然度，在以下的迭代过程中，老的模型参数集合表示为λ′_M，新的模型参数为：

在开始迭代以前，

L′设置为一个很大的负数，例如L′_k＝-10000。初始化前向因子，

初始化后向因子，

以下开始迭代运算。

步骤202：计算前向因子： ${\stackrel{\‾}{F}}_t\left(z\right)=\underset{y}{\mathrm{\Σ}}{\stackrel{\‾}{F}}_{t-1}\left(z\right){\stackrel{\‾}{a}}_{\mathrm{yz},M}b\left(x_t\right|y,{\stackrel{\‾}{\mathrm{\λ}}}_M),z\∈\left\{\mathrm{0,1}\right\}.$

步骤203：计算后向因子： ${\stackrel{\‾}{B}}_t\left(z\right)=\underset{y}{\mathrm{\Σ}}{\stackrel{\‾}{B}}_{t+1}\left(y\right){\stackrel{\‾}{a}}_{zy,M}b\left(x_{l+1}\right|y,{\stackrel{\‾}{\mathrm{\λ}}}_M),z\∈\left\{\mathrm{0,1}\right\}.$

步骤204：计算噪声和语音出现概率： $p\left(z\right|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)=\frac{{\stackrel{\‾}{F}}_t\left(z\right){\stackrel{\‾}{B}}_t\left(z\right)}{{\mathrm{\Σ}}_z{\stackrel{\‾}{F}}_t\left(z\right){\stackrel{\‾}{B}}_t\left(z\right)},z\∈\left\{\mathrm{0,1}\right\}$

步骤205：如果停止跌代。其中ζ是接近于零但大于零的小数。

步骤206：计算转移概率： $p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)=\frac{{\stackrel{\‾}{F}}_{t-1}\left(y\right){\stackrel{\‾}{B}}_t\left(z\right){\stackrel{\‾}{a}}_{\mathrm{yz},M}b\left(x_t\right|z,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{{\mathrm{\Σ}}_z{\stackrel{\‾}{F}}_{t-1}\left(y\right){\stackrel{\‾}{B}}_t\left(z\right){\stackrel{\‾}{a}}_{\mathrm{yz},M}b\left(x_t\right|z,{\stackrel{\‾}{\mathrm{\λ}}}_M)}.$

步骤207：计算新的初始化概率 ${\mathrm{\π}}_z^{\′}=p(s_1=z|x_1,{\stackrel{\‾}{\mathrm{\λ}}}_M)$

步骤208：计算新的均值 ${\mathrm{\μ}}_{z,M}^{\′}=\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{\stackrel{\‾}{\mathrm{\λ}}}_M{\stackrel{\‾}{\mathrm{\λ}}}_M{\stackrel{\‾}{\mathrm{\λ}}}_M$

步骤209：对新的均值进行约束：μ′_1，k＝max{μ′_0，k，μ′_0，k+δ}，其中δ是一个常数，取值范围在0至100之间。

步骤210：计算新方差 ${\mathrm{\κ}}_{z,M}^{\′}=\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M){(x_t-{\stackrel{\‾}{\mathrm{\μ}}}_{z,M})}^2}{{\mathrm{\Σ}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}$

步骤211：对新方差进行约束，κ′_1，M＝max{κ′_0，M，κ′_1，M}

步骤212：计算新的转移概率， $a_{\mathrm{yz},M}^{\′}=\frac{{\mathrm{\Σ}}_{t=1}^{M}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}{{\mathrm{\Σ}}_{t=1}^M{\mathrm{\Σ}}_{z}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\‾}{\mathrm{\λ}}}_M)}$

步骤213：计算新模型的似然度 $\stackrel{\‾}{L}=\log \left(p\left(x_M\right|{\mathrm{\λ}}_M^{\′})\right),$

步骤214：如果满足条件

终止迭代，其中ε是一个很小的数字，例如ε＝0.1。如果

迭代跳转至“步骤202”。

步骤300：HMM的序贯更新。在建立初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新HMM模型。迭代过程可以表述为：在每一个频点上，已知λ_l和当前的观察值x_l，推断λ_l+1。对于第l+1帧进行傅立叶变换，得到T_l+1，k，其中0≤k＜N。在每一个频点上，计算幅度值对于每一个频点，在l+1帧的参数更新步骤如下：

步骤301：计算前向因子， $F_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\underset{y}{\mathrm{\Σ}}{F}_{l|{\mathrm{\λ}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l),z\∈\left\{\mathrm{0,1}\right\}.$

步骤302：计算语音和噪声出现概率， ${\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(z\right)}{{\mathrm{\Σ}}_z{F}_{l+1|{\mathrm{\λ}}_l}\left(z\right)},z\∈\left\{\mathrm{0,1}\right\}$

步骤303：计算条件转移概率， ${\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}{{\mathrm{\Σ}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}$

步骤304：计算平均噪声语音出现概率， ${\widetilde{\mathrm{\γ}}}_{l+1}\left(z\right)=\mathrm{\α}{\widetilde{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)$

步骤305：计算时间依赖平滑因子， ${\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)=\frac{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)}{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)}$

步骤306：计算均值， ${\mathrm{\μ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\μ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]x_{l+1}$

步骤307：对新的均值进行约束：μ_1，l+1＝max{μ_1，l+1，μ_0，l+1+δ}。

步骤308：计算新的方差， ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2$

步骤309：对新方差进行约束，κ_1，l+1＝max{κ_0，l+1，κ_1，l+1}

步骤310：计算平均转移概率， ${\widetilde{\mathrm{\ξ}}}_{l+1}(y,z)=\mathrm{\α}{\widetilde{\mathrm{\ξ}}}_l(y,z)+(1-\mathrm{\α}){\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)$

步骤311：计算转移概率， $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,1-z)}$

步骤312：对新的转移概率进行约束，a_01，l＝max{a_01，l，η}，a_00，l＝1-a_01，l，a_10，l＝max{a_10，l，η}，a_11，l＝1-a_10，l。

从以上子步骤中，我们获得了λ_l+1中所有的参数，从而得到了相应的语音出现概率和噪声信号的功率谱估值μ_0，l+1。

基于上述实施例的算法，对噪声功率谱估计的性能进行评价，采用TIMIT数据库中的男女话者语音数据各8个句子，和NOISEX92噪声数据库中的高斯白噪声、F16战斗机座仓噪声、以及babble噪声按照0、5、10dB等信噪比混合。第一种评价指标线性分段误差定义如下：

${\mathrm{\ϵ}}_n=\frac{1}{L}\underset{l=1}{\overset{L}{\mathrm{\Σ}}}\{10{\log }_{10}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{[D_{k,l}-{\hat{D}}_{k,l}]}^2/\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{D}_{k,l}^2\}$

其中D(k，l)表示实际的噪声幅度谱，

表示估计的噪声幅度谱，注意误差值越小，表示估计值越接近于真实值，估计越准确。第二种评价指标对数分段误差定义如下：

${\mathrm{\ϵ}}_r=\frac{1}{L}\underset{l=1}{\overset{L}{\mathrm{\Σ}}}{\left\{\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\left[20{\log }_{10}\right|D_{k,l}|-20{\log }_{10}|{\hat{D}}_{k,l}\left|\right]}^2\right\}}^{1/2}.$

同理，误差值越小，表示估计噪声越准确。

算法分别与当前主流的三种噪声功率谱估计算法进行了比较，其中MS表示最小统计算法，MCRA表示最小控制的回归平均算法，IMCRA表示提高版本的最小控制回归平均算法，SHMM为本发明的算法。表1表示出了线谱误差SegError的结果。

表1各种环境下的噪声估计线性分段误差

表2各种环境下的噪声估计对数分段误差

从上表可以看出，本发明提出的算法对于目前主流的三种算法均具有明显的优势。

另外，图2通过一个实例比较了SHMM噪声估计算法、经典的最小统计算法(MS)、最小控制回归平均算法(MCRA)、以及它的提高版本IMCRA。本例中，带噪语音信号在3.75秒的位置，信噪比从10dB突然下降到4dB，在13.1秒的位置，又从4dB上升到10dB。(a)为带噪语音在某一子带上的功率谱；(b)为MS算法估计的噪声功率谱包络，其中虚线表示真实的噪声功率谱包络；(c)表示MCRA算法的估计结果；(d)表示IMCRA算法的估计结果；(e)表示本算法的估计结果。从图中可以看出，其它三个算法对于噪声的突然上升的反应较慢，不能快速跟上噪声的突升。而SHMM算法表现出的性能较好。

在过去的几十年中，人们发明了各种各样的算法，用于估计语音活动性和噪声功率谱。语音信号在频域的时序相关性是其中的一个重要线索，由于语音信号是非稳定信号，这种时序相关性也是随着时间变化的。然而，过去的算法对这种时序相关性没有引起足够的重视，只是简单地加以利用，没有采用自适应的方法描述时变相关性。本发明采用序贯隐马尔可夫模型SHMM描述语音在每一频率分量上的时序相关性，这种序贯估计方法是建立在一阶回归的基础之上，时序相关性和它的参数集一起，随着输入信号的变化而变化。这种统计模型准确反映了语音的时序统计特征，因此，本发明提出的估计算法在性能优于目前主流的算法(例如最小统计，最小控制的回归平均)。

此外，以往SHMM模型都是在建立在高阶回归的基础之上，本发明中提出的一阶回归SHMM相对于高阶回归的SHMM，大大节约了计算量和存储空间。这是本发明的另一个创新之处。

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (9)

1.一种噪声功率谱估计与语音活动性检测方法，该方法基于一阶回归的序贯隐马尔可夫模型SHMM描述语音在每一频率分量上的时序相关性，并采用增量学习的方式来逐步更新SHMM，最终，推演出语音的在这个频率子带上的出现概率和噪声的功率谱信息，以准确反映语音的时序统计特征，该方法包括下列步骤：

1）对于语音信号在每一个频率分量上提取对数幅度谱包络，并建立一个对应的二元隐马尔可夫模型，其中，一个分量表示语音能量的分布，另一个分量是噪声能量的分布，每一个状态均由高斯分布表示；

2）对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，采用极大似然估计算法建立一个初始化的模型；

3）在得到初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新每一频带的HMM模型，依次递推得到噪声值和语音信号的出现概率；

所述的步骤3）中的HMM的序贯更新是在建立初始化的模型λ_M之后，从第M+1帧开始，采用增量学习的方法，逐帧更新HMM模型，其迭代过程可以表述为：在每一个频点上，已知λ_l和当前的观察值x_l，推断λ_l+1；对于第l+1帧进行傅立叶变换，得到Y_l+1,k，其中，0≤k＜N；在每一个频点上，计算幅度值对于每一个频点，在l+1帧的参数更新步骤如下：

步骤301）：计算前向因子， $F_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\underset{y}{\mathrm{\Σ}}{F}_{l|{\mathrm{\λ}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l),$ z∈{0,1}；

步骤302）：计算语音和噪声出现概率， ${\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(z\right)}{{\mathrm{\Σ}}_z{F}_{l+1|{\mathrm{\λ}}_l}\left(z\right)},z\∈\left\{\mathrm{0,1}\right\};$

步骤303）：计算条件转移概率，

${\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)}{{\mathrm{\Σ}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\λ}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\λ}}_l)};$

步骤304）：计算平均噪声语音出现概率，

步骤305）：计算时间依赖平滑因子， ${\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)=\frac{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)}{\mathrm{\α}{\stackrel{\‾}{\mathrm{\γ}}}_l\left(z\right)+(1-\mathrm{\α}){\mathrm{\γ}}_{l+1|{\mathrm{\λ}}_l}\left(z\right)};$

步骤306）：计算状态均值， ${\mathrm{\μ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\μ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]x_{l+1};$

步骤307）：对新的状态均值进行约束： ${\mathrm{\μ}}_{1,l+1}=\max \{{\mathrm{\μ}}_{1,l+1},{\mathrm{\μ}}_{0,l+1}+\mathrm{\δ}\},l\≥M;$

步骤308）：计算新的状态方差， ${\mathrm{\κ}}_{z,l+1}={\widetilde{\mathrm{\α}}}_{l+1}\left(z\right){\mathrm{\κ}}_{z,l}+[1-{\widetilde{\mathrm{\α}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\μ}}_{z,l})}^2;$

步骤309）：对新状态方差进行约束，κ_1,l+1＝max{κ_0,l+1,κ_1,l+1}，l≥M；

步骤310）：计算平均转移概率， ${\mathrm{\ξ}}_{l+1}(y,z)=\mathrm{\α}{\mathrm{\ξ}}_l(y,z)+(1-\mathrm{\α}){\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z);$

步骤311）：计算状态概率， $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\ξ}}_{l+1|{\mathrm{\λ}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\‾}{\mathrm{\ξ}}}_{l+1}(y,1-z)};$

步骤312）：对新的转移概率进行约束，a_01,l＝max{a_01,l,η}，a_00,l＝1-a_01,l，a_10,l＝max{a_10,l,η}，a_11,l＝1-a_10,l，l≥M；

从以上子步骤中，获得了λ_l+1中所有的参数，从而得到了相应的语音出现概率γ_l+1|λl(1)和噪声信号的功率谱估值μ_0,l+1。

2.根据权利要求1所述的噪声功率谱估计与语音活动性检测方法，该方法的具体步骤包括：

1）对于语音信号在每一个频率分量上提取对数幅度谱包络，对于一个频率分量上的对数幅度谱时间序列x_l＝{x₁,x₂,...,x_l}，建立一个隐马尔可夫模型s_l＝{s₁,s₂,...,s_l},s_t∈{0,1}是其对应的状态序列，1表示语音出现状态，0表示噪声出现状态，λ_l表示从序列x_l中获取的模型参数估值，那么，对于给定一个参数集λ_l，对应的观察值序列x_l的概率密度函数可以表示为：

$p\left(x_l\right|{\mathrm{\λ}}_l)=\underset{s_l}{\mathrm{\Σ}}p\left(s_l\right|{\mathrm{\λ}}_l)p\left(x_l\right|{\mathrm{\λ}}_l,s_l);$

其中，p(s_l|λ_l)表示状态序列s_l出现的先验概率，高斯分量表示为：

$p\left(s_l\right|{\mathrm{\λ}}_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}{a}_{s_{t-1},s_t};$

这里的

表示状态转移概率，

表示初始状态概率，p(x_l|λ_l,s_l)表示给定状态s_l和参数集λ_l的情况下观察值序列x_l的似然度：

$p\left(x_l\right|{\mathrm{\λ}}_l,s_l)=\underset{t=1}{\overset{l}{\mathrm{\Π}}}b\left(x_t\right|s_t,{\mathrm{\λ}}_l);$

其中，

$b\left(x_t\right|s_t,{\mathrm{\λ}}_l)=\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\κ}}_{s_t,l}}}\exp \{-\frac{1}{2}{(x_t-{\mathrm{\μ}}_{s_t,l})}^2/{\mathrm{\κ}}_{s_t,l}\};$

这里表示状态s_t对应的高斯分布方差，

表示相应的均值，λ_l＝{μ_0,l,μ_1,l,κ_1,l,κ_0,l,a_01,l,a_10,l,a_00,l,a_11,l,π₀,π₁}，参数集中的初始概率π_i不随着时间变化；

这个模型中的μ_0,l就是我们想要估计的噪声，同时，我们可以推导出语音信号在第l帧的某个频点上出现的概率为

2）对于一段语音数据，设定M帧缓存，把前M帧输入信号存入缓存中，提取缓存中M帧的对数幅度谱，代入步骤1）的HMM模型对每一个频点上初始化一个隐马尔可夫模型λ_M，下标M表示初始化的时间窗口长度，l≥M；

3）在得到初始化的模型λ_M之后，从第M+1帧开始，HMM模型采用增量学习的方法，逐帧更新SHMM模型，依次递推得到λ_l；并得出噪声值μ_0,l和语音信号在第l帧的某个频点上的出现概率。

3.根据权利要求1或2所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的步骤1）中提取一帧幅度谱的步骤包括：

首先，对这一帧的数字化声音信号做预处理，设每帧长度为F点，先补零到N点，N≥F，N＝2^j，j为整数，且j≥8，进行N点离散傅里叶变换，得到离散谱

其中，y_l,n表示缓存中第l帧的第n个采样点，Y_l,k表示缓存中第i帧的第k个傅里叶变换值（k＝0,1,...,N-1）；那么，它的幅度值可以计算为 $x_l=10*{\log }_{10}\left[{\mathrm{\Σ}}_{r=-w}^{w}b\left(r\right){\left|Y_{l,k-r}\right|}^2\right],$ 式中，b(r)为加窗函数。

4.根据权利要求3所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的预处理包括加窗、或/和预加重。

5.根据权利要求3所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的加窗函数采用汉宁窗或哈宁窗。

6.根据权利要求1或2所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的步骤2）中HMM的初始化，在某个频点上具体初始化的步骤包括：

步骤201）：通过聚类的方法将M个样本分成两类：

和其中，M₀+M₁＝M，均值较大的一类用上标（1）表示，另一类用上标（0）表示；

两类的均值为

能量较小的一类的均值为

其中， ${\stackrel{\&OverBar;}{\mathrm{\&mu;}}}_{0,M}<{\stackrel{\&OverBar;}{\mathrm{\&mu;}}}_{1,M};$

两类的方差分别为： ${\stackrel{\&OverBar;}{\mathrm{\&kappa;}}}_{0,M}=\frac{1}{M_0}{\mathrm{\&Sigma;}}_{j=1}^{M_0}{(x_{i_j}-{\stackrel{\&OverBar;}{\mathrm{\&mu;}}}_{0,M})}^2,{\stackrel{\&OverBar;}{\mathrm{\&kappa;}}}_{1,M}=\frac{1}{M_1}{\mathrm{\&Sigma;}}_{j=1}^{M_1}{(x_{i_j}-{\stackrel{\&OverBar;}{\mathrm{\&mu;}}}_{1,M})}^2;$

两类的初始化权重系数为： ${\stackrel{\&OverBar;}{a}}_{00,M}={\stackrel{\&OverBar;}{a}}_{01,M}={\stackrel{\&OverBar;}{a}}_{11,M}={\stackrel{\&OverBar;}{a}}_{10,M}=0.5;$

计算新模型的似然度，

并开始迭代运算；在以下的迭代过程中，老的模型参数集合表示为λ′_M，新的模型参数为：开始迭代以前，

L′设置为一个很大的负数，初始化前向因子，

初始化后向因子， ${\stackrel{\&OverBar;}{B}}_M\left(z\right)=1$

步骤202）：计算前向因子： ${\stackrel{\&OverBar;}{F}}_l\left(z\right)=\underset{y}{\mathrm{\&Sigma;}}{\stackrel{\&OverBar;}{F}}_{l-1}\left(z\right){\stackrel{\&OverBar;}{a}}_{y,z,M}b\left(x_l\right|y,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M),$ z,y∈{0,1}；

步骤203）：计算后向因子： ${\stackrel{\&OverBar;}{B}}_l\left(z\right)=\underset{y}{\mathrm{\&Sigma;}}{\stackrel{\&OverBar;}{B}}_{l+1}\left(z\right){\stackrel{\&OverBar;}{a}}_{y,z,M}b\left(x_{l+1}\right|y,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M),$ z,y∈{0,1}；

步骤204）：计算噪声和语音出现概率：

z∈{0,1}；

步骤205）：如果 $\frac{1}{M}{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_t=1|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)<\mathrm{\&zeta;},{\mathrm{\&lambda;}}_M={\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M,$ 停止跌代，其中ζ是接近于零但大于零的小数；

步骤206）：计算转移概率：

$p(s_{l-1}=y,s_l=z|x_l,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)=\frac{{\stackrel{\&OverBar;}{F}}_{l-1}\left(y\right){\stackrel{\&OverBar;}{B}}_l\left(z\right){\stackrel{\&OverBar;}{a}}_{\mathrm{yz},M}b\left(x_l\right|z,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)}{{\mathrm{\&Sigma;}}_z{\stackrel{\&OverBar;}{F}}_{l-1}\left(y\right){\stackrel{\&OverBar;}{B}}_l\left(z\right){\stackrel{\&OverBar;}{a}}_{\mathrm{yz},M}b\left(x_l\right|z,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)};$

步骤207）：计算新的初始化概率

步骤208）：计算新的均值 ${\mathrm{\&mu;}}_{z,M}^{\text{'}}=\frac{{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)x_t}{{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)}{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M;$

步骤209）：对新的均值进行约束：μ′_1,M＝max{μ′_0,M,μ′_0,M+δ}，其中，δ是一个常数，取值范围在0至100之间；

步骤210）：计算新方差 ${\mathrm{\&kappa;}}_{z,M}^{\text{'}}=\frac{{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M){(x_t-{\stackrel{\&OverBar;}{\mathrm{\&mu;}}}_{z,M})}^2}{{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)};$

步骤211）：对新方差进行约束，κ′_1,M＝max{κ′_0,M,κ′_1,M}；

步骤212）：计算新的转移概率， $a_{\mathrm{yz},M}^{\text{'}}=\frac{{\mathrm{\&Sigma;}}_{t=1}^{M}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)}{{\mathrm{\&Sigma;}}_{t-1}^M{\mathrm{\&Sigma;}}_{z}p(s_{t-1}=y,s_t=z|x_t,{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M)};$

步骤213）：计算新模型的似然度

步骤214）：如果满足条件终止迭代，其中，ε是一个很小的数字，如果 $|L^{\&prime;}-\stackrel{\&OverBar;}{L}|>\mathrm{\&epsiv;},{\stackrel{\&OverBar;}{\mathrm{\&lambda;}}}_M={\mathrm{\&lambda;}}_M^{\&prime;},L^{\&prime;}=\stackrel{\&OverBar;}{L},$ 迭代跳转至步骤202）。

7.根据权利要求6所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的步骤201）中的聚类的方法采用LBG非监督聚类或者模糊聚类方法。

8.根据权利要求1所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的步骤3）中的HMM模型采用的增量学习方法包括：递推权重系数、递推均值和递推方差；

其中，所述的递推均值： ${\mathrm{\&mu;}}_{z,\text{'}+1}={\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right){\mathrm{\&mu;}}_{z,l}+[1-{\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right)]x_{l+1},$ 式中，

是一个依赖于语音出现概率的平滑因子，小于1但接近于1；

所述的递推方差： ${\mathrm{\&kappa;}}_{z,l+1}={\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right){\mathrm{\&kappa;}}_{z,l}+[1-{\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\&mu;}}_{z,l})}^2;$

所述的递推转移概率： $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\&OverBar;}{\mathrm{\&xi;}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\&OverBar;}{\mathrm{\&xi;}}}_{l+1}(y,1-z)};$ 或者a_yz,l+1＝βa_yz,l+(1-β)ξ_l+1|λl(y,z)；式中，β是一个小于1但接近于1的平滑因子。

9.根据权利要求1所述的噪声功率谱估计与语音活动性检测方法，其特征在于，所述的基于一阶回归的序贯隐马尔可夫模型的参数递推方法为：

计算HMM的前向因子： $F_{l+1|{\mathrm{\&lambda;}}_l}\left(z\right)=\underset{y}{\mathrm{\&Sigma;}}{F}_{l|{\mathrm{\&lambda;}}_{l-1}}\left(z\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\&lambda;}}_l),$ z∈{0,1}；

计算语音和噪声出现概率，

z∈{0,1}；

计算条件转移概率， ${\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,z)=\frac{F_{l+1|{\mathrm{\&lambda;}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\&lambda;}}_l)}{{\mathrm{\&Sigma;}}_{\mathrm{yz}}{F}_{l+1|{\mathrm{\&lambda;}}_l}\left(y\right)a_{\mathrm{yz},l}b\left(x_{l+1}\right|s_{l+1}=z,{\mathrm{\&lambda;}}_l)};$

计算平均噪声语音出现概率， ${\widetilde{\mathrm{\&gamma;}}}_{l+1}\left(z\right)=a{\widetilde{\mathrm{\&gamma;}}}_l\left(z\right)+(1-\mathrm{\&alpha;}){\mathrm{\&gamma;}}_{l+1|{\mathrm{\&lambda;}}_l}\left(z\right);$

计算时间依赖平滑因子， ${\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right)=\frac{\mathrm{\&alpha;}{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_l\left(z\right)}{\mathrm{\&alpha;}{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_l\left(z\right)+(1-\mathrm{\&alpha;}){\mathrm{\&gamma;}}_{l+1|{\mathrm{\&lambda;}}_l}\left(z\right)};$

计算均值， ${\mathrm{\&mu;}}_{z,l+1}={\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right){\mathrm{\&mu;}}_{z,l}+[1-{\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right)x_{l+1}];$

计算新的方差， ${\mathrm{\&kappa;}}_{z,l+1}={\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right){\mathrm{\&kappa;}}_{z,l}+[1-{\widetilde{\mathrm{\&alpha;}}}_{l+1}\left(z\right)]{(x_{l+1}-{\mathrm{\&mu;}}_{z,l})}^2;$

计算平均转移概率， ${\mathrm{\&xi;}}_{l+1}(y,z)=\mathrm{\&alpha;}{\mathrm{\&xi;}}_l(y,z)+(1-\mathrm{\&alpha;}){\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,z);$

计算转移概率， $a_{\mathrm{yz},l+1}=a_{\mathrm{yz},l}+\frac{\frac{{\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,z)}{a_{\mathrm{yz},l}}-\frac{{\mathrm{\&xi;}}_{l+1|{\mathrm{\&lambda;}}_l}(y,1-z)}{1-a_{\mathrm{yz},l}}}{\frac{K}{a_{\mathrm{yz},l}^2}{\stackrel{\&OverBar;}{\mathrm{\&xi;}}}_{l+1}(y,z)+\frac{K}{{(1-a_{\mathrm{yz},l})}^2}{\stackrel{\&OverBar;}{\mathrm{\&xi;}}}_{l+1}(y,1-z)}.$

Images