本发明涉及语音信号处理和机器学习技术，尤其涉及一种基于语音声学特征的声带异常检测方法。
背景技术：
语音是人类交流的重要手段之一，语音交流障碍严重影响人们的正常生活。声带病变是导致语音交流障碍的主要因素之一。声带异常情况的准确诊断是发音康复治疗的前提，在临床医学中尤为重要。传统的声带异常检测方法是采用喉镜或者电声门图仪。基于喉镜的检测方法具有入侵性，将喉镜插入测试者喉部，需要测试者极力配合，会给测试者带来难以接受的痛苦。另外，做喉镜检测时测试者没法正常发声，无法检测声带振动情况，具有一定的局限性。基于电声门图仪的检测方法需要将两个电极板紧贴在测试者颈部，也会给测试者带来不适感，不适合颈部肥胖者，特别是老人、女性和小孩。
技术实现要素：
本发明的主要目的在于克服现有技术的缺点与不足，提供一种基于语音声学特征的声带异常检测方法，本发明的方法通过提取测试者语音的声学特征，再采用gmm分类器进行判决，从而诊断测试者的声带是否异常。
为了达到上述目的，本发明采用以下技术方案：
一种基于语音声学特征的声带异常检测方法，包括下述步骤：
s1、读取语音：读入语音样本，得到语音序列s(n)；
s2、预处理：对语音数据进行预加重、分帧、加窗和去静音处理，得到语音帧st(n)，1≤t≤t，其中t表示语音帧数；
s3、提取声学特征：从语音帧st(n)中提取梅尔频率倒谱系数(mel-frequencycepstralcoefficients,mfcc)，基频(f0)，基频微扰(jitter)，振幅微扰(shimmer)，谐噪比(harmonynoiseratio,hnr)，得到特征矩阵f＝[mfcc,f0,jitter,shimmer,hnr]；
s4、训练高斯混合模型(gaussianmixturemodel,gmm)：以特征矩阵f作为输入，采用期望最大化(expectationmaximization,em)算法分别训练代表声带异常和声带正常的高斯混合模型θa和θn；
s5、声带异常判决：将测试语音的特征矩阵f分别输入高斯混合模型θa和θn，得到相应的输出概率p(f|θa)和p(f|θn)，如果p(f|θa)>p(f|θn)，则测试语音的说话人声带异常，否则正常。
优选的，步骤s2中预处理包括如下步骤：
s2.1、预加重：利用数字滤波器h(n)对语音数据进行滤波处理，h(n)的z变换h(z)表示为：
h(z)＝1-μz-1，
其中μ取0.98；
s2.2、分帧：将预加重后的语音数据进行分帧处理，设置语音帧的帧长l，帧移s；
s2.3、加窗：将每帧语音st(n)与窗函数相乘，其中窗函数为汉明窗ω(n)：
其中n表示一帧语音的采样点数，且n＝l×fs，其中fs表示语音采样频率；
s2.4、去静音：
s2.4.1、计算第t帧语音st(n)的能量et：
得到语音流的能量特征矢量e＝[e1,e2,…,et]，其中t为总帧数；
s2.4.2、用固定的能量门限判断静音与语音有很大的局限性，因为各种环境下的语音能量相差很大，但语音与静音之间能量大小的相对关系是不变的，所以定义每帧语音的能量门限te：
te＝min(e)+0.3×[mean(e)-min(e)]
其中，min(e)是各帧能量的最小值，mean(e)是各帧能量的平均值；
s2.4.3、将每帧语音的能量与能量门限比较，低于能量门限的帧为静音帧，否则为语音帧，舍弃静音帧，保留语音帧。
优选的，步骤s3提取声学特征的步骤包括：
s3.1、提取mfcc，具体步骤如下：
s3.1.1、对第t帧语音st(n)做离散傅立叶变换(discretefouriertransformation,dft)得到线性频谱st(k)：
s3.1.2、将上述线性频谱st(k)通过梅尔滤波器组滤波得到梅尔频谱，再进行对数运算得到对数频谱st(m)，其中梅尔滤波器组为若干个带通滤波器hm(k)，0≤m≤m，m为滤波器的个数，每个滤波器具有三角形滤波特性，中心频率为f(m)，m较小时相邻f(m)的间隔较小，m增加时相邻f(m)的间隔变大，每个带通滤波器的传递函数为：
其中，f(m)定义如下：
其中fl和fh分别为滤波器频率应用范围的最低和最高频率，n为一帧语音的采样点数，fs为采样频率，b-1为b的逆函数：
b-1(b)＝700(eb/1125-1)，
因此由线性谱st(k)到对数谱st(m)的函数式为：
s3.1.3、将上述对数频谱st(m)经过离散余弦变换(discretecosinetransformation,dct)变换到倒谱域，得到第t帧语音的mfcc：
s3.2、提取f0：采用自相关法提取各语音帧的基频，具体步骤如下：
s3.2.1、对语音帧st(n)做中心削波处理，得到yt(n)：
其中cl等于语音帧最大幅度的0.68倍；
s3.2.2、计算yt(n)的自相关函数ry(i)：
其中n为语音帧采样点数；
检测ry(i)的峰值点，得到峰值点集合vp，再计算vp差分的均值，得到平均基音周期t0：
其中i为峰值点个数，t0的倒数就是基频f0：
s3.3、提取jitter：语音相邻周期存在微小变化，声带异常程度不同所引起的基频变化就不同，基频微扰jitter为：
其中f0为基频，j为一个语音帧中的基音周期个数；
s3.4、提取shimmer：振幅微扰是指语音帧中各周期的语音振幅的相对变化。振幅微扰shimmer为：
其中a表示语音帧中某个周期的语音振幅，j表示振幅个数；
s3.5、提取hnr：提取语音帧谐噪比的步骤如下：
计算语音帧st(n)的自相关函数rn(i)：
其中n表示语音帧的采样点数，设rh(i)和rn(i)分别表示语音信号的谐波相关性和噪声相关性，由于谐波成分和噪声成分不相关且噪声成分没有相关性，故有：
当i＝0时，
rn(0)＝rh(0)+rn(0)，
当i＝t0时，
rn(t0)＝rh(t0)＝rh(0)，
谐噪比hnr为：
s3.6、将语音帧的mfcc、f0、jitter、shimmer和hnr依次拼接，构成特征矩阵f＝[mfcc,f0,jitter,shimmer,hnr]。
优选的，步骤s4训练高斯混合模型的步骤包括：
s4.1、包含g个高斯混合分量的d维高斯混合模型θ表示为：θ＝{wi,ui,σi}，其中wi表示混合权重系数，μi表示均值，σi表示协方差矩阵，第t帧语音的特征矩阵ft输入θ的输出概率为：
其中bi(ft)表示第i个高斯分量的概率分布，表示为：
其中d和ft分别表示特征维数和第t帧语音的特征矩阵，所有t帧语音的特征矩阵f＝[f1,…,ft,…,ft]t，t表示矩阵转置运算，f输入θ的输出概率为：
s4.2、求期望：计算函数q(θ,θ*)：
q(θ,θ*)＝e{lnp[(f,i)|θ]}，
整理得：
其中θ*表示迭代更新后的通用背景模型，ln(·)表示自然对数函数，p(i|ft,θ*)表示更新后的ubm的第i个高斯对特征矩阵ft的输出概率，p(i|ft,θ)表示更新前的ubm的第i个高斯对特征矩阵ft的输出概率，p(i|ft,θ*)与p(i|ft,θ)的表达式相同，只需用当前迭代后的ubm参数代替当前迭代前的ubm参数即可，根据贝叶斯公式，特征矩阵ft在第i个高斯的概率为：
其中pi(ft)和pj(ft)分别为：
s4.3、最大化：根据q(θ,θ*)函数，估计θ＝{wi，ui，σi}：
s4.4、em算法迭代gmm参数：当似然函数值达到最大时停止迭代，即当前迭代的输出概率p(f|θ*)值相对上次迭代时的输出概率p(f|θ)值增幅小于设定的阈值10-4，则停止迭代，得到最终的模型参数：
混合权重系数：
均值矢量：
协方差矩阵：
优选的，步骤s5中声带异常判决，对测试语音依次进行s1、s2和s3步骤的处理，得到测试语音的特征矩阵f并分别输入采用s4步骤得到的代表声带异常和正常的高斯混合模型θa和θn，得到相应的输出概率p(f|θa)和p(f|θn)，如果p(f|θa)>p(f|θn)，则测试语音的说话人声带异常，否则正常。
本发明与现有技术相比，具有如下优点和有益效果：
1、本发明采用测试者发出的语音作为分析对象，并从测试语音中提取能有效反应声带情况的多组声学特征作为高斯混合模型的输入，有效区分了声带正常和异常的语音，从而诊断测试者声带是否异常。
2、本发明方法具有非入侵、友好便利、成本低等优点，可以适合任意人群的声带异常检测，有效克服了喉镜、电声门图仪等现有检测技术的不足。
附图说明
图1为本发明的实施例的流程图。
具体实施方式
下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。
实施例
如图1所示为本发明的实施例的流程图，具体步骤如下：
首先，在步骤101中，读取采集到语音样本。
接着，在步骤102中，对读入的语音数据进行预处理，具体步骤如下：
1、预加重：利用数字滤波器h(n)对语音数据进行滤波处理，h(n)的z变换h(z)表示为：
h(z)＝1-μz-1，
其中μ取0.98；
2、分帧：将预加重后的语音数据进行分帧处理，帧长l＝30ms，帧移s＝15ms；
3、加窗：将每帧语音st(n)与窗函数相乘，其中窗函数为汉明窗ω(n)：
其中n表示一帧语音的采样点数，且n＝l×fs，其中fs表示语音采样频率；
4、去静音：
4.1、计算第t帧语音st(n)的能量et：
得到语音流的能量特征矢量e＝[e1,e2,…,et]，其中t为总帧数；
4.2、用固定的能量门限判断静音与语音有很大的局限性，因为各种环境下的语音能量相差很大，但语音与静音之间能量大小的相对关系是不变的，所以定义每帧语音的能量门限te：
te＝min(e)+0.3×[mean(e)-min(e)]
其中，min(e)是各帧能量的最小值，mean(e)是各帧能量的平均值；
4.3、将每帧语音的能量与能量门限比较，低于能量门限的帧为静音帧，否则为语音帧，舍弃静音帧，保留语音帧。
接着，在步骤103中，提取各帧语音的声学特征，包括梅尔频率倒谱系数mfcc、基频f0、基频微扰jitter、振幅微扰shimmer和谐噪比hnr，并依次将上述特征拼接构造特征矩阵f，具体步骤如下：
1、提取mfcc：具体步骤如下：
1.1、对第t帧语音st(n)做离散傅立叶变换(discretefouriertransformation,dft)得到线性频谱st(k)：
1.2、将上述线性频谱st(k)通过梅尔滤波器组滤波得到梅尔频谱，再进行对数运算得到对数频谱st(m)，其中梅尔滤波器组为若干个带通滤波器hm(k)，0≤m≤m，m为滤波器的个数，每个滤波器具有三角形滤波特性，中心频率为f(m)，m较小时相邻f(m)的间隔较小，m增加时相邻f(m)的间隔变大，每个带通滤波器的传递函数为：
其中，f(m)定义如下：
其中，fl和fh分别为滤波器频率应用范围的最低和最高频率，n为一帧语音的采样点数，fs为采样频率，b-1为b的逆函数：
b-1(b)＝700(eb/1125-1)，
因此由线性谱st(k)到对数谱st(m)的函数式为：
1.3、将上述对数频谱st(m)经过离散余弦变换(discretecosinetransformation,dct)变换到倒谱域，得到第t帧语音的mfcc：
2、提取f0：采用自相关法提取各语音帧的基频，具体步骤如下：
2.1、对语音帧st(n)做中心削波处理，得到yt(n)：
其中cl等于语音帧最大幅度的0.68倍；
2.2、计算yt(n)的自相关函数ry(i)：
其中n为语音帧采样点数；
检测ry(i)的峰值点，得到峰值点集合vp，再计算vp差分的均值，得到平均基音周期t0：
其中i为峰值点个数，t0的倒数就是基频f0：
3、提取jitter：语音相邻周期存在微小变化，声带异常程度不同所引起的基频变化就不同，基频微扰jitter为：
其中f0为基频，j为一个语音帧中的基音周期个数；
4、提取shimmer：振幅微扰是指语音帧中各周期的语音振幅的相对变化。振幅微扰shimmer为：
其中a表示语音帧中某个周期的语音振幅，j表示振幅个数；
5、提取hnr：提取语音帧谐噪比的步骤如下：
计算语音帧st(n)的自相关函数rn(i)：
其中n表示语音帧的采样点数，设rh(i)和rn(i)分别表示语音信号的谐波相关性和噪声相关性，由于谐波成分和噪声成分不相关且噪声成分没有相关性，故有：
当i＝0时，
rn(0)＝rh(0)+rn(0)，
当i＝t0时，
rn(t0)＝rh(t0)＝rh(0)，
谐噪比hnr为：
将语音帧的mfcc、f0、jitter、shimmer和hnr依次拼接，构成特征矩阵f＝[mfcc,f0,jitter,shimmer,hnr]。
接着，在步骤104中，采用em算法训练分别代表声带异常和正常的高斯混合模型θa和θn，具体步骤如下：
1、包含g个高斯混合分量的d维高斯混合模型θ表示为：θ＝{wi,ui,σi}，其中wi表示混合权重系数，μi表示均值，σi表示协方差矩阵，第t帧语音的特征矩阵ft输入θ的输出概率为：
其中bi(ft)表示第i个高斯分量的概率分布，表示为：
其中d和ft分别表示特征维数和第t帧语音的特征矩阵，所有t帧语音的特征矩阵f＝[f1,…,ft,…,ft]t，t表示矩阵转置运算，f输入θ的输出概率为：
2、求期望：计算函数q(θ,θ*)：
q(θ,θ*)＝e{lnp[(f,i)|θ]}，
整理得：
其中θ*表示迭代更新后的通用背景模型，ln(·)表示自然对数函数，p(i|ft,θ*)表示更新后的ubm的第i个高斯对特征矩阵ft的输出概率，p(i|ft,θ)表示更新前的ubm的第i个高斯对特征矩阵ft的输出概率，p(i|ft,θ*)与p(i|ft,θ)的表达式相同，只需用当前迭代后的ubm参数代替当前迭代前的ubm参数即可，根据贝叶斯公式，特征矩阵ft在第i个高斯的概率为：
其中pi(ft)和pj(ft)分别为：
3、最大化：根据q(θ,θ*)函数，估计θ＝{wi，ui，σi}：
4、em算法迭代gmm参数：当似然函数值达到最大时停止迭代，即当前迭代的输出概率p(f|θ*)值相对上次迭代时的输出概率p(f|θ)值增幅小于设定的阈值10-4，则停止迭代，得到最终的模型参数：
混合权重系数：
均值矢量：
协方差矩阵：
最后，在步骤105中，对测试语音依次进行101、102和103步骤的处理，得到测试语音的特征矩阵f并分别输入采用104步骤得到的代表声带异常和正常的高斯混合模型θa和θn，得到相应的输出概率p(f|θa)和p(f|θn)，如果p(f|θa)>p(f|θn)，则测试语音的说话人声带异常，否则正常。
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。