曾發(fā)林， 蔡嘉偉， 孫蘇民

(1.江西省汽車(chē)噪聲與振動(dòng)重點(diǎn)實(shí)驗(yàn)室，江西 南昌 330013； 2.江蘇大學(xué) 汽車(chē)工程研究院，江蘇 鎮(zhèn)江 212013)

0 引言

汽車(chē)聲品質(zhì)體現(xiàn)了人對(duì)噪聲的主觀感受，關(guān)于聲品質(zhì)的研究主要采用主觀評(píng)價(jià)試驗(yàn)，它能直接反映出聲品質(zhì)，但是評(píng)價(jià)過(guò)程比較復(fù)雜，費(fèi)時(shí)費(fèi)力。于是，很多學(xué)者基于心理聲學(xué)客觀參數(shù)建立了預(yù)測(cè)模型，但汽車(chē)噪聲多屬于非穩(wěn)態(tài)信號(hào)，只研究時(shí)域或頻域特征不能準(zhǔn)確提取非穩(wěn)態(tài)信號(hào)特征。此時(shí)，應(yīng)該從時(shí)頻聯(lián)合域上研究并提取信號(hào)特征，然而基于時(shí)頻信號(hào)處理方法(如小波分解、EMD等)建立的聲品質(zhì)預(yù)測(cè)模型并沒(méi)有得到廣泛應(yīng)用[1]。

針對(duì)非穩(wěn)態(tài)排氣噪聲聲品質(zhì)預(yù)測(cè)問(wèn)題，筆者引入CEEMD方法對(duì)非穩(wěn)態(tài)信號(hào)進(jìn)行分解，計(jì)算IMF分量樣本熵值，并使用PCA方法進(jìn)行降維處理，得到新參量矩陣，建立了基于SQP-CSP新參量的GA-BP預(yù)測(cè)模型。為了對(duì)比此模型，建立了基于Zwicker時(shí)變算法[2]計(jì)算非穩(wěn)態(tài)信號(hào)的心理聲學(xué)客觀參量的GA-BP預(yù)測(cè)模型。結(jié)果表明，新參量建立的模型預(yù)測(cè)精度較好。

1 排氣噪聲采集及主客觀評(píng)價(jià)分析

1.1 排氣噪聲樣本采集

根據(jù)GB/T 14365—2017《聲學(xué) 機(jī)動(dòng)車(chē)輛定置噪聲聲壓級(jí)測(cè)量方法》，采用LMS測(cè)試并采集了10款某國(guó)產(chǎn)車(chē)1檔急加速(1G WOT)和2檔急加速(2G WOT)排氣噪聲信號(hào)，布置如圖1所示。采集數(shù)據(jù)時(shí)參數(shù)設(shè)置如下：采樣頻率為44 100 Hz，分辨率為1 Hz，采樣時(shí)間為15 s。

圖1 試驗(yàn)示意圖Figure 1 Sketch map of the vehicle test

1.2 排氣噪聲的主觀評(píng)價(jià)及客觀參數(shù)

本次主觀試驗(yàn)截取加速工況下的特征信號(hào)10 s進(jìn)行成對(duì)比較法評(píng)價(jià)。樣本兩兩配對(duì)比較，聽(tīng)音者認(rèn)為效果好的樣本記1分，差的不計(jì)分，最終每個(gè)樣本都會(huì)有一個(gè)確定的數(shù)值用于評(píng)價(jià)其聲品質(zhì)，并用于后期模型的建立。

參與主觀評(píng)價(jià)的人員共有40名[3]，分別來(lái)自某大學(xué)車(chē)輛相關(guān)專(zhuān)業(yè)的在讀研究生、研究院的相關(guān)工作者以及普通人員。其中男性有22名，女性有18名，年齡在24～46歲。對(duì)樣本重合度[4]和一致性系數(shù)[5]進(jìn)行計(jì)算，剔除了4名評(píng)價(jià)人員的數(shù)據(jù)，最終得到樣本的平均重合度為0.779，平均一致性系數(shù)為0.916。其中一致性系數(shù)采用Kendall法[6]計(jì)算。

心理聲學(xué)參數(shù)能較好地反映人的主觀感受，所以基于Zwicker時(shí)變算法在ArtermiS中計(jì)算非穩(wěn)態(tài)信號(hào)的心理聲學(xué)客觀參量。根據(jù)研究經(jīng)驗(yàn)[7-8]，心理聲學(xué)客觀參量采用Zwicker響度、尖銳度、粗糙度、波動(dòng)度、峭度及A計(jì)權(quán)聲壓級(jí)。樣本的心理聲學(xué)客觀參量及主觀試驗(yàn)得到的統(tǒng)計(jì)結(jié)果如表1所示，其中滿(mǎn)意度歸一化公式為：

(1)

表1 心理聲學(xué)客觀參量值及滿(mǎn)意度值Table 1 Psychoacoustic objective parameter value and satisfaction value

1.3 相關(guān)性分析

為了分析主觀滿(mǎn)意度與心理聲學(xué)客觀參量之間的聯(lián)系，對(duì)主觀評(píng)價(jià)試驗(yàn)滿(mǎn)意度得分和心理聲學(xué)客觀參量采用SPSS軟件進(jìn)行相關(guān)分析，并采用spearman秩相關(guān)雙尾進(jìn)行相關(guān)分析來(lái)消除聲品質(zhì)客觀參量中存在的極端值，spearman秩相關(guān)計(jì)算公式為：

(2)

式中：Ui和Vi為兩變量的秩，主要是為了將定距型變量變?yōu)榉嵌ň嘈?，從而減小極端值對(duì)結(jié)果的影響；n為樣本數(shù)；r為spearman秩相關(guān)系數(shù)。

相關(guān)性分析的結(jié)果如表2所示。從表2中可以看出，主觀評(píng)價(jià)滿(mǎn)意度與響度的相關(guān)性系數(shù)大于0.7，相關(guān)顯著性小于0.01，可見(jiàn)主觀滿(mǎn)意度與響度具有一定的線性關(guān)系；與尖銳度的相關(guān)系數(shù)為0.571，線性相關(guān)性較弱；粗糙度、波動(dòng)度、峭度與主觀評(píng)價(jià)滿(mǎn)意度相關(guān)系數(shù)都較低，說(shuō)明相關(guān)性較差。同時(shí)也可以看出，響度、尖銳度、波動(dòng)度與主觀滿(mǎn)意度成負(fù)相關(guān)，粗糙度、峭度與主觀滿(mǎn)意度成正相關(guān)。因此，考慮到非穩(wěn)態(tài)排氣噪聲聲品質(zhì)與主觀評(píng)價(jià)之間關(guān)系的復(fù)雜性，引入具有較強(qiáng)的非線性映射能力的神經(jīng)網(wǎng)絡(luò)對(duì)聲品質(zhì)進(jìn)行預(yù)測(cè)研究。

表2 主觀評(píng)價(jià)與心理聲學(xué)客觀參量的相關(guān)系數(shù)Table 2 Correlation coefficient between subjective evaluation and objective parameters of psychoacoustics

2 基于樣本熵計(jì)算CEEMD的特征提取

2.1 CEEMD

CEEMD是在EEMD(ensemble empirical mode decomposition)基礎(chǔ)上進(jìn)行改進(jìn)后得出的一種處理非平穩(wěn)信號(hào)的方法。其原理是根據(jù)信號(hào)本身的時(shí)間尺度進(jìn)行分解，然后得到一系列不同特征尺度的數(shù)據(jù)序列IMF分量[9]。本征模態(tài)函數(shù)描述的是單分量物理意義的振動(dòng)信號(hào)。Yeh等[10]在EEMD的基礎(chǔ)上提出了CEEMD方法，將白噪聲以正負(fù)成對(duì)的形式添加到原始信號(hào)當(dāng)中，從而降低平均次數(shù)，計(jì)算效率得到明顯提高，也有效解決了模態(tài)混疊問(wèn)題，CEEMD分解步驟如下。

步驟1將n組白噪聲以正負(fù)成對(duì)的形式添加到原始信號(hào)當(dāng)中，形成兩組集合:

(3)

式中：s(t)為原始信號(hào)；n(t)為白噪聲；m1、m2分別為添加了正負(fù)成對(duì)白噪聲的合成信號(hào)。

步驟2對(duì)2n組集合信號(hào)進(jìn)行EEMD分解，設(shè)置集合平均次數(shù)，由此每一個(gè)信號(hào)都能得到一系列的IMF分量，記imfij(t)為第i個(gè)信號(hào)的第j個(gè)分量。

步驟3求解線性平均值：

(4)

式中：imfj(t)為2n組信號(hào)在分解后第j個(gè)分量的線性平均值。

經(jīng)過(guò)CEEMD分解后，原始信號(hào)可以用各層分量和殘差之和來(lái)表示：

(5)

式中：m為IMF分量的個(gè)數(shù)；s(t)為原始信號(hào);rm(t)為殘余分量。

2.2 樣本熵

樣本熵(sample entropy，SE)的物理意義是通過(guò)度量信號(hào)中產(chǎn)生新模式的概率來(lái)判斷時(shí)間序列的復(fù)雜性，并且樣本熵對(duì)數(shù)據(jù)的長(zhǎng)度沒(méi)有要求，也具備更好的一致性，因此適合對(duì)非穩(wěn)態(tài)信號(hào)的特征進(jìn)行計(jì)算。

假設(shè)一個(gè)由N個(gè)數(shù)據(jù)組成時(shí)間序列{x(n)}={x(1),x(2),…,x(N)}，樣本熵的計(jì)算方法具體如下。

步驟1按序號(hào)形成一組m維的向量組Xm(1),…,Xm(N-m+1),其中Xm(i)={x(i),x(i+1),…,x(i+m-1)},1≤i≤N-m+1。

步驟2定義兩個(gè)向量Xm(i)和Xm(j)之間的距離d[Xm(i),Xm(j)]為向量組中對(duì)應(yīng)元素差值的絕對(duì)值中的最大值。即：

d[Xm(i)，Xm(j)]=

(6)

步驟3對(duì)于給定的Xm(i)，計(jì)算Xm(i)和Xm(j)之間的距離小于等于r的數(shù)目，記為Bi(即模板匹配數(shù))，并計(jì)算Bi和距離總數(shù)N-m+1的比值，記為：

(7)

步驟4定義Bm(r)為：

(8)

步驟5將維數(shù)m增加1至(m+1)維，對(duì)(m+1)維向量重復(fù)步驟1～4，得B(m+1)(r)。

步驟6該信號(hào)序列的樣本熵值為：

(9)

式中：m為重構(gòu)維數(shù)；r為相似容限；N為數(shù)據(jù)長(zhǎng)度。m一般選擇1或2，實(shí)際使用中學(xué)者都優(yōu)先選擇m=2，r選擇原始數(shù)據(jù)標(biāo)準(zhǔn)偏差的0.1～0.25倍。

2.3 主成分分析

主成分分析(PCA)方法是一種數(shù)學(xué)方法，其主要思想是將一組相關(guān)的向量轉(zhuǎn)換成另一組不相關(guān)的向量，從而達(dá)到降維的目的，新的向量按方差進(jìn)行降序排列，將具有最大方差的向量稱(chēng)為第一主成分，第二大方差稱(chēng)為第二主成分，第K個(gè)向量稱(chēng)為第K個(gè)主成分。PCA計(jì)算步驟如下。

步驟1假設(shè)數(shù)據(jù)矩陣為X=(Xij)m×n,其中，i=1，2，…，m，其中m是樣本數(shù)量；j=1，2，…，n，其中n是特征維數(shù)；Xij是第i個(gè)樣本的第j個(gè)特征。

(10)

通過(guò)Z-score方法對(duì)矩陣X進(jìn)行標(biāo)準(zhǔn)化處理，從而得到標(biāo)準(zhǔn)化矩陣Z:

(11)

步驟2計(jì)算標(biāo)準(zhǔn)化矩陣Z的相關(guān)系數(shù)矩陣R。

(12)

式中：j=1，2，…，n；k=1，2，…，n。

步驟3求解相關(guān)系數(shù)矩陣R的特征值λ=(λ1,λ2,…,λn)和特征向量αk=(αk1,αk2,…,αkn)，k=1，2，…，n，并通過(guò)累計(jì)貢獻(xiàn)率法確定前p個(gè)主成分，累計(jì)貢獻(xiàn)率法公式為:

(13)

步驟4確定主成分得分，第l個(gè)主成分(l=1，2，…，p)計(jì)算如下:

Fl=αk1Z1+αk2Z2+…+αknZn。

(14)

2.4 非穩(wěn)態(tài)排氣噪聲信號(hào)特征構(gòu)造

圖2所示為樣本7的原始信號(hào)和處理后樣本的頻譜圖，從圖7中可以看出，非穩(wěn)態(tài)排氣噪聲能量集中分布在500 Hz以下的低頻區(qū)域，然而基于Zwicker計(jì)算的尖銳度強(qiáng)調(diào)中高頻區(qū)域的能量，即中高頻成分在整體信號(hào)中的占比，也驗(yàn)證了表2中尖銳度與滿(mǎn)意度的弱相關(guān)性,粗糙度和波動(dòng)度反映人耳對(duì)聲音調(diào)制頻率和調(diào)制幅度的感受，峭度則反映信號(hào)的平坦程度，然而它們與滿(mǎn)意度的相關(guān)性較差，對(duì)非穩(wěn)態(tài)信號(hào)特征表現(xiàn)不足?；谏鲜龇治觯隒EEMD分解及PCA處理樣本熵的方法對(duì)非穩(wěn)態(tài)排氣噪聲信號(hào)進(jìn)行特征構(gòu)造，主要步驟如下。

步驟1濾波。首先對(duì)采集的20組非穩(wěn)態(tài)排氣噪聲樣本進(jìn)行巴特沃斯高通濾波，濾除20 Hz以下的次聲波。

圖2 非穩(wěn)態(tài)樣本7的信號(hào)Figure 2 Unsteady signal of sample 7

步驟2CEEMD分解。將高通濾波后的信號(hào)進(jìn)行CEEMD分解，并將n組白噪聲以正負(fù)成對(duì)的形式添加到原始信號(hào)當(dāng)中，由于目前對(duì)添加白噪聲幅值系數(shù)的選取沒(méi)有固定規(guī)則，故參考以往學(xué)者的研究經(jīng)驗(yàn)[11]，幅值系數(shù)選擇為0.4倍原噪聲信號(hào)的標(biāo)準(zhǔn)偏差，由此計(jì)算得到16個(gè)分量，分別是15個(gè)IMF分量和1個(gè)殘余分量。分解結(jié)果如圖3所示，可以看到原始信號(hào)經(jīng)過(guò)CEEMD分解為15個(gè)IMF分量和一個(gè)殘余分量res。通過(guò)CEEMD分解，對(duì)原始信號(hào)進(jìn)行多尺度分解，每進(jìn)行一次分解，去除高頻信號(hào)，使得15個(gè)IMF分量的頻率由高頻向低頻排列。

步驟3特征構(gòu)造。通過(guò)圖3可以看出第16階分量屬于分解趨勢(shì)項(xiàng)，對(duì)信號(hào)特征的提取沒(méi)有明顯意義，所以對(duì)前15階分量進(jìn)行樣本熵計(jì)算。根據(jù)參考文獻(xiàn)[12]中計(jì)算樣本熵時(shí)參數(shù)的選擇結(jié)果顯示，當(dāng)m取2，r取原始數(shù)據(jù)標(biāo)準(zhǔn)偏差的0.25倍時(shí)，樣本熵的計(jì)算結(jié)果更具有合理性。最終得到非穩(wěn)態(tài)信號(hào)的各樣本各階分量的樣本熵值如表3所示。

步驟4降維處理。首先對(duì)計(jì)算得到的樣本熵矩陣進(jìn)行標(biāo)準(zhǔn)化處理，得到標(biāo)準(zhǔn)化矩陣Z，再計(jì)算矩陣Z的相關(guān)系數(shù)矩陣R。在MATLAB中求解相關(guān)系數(shù)矩陣R的特征值和特征向量，提取累計(jì)貢獻(xiàn)率超過(guò)84%的主成分，各主成分特征值及其貢獻(xiàn)率如表4所示。顯然前2個(gè)主成分累計(jì)貢獻(xiàn)率已經(jīng)超過(guò)84%，但是為了更好地表現(xiàn)特征，故選取3個(gè)主成分進(jìn)行分析，按式(14)計(jì)算各個(gè)樣本的SQP-CSP值如表5所示。從第一主成分的特征值和特征向量分析發(fā)現(xiàn)IMF2、IMF3分量對(duì)非穩(wěn)態(tài)信號(hào)的特征表現(xiàn)更加明顯。

表3 各樣本各階分量的樣本熵值Table 3 Sample entropy of each order component for each sample

圖3 樣本7信號(hào)CEEMD分解結(jié)果Figure 3 CEEMD decomposition result of sample 7

表4 各主成分特征值及其貢獻(xiàn)率Table 4 Principal component eigenvalues and contribution rates

表5 非穩(wěn)態(tài)樣本信號(hào)的SQP-CSP值STable 5 SQP-CSP value of the unsteady sample signal

3 排氣噪聲聲品質(zhì)預(yù)測(cè)對(duì)比分析

3.1 基于客觀參數(shù)的聲品質(zhì)預(yù)測(cè)

圖4 GA-BP神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)Figure 4 Neural network topology structure of GA-BP

模型將tansig和purelin函數(shù)分別作為隱含層和輸出層的傳遞函數(shù)，選取梯度下降算法為網(wǎng)絡(luò)學(xué)習(xí)算法，學(xué)習(xí)效率為0.1，動(dòng)量系數(shù)設(shè)為0.9，訓(xùn)練目標(biāo)均方誤差設(shè)為0.001。遺傳算法參數(shù)如表6所示。

模型的訓(xùn)練過(guò)程如圖5所示，由訓(xùn)練結(jié)果可以看到，隨著遺傳算法迭代次數(shù)的增加，種群均值對(duì)應(yīng)的目標(biāo)值不斷減小，即適應(yīng)度不斷增大，大約在第50次迭代時(shí)趨于穩(wěn)定。遺傳算法優(yōu)化后，網(wǎng)絡(luò)訓(xùn)練的目標(biāo)函數(shù)值隨著訓(xùn)練次數(shù)的不斷增加而越來(lái)越小，在100次迭代后模型目標(biāo)值趨于穩(wěn)定，誤差達(dá)到設(shè)定目標(biāo)。訓(xùn)練結(jié)果擬合校驗(yàn)值R2為0.923，樣本期望值和訓(xùn)練值誤差較小。

表6 遺傳算法主要參數(shù)Table 6 Genetic algorithm main parameter

圖5 最佳適應(yīng)度曲線和GA-BP神經(jīng)網(wǎng)絡(luò)收斂曲線Figure 5 Best-fitness curve and GA-BP neural network convergence curve

訓(xùn)練完成后，將剩余的兩個(gè)非穩(wěn)態(tài)排氣噪聲信號(hào)的心理聲學(xué)客觀參量值作為訓(xùn)練GA-BP模型的輸入，兩種樣本預(yù)測(cè)的結(jié)果對(duì)比如圖6所示，驗(yàn)證誤差分別為12.19%、9.00%，平均驗(yàn)證誤差為10.595%。

圖6 兩種樣本預(yù)測(cè)結(jié)果Figure 6 Two prediction results of the model

3.2 基于SQP-CSP的聲品質(zhì)預(yù)測(cè)

將計(jì)算得到的SQP-CSP值S1、S2、S3歸一化后作為GA-BP模型的輸入，根據(jù)隱含層節(jié)點(diǎn)數(shù)的選取規(guī)則，建立一個(gè)3-5-1的聲品質(zhì)模型，遺傳算法的主要參數(shù)及目標(biāo)函數(shù)與前文訓(xùn)練心理聲學(xué)客觀參量模型保持一致。將前18組樣本的SQP-CSP值作為模型的輸入，主觀評(píng)價(jià)滿(mǎn)意度值作為模型輸出。訓(xùn)練后預(yù)測(cè)值與目標(biāo)值的均方根誤差為0.759，預(yù)測(cè)結(jié)果的相關(guān)系數(shù)R2為0.976，效果均優(yōu)于將心理聲學(xué)客觀參量作為輸入的預(yù)測(cè)。

將樣本19和20的SQP-CSP值S代入已訓(xùn)練好的模型中，得到預(yù)測(cè)結(jié)果。表7是兩種模型的預(yù)測(cè)結(jié)果對(duì)比。由此可以發(fā)現(xiàn)，IMF分量和樣本熵構(gòu)造出的聲信號(hào)特征再經(jīng)PCA降維后，較心理聲學(xué)客觀參量更適合用作非穩(wěn)態(tài)排氣噪聲聲品質(zhì)預(yù)測(cè)。

表7 兩種模型預(yù)測(cè)結(jié)果對(duì)比Table 7 Comparison of prediction results of two models

4 結(jié)論

(1)基于Zwicker計(jì)算了非穩(wěn)態(tài)排氣噪聲的心理聲學(xué)客觀參量，進(jìn)行了定性分析，結(jié)果表明主觀滿(mǎn)意度與響度與具有一定的線性關(guān)系，與尖銳度的線性相關(guān)性較弱，粗糙度、波動(dòng)度、峭度與主觀評(píng)價(jià)滿(mǎn)意度的相關(guān)性較差。同時(shí)，建立了GA-BP聲品質(zhì)預(yù)測(cè)模型，對(duì)非穩(wěn)態(tài)排氣噪聲進(jìn)行預(yù)測(cè)。

(2)利用CEEMD分解，獲得相應(yīng)的IMF分量，計(jì)算各階分量對(duì)應(yīng)的樣本熵值并得到新的特征向量，并通過(guò)PCA對(duì)數(shù)據(jù)進(jìn)行降維處理，得到新參量SQP-CSP值S與心理聲學(xué)客觀參量進(jìn)行對(duì)比分析，驗(yàn)證所構(gòu)造特征的可行性，預(yù)測(cè)結(jié)果表明新參量能夠有效提取信號(hào)特征，可以用于非穩(wěn)態(tài)信號(hào)的研究。