单通道语音增强之统计信号模型

1. 信号估计理论简述

信号估计理论是现代统计处理的基础课题[@ZhangXianDa2002ModernSP],在通信、语音、图像领域均有广泛应用。语音增强,就是从带噪的语音测量信号中估计原始的无噪语音,这是典型的信号估计问题。
《语音增强–理论与实践》[@loizou2007speech]一书中列举了用于语音增强的一系列统计模型。

假设麦克风采集到的带噪语音序列为y[n],并且噪声都是加性噪声。则带噪语音序列为无噪语音序列与噪声序列的和。原始语音信号与噪声均可视为随机信号。

y[n] = x[n] + d[n]

在时域对x[n]进行估计是非常困难的,通过傅立叶变换,我们可以将信号分解为频域上互相独立的系数。信号估计模型转变为对每一个频点的系数进行估计的模型,不同频点之间的参数是相互独立的。

Y(\omega_{k}) = X(\omega_{k}) + D(\omega_{k})

这个方法就叫做统计信号谱分析(Statistical Spectral
Analysis)。显然地,纯净信号谱X(\omega_{k})带有幅度与相位两参数,我们实际上是对幅度X_{k}和相位\theta_{y}(k)进行参数估计。

Y_{k}e^{j\theta_{y}(k)} = X_{k}e^{j\theta_{x}(k)} + D_{k}e^{j\theta_{d}(k)}

重组谱幅度和谱相位估计值即可恢复纯净语音谱,估计值用上标来表示:\hat{X}(\omega_{k})=\hat{X}(k)e^{j\hat{\theta}_{x}(k)}

实际信号的幅度和相位是不方便直接用在运算过程中的,因为信号取值范围不定,且瞬时变化。在噪声抑制领域,更常用的语音谱估计方法是对抑制增益(Suppression
Gain)进行估计,不同的估计准则称为抑制准则(Supression Rule)。

\hat{X}(\omega_{k}) = H_{k}Y(\omega_{k})

通常会根据先验信噪比、后验信噪比来估计抑制增益H_{k}
并且可以在只有噪声出现的时刻更新H_{k}
在语音存在的时刻进行抑制,无须每帧去调用噪声抑制算法,计算过程比直接估计信号谱灵活。

综上,语音增强的典型流程就是:

  1. 对带噪语音y[n]分帧, 每一帧进行DFT得到Y(\omega_{k})

  2. 估计或者沿用上一帧的抑制增益H_{k},得到纯净语音谱\hat{X}(\omega_{k})

  3. \hat{X}(\omega_{k})进行IDFT,得到纯净语音序列的估计x[n]

为了估计模型的建模,对测量信号、估计信号、噪声信号都需要作一些数学上的假设和简化。其中对噪声一般会作以下假设:

  1. 噪声是与语音独立的加性噪声;

  2. 每一帧噪声的统计分布是稳态的;

  3. 噪声的傅立叶级数是零均值复高斯分布。

\left[D(\omega_{k})\right] = \frac{1}{\pi\lambda_{d}(k)} \exp\left[-\frac{|D(\omega_{k})|^2}{\lambda_{d}(k)}\right]

2. 最大似然估计ML

如果不考虑信号的先验分布,即认为信号值是确定信号,而不是随机信号,我们只需要分析含有信号x为参数的带噪信号y的概率分布p(y;x),并使之最大。这种估计方法叫最大似然估计器(Maximum
Likelihood
Estimation)。将y的带参概率分布p(y;x)称为似然函数(Likelihood
Function)。对纯净信号x的估计,表达为求解合适的x值,使得似然函数p(y;x)最大。

\hat{x} = \arg max_{x}p(y;x)

文献[@mcaulay1980speech]最早将最大似然估计法用在语音增强领域。对于纯净语音,可以假设纯净语音幅度X_{k}和相位\theta_{y}(k)是未知但确定,无需考虑其先验概率分布。最大似然语音增强模型表达为:

\hat{X}_{k},\hat{\theta}_{k} = \arg max_{X_{k},\theta_{k}} p\left[Y(\omega_{k});X_{k},\theta_{k}\right]

D_{k}e^{\theta_{d}(k)}=Y_{k}e^{\theta_{y}(k)} - X_{k}e^{\theta_{x}(k)}代入噪声零均值复高斯分布公式中,得到:

p\left[Y(\omega_{k});X_{k},\theta_{k}\right] = \frac{1}{\pi \lambda_{d}(k)}\exp\left[-\frac{|Y(\omega_{k}) - X_{k}e^{j\theta_{x}(k)}|^2}{\lambda_{d}(k)}\right]

H_{k} = \frac{1}{2}+\frac{1}{2} \sqrt{\frac{\xi_{k}}{1+\xi_{k}}}

3.贝叶斯估计

如果比最大似然估计更进一步,考虑待估计量x也是随机变量,且x的先验分布为p(x),这种假设下的估计方法叫做贝叶斯估计[@ZhangXianDa2002ModernSP]。定义估计值\hat{x}与实际值x之间的误差函数为c(\hat{x},x),贝叶斯估计器的目标即为找出是平均误差E[c(\hat{x},x)]最小的估计值x

\hat{x} = \arg min_{\hat{x}} E\left[ c(\hat{x}, x) \right]

对于待估计的纯净语音谱,贝叶斯估计器可以表达为:

\hat{X}(\omega_{k}) = \arg min_{\hat{X_{k}}} E[c(\hat{X}(\omega_{k}), X(\omega_{k}))]

误差c(\hat{x},x)的期望值取决于待测信号与测量信号的联合概率分布。

E[c(\hat{x}, x)] = \int_{x}\int_{y}c(\hat{x}, x)p(x, y)dxdy

根据条件概率密度公式对联合概率密度进行分解:

E[c(\hat{x}, x)] = \int_{x}\int_{y}c(\hat{x}, x)p(x|y)p(y)dxdy = \int_{y} \left[ \int_{x}c(\hat{x}, x)p(x|y)dx\right] p(y)dy

因为估计值\hat{x}p(y)相互独立,可以把外层对y的积分消除,也就是:

\hat{x} = \arg min_{\hat{x}} E[c(\hat{x}, x)|y]

误差函数c(\hat{x}, x)的计算模型,会引出不同种类的估计器。典型的误差函数有几种类型[@ZhangXianDa2002ModernSP]:

  1. 平方误差函数,对应最小均方估计。

    c(\hat{x}, x) = (\hat{x} - x)^2

  2. 绝对值误差函数,对应条件中位数估计。

    c(\hat{x}, x) = |\hat{x} - x|

  3. 均匀误差函数,对应最大后验估计。

c(\hat{x}, x) = \left\{ \begin{array}{lr} 1 \quad (|\hat{x}-x|\geq\Delta/2) & \\ &\Delta>0\\ 0 \quad (|\hat{x}-x|<\Delta/2) & \end{array} \right.

3.1 最小均方估计(MMSE)

最小均方估计(Minimum Mean Square Error
Estimation)[@ephraim1985speech]使用平方误差函数,平均误差为:

E[c(\hat{x}, x)|y]=\int_{x}(\hat{x} - x)^2 p(x|y)dx

为求平均误差的极大值,可对估计量\hat{x}求导,并求极值点。

\frac{\partial E[c(\hat{x}, x)|y]}{\partial \hat{x}}=\int_{x}2(\hat{x} - x) p(x|y)dx =0

显然极值点的\hat{x}为被估计量x的条件均值:

\hat{x} =\int_{x}x p(x|y)dx = E[ x| y]

亦可以使用贝叶斯公式展开,得到:

\begin{aligned} \hat{x} = \int_{x}x \frac{p(y|x)p(x)}{p(y)} dx \\ = \frac{\int_{x}xp(y|x)p(x)dx}{\int_{x}p(y|x)p(x)dx} \end{aligned}

在语音增强的模型里,纯净语音谱的估计为其在带噪语音谱下的条件均值。

\hat{X}(\omega_{k}) = E[X(\omega_{k})| Y(\omega_{k})]

上式中,为得到纯净语音谱需要分别估计谱幅度和谱相位: \hat{X}, \hat{\theta}_{x}

\hat{X}_{k}e^{j\hat{\theta}_{x}(k)} = E[ X_{k}e^{j\theta_{x}(k)}| Y(\omega_{k})]

同时估计谱幅度和谱相位是很难的,研究者提出了许多分别估计谱幅度和谱相位的方法,估计完成后再用两者重组复语音信号。

\hat{X}_{k} = E[ X_{k}| Y(\omega_{k})]

e^{j\hat{\theta}_{x}(k)} = E[ e^{j\theta_{x}(k)}| Y(\omega_{k})]

3.1.1 MMSE谱幅度估计

最小均方根估计器MMSE short-time spectral amplitude

\hat{X_{k}} = E[ X_{k}| Y(\omega_{k})]

v_{k} = \frac{\xi_{k}}{1+\xi_{k}}\gamma_{k}

H_{k} =\frac{\sqrt{\pi v_{k}}}{ 2\gamma_{k}}[(1+v_{k})I_{0}(\frac{v_k}{2})+v_{1}I_{1}(\frac{v_{k}}{2})]\exp(\frac{-v_{k}}{2})

    def mmse_stsa_gain(parameters=None):
        gamma = parameters['gamma']
        ksi = parameters['ksi']
            
        vk = ksi * gamma / (1 + ksi)
        j0 = np.i0(vk / 2)
        j1 = np.i1(vk / 2)
            
        A = sqrt(pi * vk) / 2 / gamma
        B = (1 + vk) * j0 + vk * j1
        C = np.exp(-0.5 * vk)
        gain = A * B * C
        return gain

3.1.2 MMSE对数谱幅度估计

对数最小均方根估计器The MMSE log spectral amplitude (MMSE-LSA),
或者缩写为LogMMSE估计器。

c(\hat{X_{k}}, X_{k}) = (\log{X_{k}}- \log{X_{k}})^2

\log{\hat{X_{k}}} = E[ \log{X_{k}}| Y(\omega_{k})]

H_{k} = \frac{\xi_{k}}{1+\xi_{k}}\exp(\frac{1}{2}\int_{v_{k}}^{\infty}\frac{e^{-t}}{t}dt)

    def logmmse_gain(parameters=None):
        gamma = parameters['gamma']
        ksi = parameters['ksi']
        A = ksi / (1 + ksi)
        vk = A * gamma
        ei_vk = 0.5 * expn(1, vk)
        gain = A * np.exp(ei_vk)
        return gain

3.1.3 MMSE平方谱幅度估计

频谱幅度平方估计器MMSE magnitude squared[@wolfe2003efficient]

\hat{X_{k}^2} = E[ X_{k}^2| Y_{k}]

\hat{X_{k}^2} = \frac{\xi_{k}}{1+\xi_{k}} (\frac{1+v_{k}}{\gamma_{k}})Y_{k}^2

H_{k} = \sqrt{\frac{\xi_{k}}{1+\xi_{k}} (\frac{1+v_{k}}{\gamma_{k}})}

    def mmse_sqr_gain(parameters=None):
        gamma = parameters['gamma']
        ksi = parameters['ksi']
            
        vk = ksi * gamma / (1 + ksi)
        j0 = np.i0(vk / 2)
        j1 = np.i1(vk / 2)
            
        A = ksi / (1 + ksi)
        B = (1 + vk) / gamma
        gain = sqrt(A * B)
        return gain

3.2 最大后验估计 MAP

当贝叶斯估计器采用均匀误差函数时,平均误差为:

\begin{aligned} E[c(\hat{x}, x)|y] &= \int_{-\infty}^{\hat{x}-\Delta/2} p(x|y)dx + \int_{\hat{x}+\Delta/2}^{\infty} p(x|y)dx \\ &= 1-\int_{\hat{x}\Delta/2}^{\hat{x}+\Delta/2} p(x|y)dx \end{aligned}

显然要使得平均误差最小,就是要求目标估计\hat{x},使得
p(x|y)最大。这种估计模型称作最大后验估计(Maximum A Posteriori
Estimation)。这个模型的意思是只有估计值\hat{x}与原始值x相等,误差才为0,其他时候误差均匀为1。估计值可以表达为:

\hat{x} =\arg min_{x} p(x|y)

在语音增强的模型里,纯净语音谱的估计为其在带噪语音谱下的条件均值。

\hat{X}(\omega_{k}) = \arg min_{ X(\omega_{k})}\quad p\left[X(\omega_{k})| Y(\omega_{k})\right]

文献[@wolfe2003efficient]提出了两种基于最大后验估计(Maximum A Posteriori
Estimation)的语音增强算法。一种是同时求解幅度和相位的混合最大后验估计:
另外一种是单纯估计幅度的方法,两种估计器最后的噪声抑制增益略有不同。

3.2.1 幅度和相位混合最大后验估计

[@wolfe2003efficient]

H_{k} = \frac{\xi_{k}+\sqrt{\xi_{k}^2+2(1+\xi_{k})\frac{\xi_{k}}{\gamma_{k}}}}{2(1+\xi_{k})}

def map_joint_gain(parameters=None):
    gamma = parameters['gamma']
    ksi = parameters['ksi']
    
    eps = 1e-6
    gain = (ksi + sqrt(ksi^ 2 + 2 * (1.0 + ksi)* ksi/ (gamma + eps))) / 2.0/ (1.0 + ksi)
    return gain

3.2.2 纯幅度最大后验估计

[@wolfe2003efficient]

H_{k} = \frac{\xi_{k}+\sqrt{\xi_{k}^2+(1+\xi_{k})\frac{\xi_{k}}{\gamma_{k}}}}{2(1+\xi_{k})}

def map_sa_gain(parameters=None):
    gamma = parameters['gamma']
    ksi = parameters['ksi']
    
    eps = 1e-6
    gain = (ksi + sqrt(ksi^ 2 + (1.0 + ksi)* ksi/ (gamma + eps))) / 2.0/ (1.0 + ksi)
    return gain


参考文献

[1] 张贤达, 现代信号处理. 清华大学出版社有限公司, 2002.

[2] P. C. Loizou, Speech enhancement: theory and practice. CRC press, 2007.

[3] R. McAulay and M. Malpass, Speech enhancement using a soft-decision noise suppression filter," IEEE Transactions on Acoustics, Speech, and Signal Processing, vol. 28, no. 2, pp. 137-145,1980.

[4] Y. Ephraim and D. Malah, Speech enhancement using a minimum mean-square error log-spectral amplitude estimator," IEEE transactions on acoustics, speech, and signal processing,vol. 33, no. 2, pp. 443-445, 1985.

[5] P. J. Wolfe and S. J. Godsill, Efficient alternatives to the ephraim and malah suppression rule for audio signal enhancement," EURASIP Journal on Applied Signal Processing, vol. 2003, pp.1043-1051, 2003.

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 200,302评论 5 470
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 84,232评论 2 377
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 147,337评论 0 332
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 53,977评论 1 272
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 62,920评论 5 360
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,194评论 1 277
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,638评论 3 390
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,319评论 0 254
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,455评论 1 294
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,379评论 2 317
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,426评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,106评论 3 315
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,696评论 3 303
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,786评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,996评论 1 255
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,467评论 2 346
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,043评论 2 341

推荐阅读更多精彩内容