噪声估计之MCRA

$MCRA$¹，全称为最小值控制的递归平均，是cohen提出的一种常用的噪声估计方法，处理流程框图 ² 如下

从命名上以及以上框图能看出来着个方法主要包含两个部分，噪声谱递归平均和最小值控制(跟踪)，下面分别看看这两个部分

1. 噪声谱估计（递归平均）

还是老套路先定义信号表达式
$$Y(k,\ell )=\sum _{n=0}^{N-1}y(n+\ell M)h(n)e^{-j\frac{2\pi }{N}nk}\text{\hspace{1em}(1)}$$

$k$:frequency bin index

$\ell$:time frame index

$h$:analysis window of size N

$M$:frame update step in time

定义两个假设，在第 $\ell$ 帧第$k$ 个频率点上语音不存在和语音存在$H_0(k,\ell ),H_1(k,\ell )$概率如下：
$$\begin{array}{rl}& H_0(k,\ell ):Y(k,\ell )=D(k,\ell )\\ & H_1(k,\ell ):Y(k,\ell )=X(k,\ell )+D(k,\ell )\end{array}\text{\hspace{1em}(2)}$$
其中$X(k,\ell )和D(k,\ell )$分别为纯净语音和噪声的$STFT$表示

定义${\lambda }_d(k,\ell )=E[{\left|\begin{array}{c}D(k,\ell )\end{array}\right|}^2]$为$k帧、\ell 子带$处噪声谱，那么就可以在无语音段用一个时间维度上的递归平滑来更新噪声，用公式表示如下
$$\begin{array}{l}{H}_0^{\prime }(k,\ell ):{\hat{\lambda }}_d(k,\ell +1)={\alpha }_d{\hat{\lambda }}_d(k,\ell )+\left(1-{\alpha }_d\right)|Y(k,\ell ){|}^2\\ H_1^{\prime }(k,\ell ):{\hat{\lambda }}_d(k,\ell +1)={\hat{\lambda }}_d(k,\ell )\end{array}\text{\hspace{1em}(3)}$$
其中${\alpha }_d(0<{\alpha }_d<1)$为平滑因子

用$p^{\prime }\triangleq P(H_1^{\prime }|Y(k,\ell ))$表示语音存在的条件概率，则(3)式可以写成如下形式
$${\tilde{\lambda }}_d(k,\ell +1)={\tilde{\alpha }}_d(k,\ell ){\hat{\lambda }}_d(k,\ell )+\left[1-{\tilde{\alpha }}_d(k,\ell )\right]|Y(k,\ell ){|}^2\text{\hspace{1em}(4)}$$
其中${\tilde{\alpha }}_d(k,\ell )$为时变的平滑参数
$${\tilde{\alpha }}_d(k,\ell )\triangleq {\alpha }_d+(1-{\alpha }_d)p^{\prime }(k,\ell )\text{\hspace{1em}(5)}$$
上面的(4)、(5)两式就是递归平均更新噪声谱的核心内容，现在的问题就是要求出时变参数${\tilde{\alpha }}_d(k,\ell )$，也就是要求出语音存在概率这个关键变量$p^{\prime }(k,\ell )$

2. 语音存在概率（最小值控制）

2.1. 最小值跟踪

​ 语音存在概率由带噪语音当前的能量和指定长度窗内谱最小值的比值来计算，先对带噪语音谱分别做时间、频率两个维度上的平滑

频率平滑：
$$S_f(k,\ell )=\sum _{i=-w}^{w}b(i)|Y(k-i,\ell ){|}^2\text{\hspace{1em}(5)}$$
时间平滑:
$$S(k,\ell )={\alpha }_s(k,\ell )S(k,\ell -1)+S_f(k,\ell )\text{\hspace{1em}(6)}$$
其中${\alpha }_s(0<{\alpha }_s<1)$为平滑常数

谱最小值$S_{min}(k,\ell )$搜索过程如下：

初始化:

$S_{min}(k,\ell )=S(k,0)$

$S_{tmp}(k,\ell )=S(k,0)$

然后按样本点(频谱)比较
$$\begin{gathered} S_{min}(k,\ell )=min\left\{\begin{array}{c}{S}_{min}(k,\ell -1),S(k,\ell )\end{array}\right\} \\ {S}_{tmp}(k,\ell )=min\left\{\begin{array}{c}{S}_{tmp}(k,\ell -1),S(k,\ell )\end{array}\right\}\text{\hspace{1em}(7)} \end{gathered}$$
这个时候$S_{min}和S_{tmp}$都还是相等的，当比较了L帧（mod($\ell$,L)=0）后
$$\begin{array}{rl}& S_{min}(k,\ell )=min\left\{\begin{array}{c}{S}_{tmp}(k,\ell -1),S(k,\ell )\end{array}\right\}\\ & S_{tmp}(k,\ell )=S(k,\ell )\end{array}\text{\hspace{1em}(8)}$$
重复(7)、(8)过程得到最小值谱，其中，搜索窗的长度L会影响到噪声的跟踪速度，一般按照经验选0.5s~1.5s左右。

2.2. 语音存在概率计算

​ 定义带噪语音能量与局部最小能量与比$S_r(k,\ell )$如下
$$S_r(k,\ell )\triangleq \frac{S(k,l)}{S_{min}(k,\ell )}\text{\hspace{1em}(9)}$$
定义二值$I(k,\ell )$如下
$$I(k,\ell )=\{\begin{array}{c}\begin{array}{rl}& 1,S_r(k,\ell )>\delta \\ & 0,otherwise\end{array}\end{array}\text{\hspace{1em}(10)}$$
最后，语音存在概率更新如下：
$$\widehat{p^{\prime }}(k,\ell )={\alpha }_p\widehat{p^{\prime }}(k,\ell -1)+\left(1-{\alpha }_p\right)I(k,\ell )\text{\hspace{1em}(11)}$$
其中$\delta 为预先设定的门限，{\alpha }_p(0<{\alpha }_p<1)$为平滑常数

3.code & result

整个过程比较简单，参照上述过程编写程序估计噪音，结合谱减法降噪效果如下


对比处理前后可以看到背景噪声有了明显的消除，当然因为谱减法的原因引入了较多的音乐噪声，这个时候还可以尝试其它的谱修改方法如维纳滤波或者OMLSA等
References:

---

1. Cohen, I., & Berdugo, B. (2002). Noise estimation by minima controlled recursive averaging for robust speech enhancement. IEEE Signal Processing Letters, 9(1), 12–15 ↩︎

2. Cohen, I., & Berdugo, B. (2001). Speech enhancement for non-stationary noise environments. Signal Processing, 81(11), 2403–2418 ↩︎