CN108899052B - 一种基于多带谱减法的帕金森语音增强方法 - Google Patents

Abstract

本发明揭示了一种多带谱减法的帕金森语音增强方法，包括以下步骤：通过非接触式方式对帕金森患者进行语音信号采集；采集到的帕金森病患者语音进行采样和短时帧划分；对每一帧语音进行短时能量计算，并根据短时能量进行自适应频带划分；计算每个窗长内含噪声语音的短时功率谱、时变的平滑参数及平滑功率谱；计算偏差修正因子；搜索功率谱密度的最小值；计算并更新噪声的功率谱密度；及采用多带谱减法进行语音增强；本发明在非实验环境下，对录制的帕金森病患者语音进行有效的去噪处理，提高帕金森病患者语音质量，最终有效提高帕金森病的诊断或帕金森病情严重程度预测的准确性。

Description

一种基于多带谱减法的帕金森语音增强方法

技术领域

本发明涉及一种语音增强方法，尤其涉及一种基于多带谱减法的帕金森语音增强方法，属于语音增强领域。

背景技术

帕金森病主要是由于中脑黑质致密部多巴胺能神经元变异，以及残存神经元细胞多巴胺生物合成能力下降，导致纹状体区多巴胺逐渐缺乏而引起的。目前，纹状体区多巴胺递质降低的原因尚不知晓，这也是导致帕金森病只能进行康复诊疗，却无法彻底根治的重要原因。帕金森病将给老年人群及其家庭生活带来极大不便。特别是晚期帕金森病患者，通常出现肌肉僵直症状，只能长期卧床。这不仅给患者本人带来极大痛苦，也给患者家庭造成极大不便。研究表明，帕金森病患者如能尽早发现病情，通过及时的治疗，能够有效延缓病情进展。因此，帕金森病的早期诊断和康复治疗的研究具有十分重要的意义。

研究表明，90％的帕金森病患者都有一定程度的语音障碍。利用语音障碍来诊断和治疗帕金森病，不但费用低廉，而且患者可通过非接触式方法自助测量，简单方便，易于实现远程诊断。通过语音信号进行帕金森病的诊断和康复治疗已经得到国内外研究学者的广泛关注。Little等人使用支持向量机(Support Vector Machine,SVM)对帕金森病患者的语音进行检测，判断用户是否患有帕金森病。叶晓江等人基于该研究，利用安卓技术开发出了基于语音的帕金森检测系统，用户在室内环境下自行采集持续的长元音/a/，然后将采集到的语音上传至该系统，由系统进行帕金森病的诊断。然而，在实际应用场景下，室内广泛存在的环境噪声无疑会对语音的平稳性及其它质量产生影响，进而影响语音特征提取的准确性，最终会影响基于语音和机器学习方法进行的帕金森病相关研究结果的准确性。因此，对帕金森病患者使用非接触方法采集到的连续长元音/a/进行语音增强，是提高帕金森症检测系统性能的关键之一。

语音增强的目的是从带噪声的语音中提取出纯净的语音，同时保证语音的质量(quality)以及可懂度(intelligibility)。谱减法的基本思想是假设噪声具有统计平稳性，用语音的“寂静段”来估计噪声的功率谱；然后，在发声段，用带噪声的语音功率谱减去噪声功率谱，就可以直接得到纯净的语音谱。研究表明，当频谱的信噪比大于8dB时，可直接采用带噪声语音的相位进行增强语音重构。由于室内环境信噪比较高，所以该方法适用于加性噪声语音的估计，且在室内环境下展示了良好的语音增强性能。虽然谱减法能够抑制背景噪声的影响，但噪声的非平稳性以及分布不均匀等特点，将使得寂静段噪声与发声段噪声并不匹配，使用寂静段噪声估计发生段噪声的方法将导致增强后的语音存在较多的“残留噪声”。为此，Berouti等假设噪声对所有语音频谱具有同等影响，提出将噪声估计功率谱乘以一个修正因子，来削弱噪声对频谱的影响。但现实生活中，不同噪声对语音的高频和低频的影响不同。Lockwood提出非线性的谱减法，即根据频率点设计语音每一帧的修正因子，但这会导致增强后的语音存在较多的“音乐噪声”。为了克服这个缺点，Singh提出多带谱减法，将语音频谱分为多个频带，假设噪声对每个频带内的每一帧语音具有同等的影响，将问题转化为设计每个频带的噪声谱修正系数。然而，针对帕金森病患者语音的语音增强方法尚未见报道。

综上所述，如何提供一种基于多带谱减法的帕金森语音增强方法，就成为本领域技术人员亟待解决的问题。

发明内容

本发明的在基于语音进行帕金森病诊断和病情预测的可行性分析基础上，结合采集到的帕金森病患者语音的特点，提出一种基于多带谱减法的帕金森语音增强方法。该方法采用多带谱减法对帕金森病患者的语音进行语音增强，采用基于最小值统计的噪声估计算法作为帕金森语音增强中的噪声估计算法，不仅能够有效降低语音增强算法的复杂度，而且能够有效降低语音增强后的失真度。

本发明的技术解决方案是：

一种多带谱减法的帕金森语音增强方法，包括以下步骤：

S1：通过非接触式方式对帕金森患者进行语音信号采集；

S2：将步骤S1中采集到的帕金森病患者语音进行采样和短时帧划分；

S3：对每一帧语音进行短时能量计算，并根据短时能量进行自适应频带划分；

S4：计算每个窗长内含噪声语音的短时功率谱、时变的平滑参数及平滑功率谱；

S5：计算偏差修正因子；

S6：寻找搜索功率谱密度的最小值；

S7：计算并更新噪声的功率谱密度；

S8：采用多带谱减法进行语音增强；

S81：将自适应频带划分为N个互不重叠的子带，每个子带独立且同步运行谱减法，假设x(t)与d(t)独立，其中，x(t)为纯净语音信号建模，d(t)为噪声建模，且d(t)服从零均值的高斯分布，用Y_i(W)表示第i个子带上的含噪语音的傅里叶变换、X_i(W)表示第i个子带上的纯净语音的傅里叶变换、D_i(W)表示第i个子带上的噪声信号的傅里叶变换；

S82：对接收到的带噪声信号幅度谱进行平滑处理，并计算平滑后的含噪声语音谱；

S83：根据多带谱减法，利用噪声修正系数对噪声谱进行过减处理，即利用多个子带的权重控制每个子带，并估算第i个子带的纯净语音信号谱；

S84：根据修饰公式对增强后的语音谱进行进一步修饰，所述修饰公式为：

其中，

为最终得到的增强语音，

为第i个子带的纯净语音信号谱，

为平滑后的含噪声语音谱；

噪声修正系数α_i是第i个频率子带的信噪比函数，计算如下：

其中子带的信噪比SNR_i的计算如下：

其中，b_i和e_i是第i个子带上频点的起始点和结束点。

优选地，所述步骤S1中的语音信号采集是用户通过非接触方式，在室内安静环境下采集一段时间的语音信号，并将帕金森病患者的纯净语音信号建模为x(t)，噪声建模为d(t)，且d(t)服从零均值的高斯分布，则采集的帕金森病患者的含噪语音为：

y(t)＝x(t)+d(t)。

优选地，所述步骤S2中采样包括采样频率、采样精度及采样后的含噪语音信号，并设采样频率为F_s kHz，采样精度为B bit，采样后的含噪语音信号为y(m)，其中，m为离散时间度量单位；

所述步骤S2中短时帧划分的帧长为L，相邻的帧为R点重叠。

优选地，所述短时帧划分时所用的窗函数类型为汉宁窗。

优选地，所述步骤S3中计算每一帧语音进行短时能量为计算每一帧中所有样本值的加权平方和，计算公式为：

取所述汉宁窗的窗函数为w(m)，其中0≤m≤L-1，则第λ帧的短时能量计算公式如下：

其中，0≤λ≤J-1，J为总帧数；

所述步骤S3中自适应频带划分是在时域将语音信号分离为多个子带信号，并设置阈值e_low和e_high，取短时能量范围在e_low至e_high的频带，若连续多个帧的短时能量E_λ满足E_λ＞e_high，则以此连续帧的第一帧作为一个划分的起点；若连续多帧的短时能量E_λ满足E_λ＜e_low，则以此连续帧的最后一帧作为一个划分的终点。

优选地，所述步骤S4中短时功率谱的计算方法为：采用步骤S3中经过自适应频带划分所获得的每段频带，针对其中的每一帧，计算出Y(λ,κ)，再计算含噪声语音的短时功率谱|Y(λ,κ)|²，计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，L是帧长，M为各子带内的频点数；

所述步骤S4中时变的平滑参数的计算方法为：

初始化平滑功率谱P(λ，κ)＝0，其中λ＜0，噪声的功率谱密度

时变的平滑参数δ(λ,κ)计算公式为：

其中，P(λ-1，κ)为前一帧的平滑功率谱,λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，M为各子带内的频点数；

为λ-1帧的噪声的功率谱密度；

所述步骤S4中平滑功率谱的计算公式为：

P(λ，κ)＝δ(λ，κ)P(λ-1，κ)+(1-δ(λ，κ))|Y(λ，κ)|²

其中，P(λ,κ)为更新平滑功率谱，δ(λ,κ)为平滑参数，P(λ-1,κ)为前一帧的平滑功率谱。

优选地，所述步骤S5中偏差修正因子B_min(λ,κ)的计算步骤为：

S51：计算

与

其中

是E{P(λ,κ)}的一阶平滑估计，

采用的是E{P²(λ,κ)}的一阶平滑估计，

与

的计算公式分别为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}；P(λ,κ)为更新平滑功率谱；δ(λ,κ)为平滑参数；M为各子带内的频点数；

S52：计算P(λ,κ)方差的估值

S53：计算偏差修正因子B_min(λ,κ)：

其中，D为搜索最小值的窗长或帧长，Q_eq(λ,κ)为：

S54：将B_min(λ,κ)用系数B_c(λ)进行修正，B_c(λ)取值如下：

其中，

优选地，所述步骤S6中搜索功率谱密度的最小值是在长度为D帧的窗内，搜索最小功率谱密度P_min(λ,κ)，其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，且当处理第V帧时，其中V＜D，将最小功率谱密度P_min(λ,κ)更新。

优选地，所述步骤S7中噪声的功率谱密度是当取到更小值时，就更新噪声的功率谱密度

噪声的功率谱密度

的计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，P_min(λ,κ)为最小功率谱密度，B_min(λ,κ)为偏差修正因子。

优选地，所述步骤S81中第i个子带上含噪语音、纯净语音、噪声信号的傅里叶变换公式为：

|Yi(w)|²＝|Xi(w)|²+|Wi(w)|²

其中，Y_i(W)表示第i个子带上的含噪语音、X_i(W)表示第i个子带上的纯净语音、D_i(W)表示第i个子带上的噪声信号；采集的帕金森病患者语音为持续长元音/a/，用户在发声前和发声期间的平稳噪声功率谱相同，并通过基于最小值统计的噪声估计算法估算出噪声的功率谱

所述步骤S82中平滑后的含噪声语音谱的计算公式为：

其中，取P＝2；w_k＝2πk/M是离散频率,k＝0，1，...，M-1；W_j是分配给每一子带的权重，0＜W_j＜1；

所述步骤S83中第i个子带的纯净语音信号谱的估算公式为：

其中，α_i为第i个子带的噪声修正系数，亦称过减因子，且α_i≥1；δ_i为第i个子带的权重因子或减法因子；

为所估计的噪声功率谱；b_i和e_i是第i个子带上频点的起始点和结束点；

减法过程产生的负值按带噪信号谱取下限，则第i个子带的纯净语音信号谱转换为：

其中，β为谱下限参数；

所述步骤S83利用多个子带的权重控制每个子带的权重因子δ_i的计算公式为：

其中，f_i为第i个子带的频率上界，F_S为采样频率。

本发明提供了一种基于多带谱减法的帕金森语音增强方法，其优点主要体现在：

(1)本发明对帕金森病患者在非实验室场景下采集到的语音进行语音增强，同时尽可能降低对语音造成的失真。

(2)本发明在非实验环境下，如室内环境、复杂的噪声环境等，对录制的帕金森病患者语音进行有效的去噪处理，提高帕金森病患者语音质量，最终有效提高帕金森病的诊断或帕金森病情严重程度预测的准确性。

以下便结合实施例附图，对本发明的具体实施方式作进一步的详述，以使本发明技术方案更易于理解、掌握。

附图说明

图1是本发明的流程图；

图2是本发明中自适应频带划分的流程图。

具体实施方式

一种多带谱减法的帕金森语音增强方法，如图1所示：包括以下步骤：

S1：通过非接触式方式对帕金森患者进行语音信号采集，且在室内安静环境下采集一段时间的语音信号，在本实施例中，采集的时间长为4s至8s，并将帕金森病患者的纯净语音信号建模为x(t)，噪声建模为d(t)，且d(t)服从零均值的高斯分布，则采集的帕金森病患者的含噪语音为：

y(t)＝x(t)+d(t)

S2：将步骤S1中采集到的帕金森病患者语音进行采样和短时帧划分；其中，采样包括采样频率、采样精度及采样后的含噪语音信号，并设采样频率为F_skHz，采样精度为Bbit，采样后的含噪语音信号为y(m)，其中，m为离散时间度量单位；短时帧划分的帧长为L，为保持信号的连贯，相邻的帧为R点重叠，且短时帧划分时所用的窗函数类型可根据需要设定为矩形窗、汉宁窗或汉明窗。在本发明的技术方案中，短时帧划分时所用的窗函数类型为汉宁窗。

S3：对每一帧语音进行短时能量计算，并根据短时能量进行自适应频带划分；在每一帧中，计算每一帧语音进行短时能量为计算每一帧中所有样本值的加权平方和，计算公式为：

取所述汉宁窗的窗函数为w(m)，其中0≤m≤L-1，则第λ帧的短时能量计算公式如下：

其中，0≤λ≤J-1，J为总帧数；

另外，自适应频带划分，如图2所示，是在时域将语音信号分离为多个子带信号，并设置阈值e_low和e_high，取短时能量范围在e_low至e_high的频带，若连续多个帧的短时能量E_λ满足E_λ＞e_high，则以此连续帧的第一帧作为一个划分的起点；若连续多帧的短时能量E_λ满足E_λ＜e_low，则以此连续帧的最后一帧作为一个划分的终点。

S4：计算每个窗长内含噪声语音的短时功率谱、时变的平滑参数及平滑功率谱；

其中，短时功率谱的计算方法为：采用步骤S3中经过自适应频带划分所获得的每段频带，针对其中的每一帧，计算出Y(λ,κ)，再计算含噪声语音的短时功率谱|Y(λ,κ)|²，计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，L是帧长，M为各子带内的频点数；

另外，时变的平滑参数的计算方法为：

初始化平滑功率谱P(λ，κ)＝0(λ＜0)，其中λ＜0，噪声的功率谱密度

时变的平滑参数δ(λ,κ)计算公式为：

其中，P(λ-1，κ)为前一帧的平滑功率谱，

为λ-1帧的噪声的功率谱密度；

进一步地，平滑功率谱的计算公式为：

P(λ，κ)＝δ(λ，κ)P(λ-1，κ)+(1-δ(λ，κ))|Y(λ，κ)|²

其中，P(λ,κ)为更新平滑功率谱，δ(λ,κ)为平滑参数，P(λ-1,κ)为前一帧的平滑功率谱。

S5：计算偏差修正因子；偏差修正因子B_min(λ,κ)的计算步骤为：

S51：计算

与

其中

是E{P(λ,κ)}的一阶平滑估计，

采用的是E{P²(λ,κ)}的一阶平滑估计，

与

的计算公式分别为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}；P(λ,κ)为更新平滑功率谱；δ(λ,κ)为平滑参数；M为各子带内的频点数；

S52：计算P(λ,κ)方差的估值

S53：计算偏差修正因子B_min(λ,κ)：

其中，D为搜索最小值的窗长或帧长，Q_eq(λ,κ)为：

S54：将B_min(λ,κ)用系数B_c(λ)进行修正，B_c(λ)取值如下：

其中，

S6：寻找搜索功率谱密度的最小值，是在长度为D帧的窗内，搜索最小功率谱密度P_min(λ,κ)，其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，且当处理第V帧时，其中V＜D，将最小功率谱密度P_min(λ,κ)更新。

S7：计算并更新噪声的功率谱密度；当取到更小值时，就更新噪声的功率谱密度

噪声的功率谱密度

的计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，P_min(λ,κ)为最小功率谱密度，B_min(λ,κ)为偏差修正因子；

S8：采用多带谱减法进行语音增强；

S81：将自适应频带划分为N个互不重叠的子带，每个子带独立且同步运行谱减法，假设x(t)与d(t)独立，其中，x(t)为纯净语音信号建模，d(t)为噪声建模，且d(t)服从零均值的高斯分布，用Y_i(W)表示第i个子带上的含噪语音、X_i(W)表示第i个子带上的纯净语音、D_i(W)表示第i个子带上的噪声信号，其傅里叶变换公式为：

|Y_i(w)|²＝|X_i(w)|²+|D_i(w)|²

其中，Y_i(W)表示第i个子带上的含噪语音、X_i(W)表示第i个子带上的纯净语音、D_i(W)表示第i个子带上的噪声信号；采集的帕金森病患者语音为持续长元音/a/，用户在发声前和发声期间的平稳噪声功率谱相同，并通过基于最小值统计的噪声估计算法估算出噪声的功率谱

S82：对接收到的带噪声信号幅度谱进行平滑处理，并计算平滑后的含噪声语音谱，平滑后的含噪声语音谱的计算公式为：

其中，取P＝2；w_k＝2πk/M是离散频率,k＝0，1，...，M-1；W_j是分配给每一子带的权重，0＜W_j＜1；

S83：根据多带谱减法，利用噪声修正系数对噪声谱进行过减处理，即利用多个子带的权重控制每个子带，并估算第i个子带的纯净语音信号谱；

第i个子带的纯净语音信号谱的估算公式为：

其中，α_i为第i个子带的噪声修正系数，亦称过减因子，且α_i≥1；δ_i为第i个子带的权重因子或减法因子；

为所估计的噪声功率谱；b_i和e_i是第i个子带上频点的起始点和结束点；

减法过程产生的负值按带噪信号谱取下限，则第i个子带的纯净语音信号谱转换为：

其中，β为谱下限参数；

另外，利用多个子带的权重控制每个子带的权重因子δ_i的计算公式为：

其中，f_i为第i个子带的频率上界，F_S为采样频率；

S84：根据修饰公式对增强后的语音谱进行进一步修饰。

修饰公式为：

其中，

为最终得到的增强语音；

噪声修正系数α_i是第i个频率子带的信噪比函数，计算如下：

其中子带的信噪比SNR_i的计算如下：

以下用具体实施例说明多带谱减法的帕金森语音增强方法：

步骤S1：通过非接触式方式对帕金森患者进行语音信号采集；

在本发明的技术方案中，用户通过手机等非接触方式，在室内安静环境下采集4至8秒长的语音信号，其中第1秒不发声，从第2秒开始持续发元音/a/。为了后续步骤表述上的方便，将帕金森病患者的纯净语音信号建模为x(t)，加性噪声建模为d(t)，d(t)服从零均值的高斯分布，则采集的含噪帕金森病患者语音可以表示为：

y(t)＝x(t)+d(t)

步骤S2：将步骤S1中采集到的帕金森病患者语音进行采样和短时帧划分。

采样频率设为F_s kHz，采样精度设为Bbit。采样后的含噪语音信号记为y(m)，其中，m为离散时间度量单位。短时帧划分的帧长为L，为保持信号的连贯，相邻的帧重叠R点。短时帧划分时所用的窗函数类型可根据需要设定矩形窗、汉宁窗、汉明窗等。本发明中，可设采样频率8kHz，采样精度16bit；短时帧划分的帧长L为256(每帧时长32毫秒)，窗函数选汉明窗；帧移R设为128，即帧之间的重叠度为50％。

步骤S3：对每一帧语音进行短时平均能量计算，并根据短时能量进行自适应频带划分；

在每一帧中，计算该帧中所有样本值的加权平方和作为该帧的短时能量。取窗函数w(m)，其中0≤m≤L-1，则第λ帧的短时平均能量计算公式如下：

其中，0≤λ≤J-1，J为总帧数；

根据短时平均能量进行自适应频带划分，具体过程为：

设置阈值e_low和e_high，并根据计算所得的每一帧语音的短时平均能量E_λ，在时域将语音信号分离为多个子带信号，若连续多个帧的短时平均能量E_λ满足E_λ＞e_high，则以此连续帧的第一帧作为一个划分的起点；若连续多帧的短时平均能量E_λ满足E_λ＜e_low，则以此连续帧的最后一帧作为一个划分的终点。由于采集到的语音是在第2秒开始发声，第1秒未发声段语音包含的只有噪声。设第1秒内未发声段的短时平均能量为阈值e_low，由于噪声的随机性，第1秒未发声段语音的短时平均能量不同。因此，对于不同语音样本，e_low的取值不同，即可以根据不同的语音样本自适应调整边界条件或阈值。通过步骤S3，最终实现所采集语音信号的自适应频带划分。取短时能量范围在e_low至e_high的频带，对其中每一帧执行接下来的步骤。

步骤S4：针对自适应划分的频带，计算每个窗长内含噪声语音的短时功率谱、时变的平滑参数及平滑功率谱。

经过自适应频带划分所获得的每段频带，针对其中的每一帧，计算出Y(λ,κ)，再计算含噪声语音的短时功率谱|Y(λ,κ)|²，具体计算过程如下：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，M表示各子带内的频点数，L是帧长。

计算时变的平滑参数。

初始化平滑功率谱P(λ，κ)＝0，其中λ＜0，噪声的功率谱密度

时变的平滑参数δ(λ,κ)计算公式如下：

其中，P(λ-1,κ)为前一帧的平滑功率谱，

为λ-1帧的噪声的功率谱密度。

计算平滑功率谱过程为：更新平滑功率谱P(λ,κ)，计算公式如下：

P(λ，κ)＝δ(λ，κ)P(λ-1，κ)+(1-δ(λ，κ))|Y(λ，κ)|²

其中，δ(λ,κ)为平滑参数，P(λ-1,κ)为前一帧的平滑功率谱。

步骤S5：计算偏差修正因子；计算偏差修正因子B_min(λ,κ)，计算公式如下：

其中，D为搜索最小值的窗长或帧长，本实施例中最小值的窗长取值为96，Q_eq(λ,κ)的计算如下：

其中，

为P(λ,κ)方差的估值，计算公式如下：

其中，

采用的是E{P(λ,κ)}的一阶平滑估计；

采用的是E{P²(λ,κ)}的一阶平滑估计，计算公式如下：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}；P(λ,κ)为更新平滑功率谱；δ(λ,κ)为平滑参数；M为各子带内的频点数；

为了增强在随机噪声环境下噪声估计的稳定性，将B_min(λ,κ)用系数B_c(λ)进行修正，B_c(λ)取值如下：

其中，

步骤S6：寻找搜索功率谱密度的最小值；

在长度为D帧的窗内，搜索最小功率谱密度P_min(λ,κ)。λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，且当处理第V帧时，其中V＜D，将该最小值更新。

步骤S7：计算并更新噪声的功率谱密度；

计算噪声的功率谱密度，若取到更小值，就更新噪声的功率谱密度

计算公式如下：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0,1,...M-1}，P_min(λ,κ)为最小功率谱密度，B_min(λ,κ)为偏差修正因子。

步骤S8：采用多带谱减法进行语音增强；

针对前面划分的N个互不重叠的子带，谱减法在每个子带独立且同步运行。第i(1≤i≤N)个子带的纯净语音信号谱的估计如下：

假设x(t)与d(t)独立，用Y_i(W)表示第i个子带上的含噪语音的傅里叶变换、X_i(W)表示第i个子带上的纯净语音的傅里叶变换、D_i(W)表示第i个子带上的噪声信号的傅里叶变换，且满足：

|Y_i(w)|²＝|X_i(w)|²+|D_i(w)|²

其中，Y_i(W)表示第i个子带上的含噪语音、X_i(W)表示第i个子带上的纯净语音、D_i(W)表示第i个子带上的噪声信号；采集的帕金森病患者语音为持续长元音/a/,可认为用户在发声前和发声期间的平稳噪声功率谱几乎不产生变化，通过基于最小值统计的噪声估计算法可估计出噪声的功率谱

为防止对语音谱信息造成损伤，首先对接收到的带噪声信号幅度谱Y_i(w_k)进行平滑处理。平滑后的带噪声语音谱为

其中，一般取P＝2；w_k＝2πk/M是离散频率,k＝0，1，...，M-1；W_j是分配给每一子带的权重，0＜W_j＜1，本实施例中取W_j＝[0.09,0.25,0.32,0.25,0.09]。

根据多带谱减法，利用噪声修正系数对噪声谱进行过减处理，还可利用多个子带的权重更加灵活地控制每个子带。则第i个子带的纯净语音信号谱的估计如下：

其中，α_i为第i个子带的噪声修正系数，亦称过减因子，且α_i≥1；δ_i为第i个子带的权重因子或减法因子；

为所估计的噪声功率谱；b_i和e_i是第i个子带上频点的起始点和结束点，它可以按子带独立设置以满足不同的噪声抑制需求。减法过程产生的负值按带噪信号谱取下限，则上式转换为：

其中，谱下限参数β的值为0.002。为进一步消除“音乐噪声”，可对增强后的语音谱进行进一步修饰，公式如下：

其中，

为最终得到的增强语音。

噪声修正系数α_i是第i个频率子带的信噪比函数，根据经验确定其计算如下：

其中子带的信噪比SNR_i的计算如下：

权重因子(减法因子)δ_i的计算公式如下：

其中，f_i为第i个子带的频率上界，F_S为采样频率。

本发明提供的基于多带谱减法的帕金森语音增强方法，主要基于采集患者的持续长元音/a/，持续长元音/a/除了语音的起始和终止，发声期间并没有“寂静段”，通过计算语音的短时能量，并根据短时能量划分频带，即自适应频带划分，该方法可有效降低语音增强算法的复杂度。此外，对于一般语音而言，在语音的“寂静段”，单个频带的带噪语音信号的功率也会衰减到噪声的平均功率，基于最小值统计的噪声估计算法无需对语音的“寂静段”和“发声段”的噪声分别进行估计，因而该噪声估计算法对低频带的语音增强具有更好的质量。因此，采用基于最小值统计的噪声估计算法作为帕金森语音增强中的噪声估计算法，不仅能够有效降低语音增强算法的复杂度，而且能够有效降低语音增强后的失真度。

应该注意的是，上述实施例对本发明进行说明而不是对本发明进行限制，并且本领域技术人员在不脱离所附权利要求的范围的情况下可设计出替换实施例。

Claims (10)

1.一种多带谱减法的帕金森语音增强方法，其特征在于，包括以下步骤：

S1：通过非接触式方式对帕金森患者进行语音信号采集；

S2：将步骤S1中采集到的帕金森病患者语音进行采样和短时帧划分；

S3：对每一帧语音进行短时能量计算，并根据短时能量进行自适应频带划分；

S4：计算每个窗长内含噪声语音的短时功率谱、时变的平滑参数及平滑功率谱；

S5：计算偏差修正因子；

S6：寻找搜索功率谱密度的最小值；

S7：计算并更新噪声的功率谱密度；

S8：采用多带谱减法进行语音增强；

S81：将自适应频带划分为N个互不重叠的子带，每个子带独立且同步运行谱减法，假设x(t)与d(t)独立，其中，x(t)为纯净语音信号建模，d(t)为噪声建模，且d(t)服从零均值的高斯分布，用Y_i(W)表示第i个子带上的含噪语音的傅里叶变换、X_i(W)表示第i个子带上的纯净语音的傅里叶变换、D_i(W)表示第i个子带上的噪声信号的傅里叶变换；

S82：对接收到的带噪声信号幅度谱进行平滑处理，并计算平滑后的含噪声语音谱；

S83：根据多带谱减法，利用噪声修正系数对噪声谱进行过减处理，即利用多个子带的权重控制每个子带，并估算第i个子带的纯净语音信号谱；

S84：根据修饰公式对增强后的语音谱进行进一步修饰，所述修饰公式为：

其中，

为最终得到的增强语音，

为第i个子带的纯净语音信号谱，

为平滑后的含噪声语音谱；

噪声修正系数α_i是第i个频率子带的信噪比函数，计算如下：

其中子带的信噪比SNR_i的计算如下：

其中，b_i和e_i是第i个子带上频点的起始点和结束点。

2.根据权利要求1所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S1中的语音信号采集是用户通过非接触方式，在室内安静环境下采集一段时间的语音信号，并将帕金森病患者的纯净语音信号建模为x(t)，噪声建模为d(t)，且d(t)服从零均值的高斯分布，则采集的帕金森病患者的含噪语音为：

y(t)＝x(t)+d(t)。

3.根据权利要求2所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S2中采样包括采样频率、采样精度及采样后的含噪语音信号，并设采样频率为F_s kHz，采样精度为B bit，采样后的含噪语音信号为y(m)，其中，m为离散时间度量单位；

所述步骤S2中短时帧划分的帧长为L，相邻的帧为R点重叠。

4.根据权利要求3所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述短时帧划分时所用的窗函数类型为汉宁窗。

5.根据权利要求4所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S3中计算每一帧语音进行短时能量为计算每一帧中所有样本值的加权平方和，计算公式为：

取所述汉宁窗的窗函数为w(m)，其中0≤m≤L-1，则第λ帧的短时能量计算公式如下：

其中，0≤λ≤J-1，J为总帧数；

所述步骤S3中自适应频带划分是在时域将语音信号分离为多个子带信号，并设置阈值e_low和e_high，取短时能量范围在e_low至e_high的频带，若连续多个帧的短时能量E_λ满足E_λ＞e_high，则以此连续帧的第一帧作为一个划分的起点；若连续多帧的短时能量E_λ满足E_λ＜e_low，则以此连续帧的最后一帧作为一个划分的终点。

6.根据权利要求5所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S4中短时功率谱的计算方法为：采用步骤S3中经过自适应频带划分所获得的每段频带，针对其中的每一帧，计算出Y(λ,κ)，再计算含噪声语音的短时功率谱|Y(λ,κ)|²，计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，L是帧长，M为各子带内的频点数；

所述步骤S4中时变的平滑参数的计算方法为：

初始化平滑功率谱P(λ，κ)＝0，其中λ＜0，噪声的功率谱密度

时变的平滑参数δ(λ,κ)计算公式为：

其中，P(λ-1，κ)为前一帧的平滑功率谱,λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，M为各子带内的频点数；

为λ-1帧的噪声的功率谱密度；

所述步骤S4中平滑功率谱的计算公式为：

P(λ，κ)＝δ(λ，κ)P(λ-1，κ)+(1-δ(λ，κ))|Y(λ，κ)|²

其中，P(λ,κ)为更新平滑功率谱，δ(λ,κ)为平滑参数，P(λ-1,κ)为前一帧的平滑功率谱。

7.根据权利要求6所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S5中偏差修正因子B_min(λ,κ)的计算步骤为：

S51：计算

与

其中

是E{P(λ,κ)}的一阶平滑估计，

采用的是E{P²(λ,κ)}的一阶平滑估计，

与

的计算公式分别为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}；P(λ,κ)为更新平滑功率谱；δ(λ,κ)为平滑参数；M为各子带内的频点数；

S52：计算P(λ,κ)方差的估值

S53：计算偏差修正因子B_min(λ,κ)：

其中，D为搜索最小值的窗长或帧长，Q_eq(λ,κ)为：

S54：将B_min(λ,κ)用系数B_c(λ)进行修正，B_c(λ)取值如下：

其中，

8.根据权利要求7所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S6中搜索功率谱密度的最小值是在长度为D帧的窗内，搜索最小功率谱密度P_min(λ,κ)，其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，且当处理第V帧时，其中V＜D，将最小功率谱密度P_min(λ,κ)更新。

9.根据权利要求8所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S7中噪声的功率谱密度是当取到更小值时，就更新噪声的功率谱密度

噪声的功率谱密度

的计算公式为：

其中，λ为帧数标记，λ∈Z；κ表示频率点，κ∈{0，1，...M-1}，P_min(λ,κ)为最小功率谱密度，B_min(λ,κ)为偏差修正因子。

10.根据权利要求9所述一种基于多带谱减法的帕金森语音增强方法，其特征在于：所述步骤S81中第i个子带上含噪语音、纯净语音、噪声信号的傅里叶变换公式为：

|Y_i(w)|²＝|X_i(w)|²+|D_i(w)|²

其中，Y_i(W)表示第i个子带上的含噪语音、X_i(W)表示第i个子带上的纯净语音、D_i(W)表示第i个子带上的噪声信号；采集的帕金森病患者语音为持续长元音/a/，用户在发声前和发声期间的平稳噪声功率谱相同，并通过基于最小值统计的噪声估计算法估算出噪声的功率谱

所述步骤S82中平滑后的含噪声语音谱的计算公式为：

其中，取P＝2；w_k＝2πk/M是离散频率,k＝0，1，...，M-1；W_j是分配给每一子带的权重，0＜W_j＜1；

所述步骤S83中第i个子带的纯净语音信号谱的估算公式为：

其中，α_i为第i个子带的噪声修正系数，亦称过减因子，且α_i≥1；δ_i为第i个子带的权重因子或减法因子；

为所估计的噪声功率谱；b_i和e_i是第i个子带上频点的起始点和结束点；

减法过程产生的负值按带噪信号谱取下限，则第i个子带的纯净语音信号谱转换为：

其中，β为谱下限参数；

所述步骤S83利用多个子带的权重控制每个子带的权重因子δ_i的计算公式为：

其中，f_i为第i个子带的频率上界，F_S为采样频率。

Images