GWPE - Multi-channel Dereverberation

最近在一批实录的数据上尝试了一下GWPE算法，发现有一定的效果，随之研究了一下原始论文，并在这里对其原理做一下简要介绍。

GWPE算法在提出的时候，对应用的声学环境没有建立一些常见的约束条件，比如speaker数量（原始WPE算法要求single source），noise环境等等，因此在一些复杂条件下的鲁棒性会相对强一些。原论文中将一些随机向量用黑体大写字母表示，实际观测向量用小写字母，在这篇文章中，我统一用后者表示。

从论文题目可以看出，GWPE解决的是MIMO下，impluse response的shorten问题，也可以理解为multi-channel的dereverbration。主要贡献在于使用了一种新的相关性度量作为predict filter估计的代价函数，并导出了对应的估计方法。

New Cost Function

先介绍一下是PE（predict error），$\mathbf{y}_{f,t}$表示$N$路麦克风的观测信号，滤波配置（$\Delta, K, \mathcal{G}_f=\{\boldsymbol{G}_{f,\tau}\}_{\tau = \Delta}^{\Delta + K -1}$）下对$\mathbf{y}_{f,t}$的PE表示为：

$$\hat{\mathbf{x}_{f,t}} = \mathbf{y}_{f,t} - \sum_{\tau = \Delta}^{\Delta + K - 1} \boldsymbol{G}_{f,t}^H \mathbf{y}_{f,\tau} \tag{1}$$

其中第二项

$$\mathbf{x}_{f,t} = \sum_{\tau = \Delta}^{\Delta + K - 1} \boldsymbol{G}_{f,t}^H \mathbf{y}_{f,\tau}$$

称为在该滤波配置下对$\mathbf{y}_{f,t}$的线性预测。$K$表示滤波器阶数，$\Delta$表示延迟，$\boldsymbol{G}_{f,\tau} \in \mathbf{C}^{N \times N}, \mathbf{x}_{f,t}, \mathbf{y}_{f,t} \in \mathbf{C}^{N \times 1}$。

WPE算法中，在PE中引入权重因子作为cost function，GWPE使用新的Hadamard-Fischer mutual correlation作为cost function，定义如下：

$$C_{HF}(\boldsymbol{U}_1, \dots, \boldsymbol{U}_N) = \frac{1}{N} \sum_n \log (\det E(\boldsymbol{U}_n\boldsymbol{U}_n^H)) - \log (\det E(\boldsymbol{U}\boldsymbol{U}^H)) \tag{2}$$

$\boldsymbol{U}_{1\cdots N}$表示复值随机向量。$C_{HF}$具有非负性，当且仅当$\boldsymbol{U}_{1\cdots N}$两两完全不相关时取0值。

借用$C_{HF}$来衡量在滤波矩阵$\mathcal{G}_f = \{\boldsymbol{G}_{f,\tau}\}_{\tau = \Delta}^{\Delta + K -1}$下线性预测$\mathbf{x}_{f,t} = \sum_{\tau = \Delta}^{\Delta + K - 1} \boldsymbol{G}_{f,t}^H \mathbf{y}_{f,\tau}$的相关性：

$$\mathcal{J}(\mathcal{G}_f) = \frac{1}{T}\sum_t \log \left( \det E(\mathbf{x}_{f,t}\mathbf{x}_{f,t}^H) \right) - \log (\det E(\boldsymbol{X}_f\boldsymbol{X}_f^H)) \tag{3}$$

可以证明，$(3)$式的第二项是一个常数，因此，$\mathcal{J}(\mathcal{G}_f)$可以简化为：

$$\mathcal{J}(\mathcal{G}_f) = \frac{1}{T}\sum_t \log \left( \det E(\mathbf{x}_{f,t}\mathbf{x}_{f,t}^H) \right) \tag{4}$$

上式便是GWPE的cost function。到此为止，MIMO的response shorten问题转化成了：

$$\mathcal{G}_f = \arg \min_{\mathcal{G}_f} \mathcal{J}(\mathcal{G}_f)$$

的优化问题。

Estimate Filters

$(4)$式的最小值没有解析解，因此，需要导出一种稳定的求解方法。论文中给出的辅助函数方法，即构造辅助函数$\mathcal{\hat{J}}(\mathcal{G}_f, \mathcal{L}_f)$，将原始问题分解成两个子问题依次优化：

$$\mathcal{\hat{L}}_f = \arg \min_{\mathcal{L}_f} \mathcal{\hat{J}}(\mathcal{\hat{G}}_f, \mathcal{L}_f) \tag{5} \\ \mathcal{\hat{G}}_f = \arg \min_{\mathcal{G}_f} \mathcal{\hat{J}}(\mathcal{G}_f, \mathcal{\hat{L}}_f)$$

其中$\mathcal{L}_f = \{\boldsymbol{\Lambda}_{f,t} \}_{t = 0}^{T - 1}$。

辅助函数定义为：

$$\mathcal{\hat{J}}(\mathcal{G}_f, \mathcal{L}_f) = \frac{1}{T}\sum_t \left( E(\mathbf{x}_{f,t}^H \boldsymbol{\Lambda}_{f,t}^{-1} \mathbf{x}_{f,t}) - N + \log(\det(\boldsymbol{\Lambda}_{f,t})) \right)$$

依据性质：

$$\log \left| \det E(\boldsymbol{U}\boldsymbol{U}^H) \right| \leqslant E(\boldsymbol{U}^H\boldsymbol{\Lambda}^{-1}\boldsymbol{U}) -N + \log(\det \boldsymbol{\Lambda})$$

当且仅当$\boldsymbol{\Lambda} = E(\boldsymbol{U}\boldsymbol{U}^H)$时取等号。因此，对于$\mathcal{\hat{J}}(\mathcal{G}_f, \mathcal{L}_f)$：

$$\mathcal{J}(\mathcal{G}_f) \leqslant \mathcal{\hat{J}}(\mathcal{G}_f, \mathcal{L}_f)$$

取得等号时，$(5)$式的解为：

$$\mathcal{\hat{L}}_f = \{\boldsymbol{\hat{\Lambda}}_{f,t} \}_{t = 0}^{T - 1} = \{E(\mathbf{\hat{x}}_{f,t}\mathbf{\hat{x}}_{f,t}^H)\}_{t=0}^{T-1}$$

这一步的核心在于$E(\mathbf{\hat{x}}_{f,t}\mathbf{\hat{x}}_{f,t}^H)$（spatial correlation matrix）的估计，最常见的方法是做time average：

$$\boldsymbol{\hat{\Lambda}}_{f,t} = E(\mathbf{\hat{x}}_{f,t}\mathbf{\hat{x}}_{f,t}^H) = \sum_{k=t-\delta}^{t+\delta} = \frac{1}{2\delta+1} \mathbf{\hat{x}}_{f,k}\mathbf{\hat{x}}_{f,k}^H$$

论文中也提供了其他四种可选方法，可用于简化计算。

$\mathcal{\hat{G}}_f$的解相推导对麻烦一些，我只说一下计算方法，首先构造$\boldsymbol{\psi}_{f,t}$矩阵：

$$\boldsymbol{\psi}_{f,t} = [\boldsymbol{Y}_{f,t}^H, \cdots, \boldsymbol{Y}_{f,t-K+1}^H]^T \in \mathbf{C}^{KN^2 \times N}$$

其中：

$$\boldsymbol{Y}_{f,t} = \begin{bmatrix} \mathbf{y}_{f,t} & & \mathbf{o}\\ & \ddots & \\ \mathbf{o} & & \mathbf{y}_{f,t} \\ \end{bmatrix} \in \mathbf{C}^{N^2 \times N}$$

其次，计算$\boldsymbol{g}_f = \overline{\boldsymbol{R}_f^{-1}\boldsymbol{r}_f}$，其中：

$$\begin{align} \boldsymbol{R}_f &= \sum_t \boldsymbol{\psi}_{f,t - \Delta} \boldsymbol{\hat{\Lambda}}_{f,t} \boldsymbol{\psi}_{f,t - \Delta}^H \in \mathbf{C}^{KN^2 \times KN^2} \\ \boldsymbol{r}_f &= \sum_t \boldsymbol{\psi}_{f,t - \Delta} \boldsymbol{\hat{\Lambda}}_{f,t} \mathbf{y}_{f,t} \in \mathbf{C}^{KN^2 \times 1} \end{align}$$

$\boldsymbol{g}_f$和$\boldsymbol{G}_{f,\tau}$的对应关系为：

$$\boldsymbol{g}_l = \begin{bmatrix} \boldsymbol{G}_{f,\Delta,:, 1} \\ \vdots \\ \boldsymbol{G}_{f,\Delta,:, N} \\ \vdots \\ \boldsymbol{G}_{f,\Delta + K - 1,:, N} \end{bmatrix} \\$$

最后reshape一下成为$\mathcal{\hat{G}}_f$的更新值。

Application

目前GWPE由fgnt开源了一版python的实现（https://github.com/fgnt/nara_wpe），效果挺好，可以参考一下实现细节，有空的话，我自己也想实现一版。

GWPE最直接的用处就是multi-channel的dereverbration，如下图（这是4麦的实录数据，混响不是非常严重），对比来看，还是可以比较明显的看出GWPE处理之后的音频扫尾现象少了很多。由于混响对识别任务而言是一个影响较大的因素，因此，在WER上往往可以比较明显的对比出差距。

其次就是和前端的一些mask估计，beamforming结合起来用，因为增强/分离/定位这些任务也是混响敏感的，用GWPE做一遍数据预处理通常会有ASR上的增益。fgnt在interspeech2018上有一篇文章[2]专门分析WPE和beamforming的结合，有兴趣的同学可以参考一下。目前我手上的实验也验证了这一点，参考结果如下（WER绝对提升）：

RAW-CH1	GWPE-CH1	DS	CGMM-MVDR	GWPE-CGMM-MVDR
0%	2.07%	1.98%	3.37%	4.35%

可以看出，单独过一遍GWPE就可以获得2%的绝对提升，在CGMM-MVDR基础上，替换输入为GWPE结果之后，可以继续获得一个点的绝对提升。

Reference

[1]. Yoshioka T, Nakatani T. Generalization of multi-channel linear prediction methods for blind MIMO impulse response shortening[J]. IEEE Transactions on Audio, Speech and Language Processing, 2012, 20(10): 2707-2720.
[2]. Drude L, Boeddeker C, Heymann J, et al. Integrating Neural Network Based Beamforming and Weighted Prediction Error Dereverberation[J]. Proc. Interspeech 2018, 2018: 3043-3047.