應啟帆，謝代梁，徐志鵬，徐 雅，劉鐵軍，黃震威

(中國計量大學浙江省流量計量技術研究重點實驗室，浙江 杭州 310018)

0 引 言

水體懸移質(zhì)粒徑分布測量在河流、湖泊、海洋和生態(tài)系統(tǒng)等領域中有著重要價值，其準確測量有助于認識懸移質(zhì)沉降、擴散規(guī)律[1]，對水土流失規(guī)律研究[2]和水污染治理[3]等領域的發(fā)展有重大意義。

粒徑測量方法按照工作原理[4]，可分為篩分法、顯微鏡法、沉降法、光散射法和超聲法等。其中，超聲法測量具有非浸入式、穿透性強、操作簡單等優(yōu)點，受到了廣泛關注。Shukla等[5]使用超聲波譜在線監(jiān)測乙酰氨基酚結晶過程中的顆粒生長；章維等使用超聲波衰減譜和相速度譜及Epstein-Carhart-Allegra-Hawley (ECAH)模型反演計算懸濁液顆粒粒徑分布[6]；Ding等[7]使用超聲頻譜測量了水分散體系中重油的粒徑；李燁明等采用人工蜂群算法優(yōu)化超聲粒徑分布計算結果[8]。

對粒徑分布的測量，傳統(tǒng)學者多聚焦在選用更精準的模型和更加優(yōu)化的算法來提升粒徑分布結果的準確度。蘇明旭等[9]通過構造與貝塞爾(Bessel)函數(shù)有關的變量改造系數(shù)矩陣，降低了矩陣求解的條件數(shù)，拓展了模型適用的頻率和顆粒粒徑。Silva等[10]使用 6～14 MHz頻率光譜法和擴展多重散射模型測量體積分數(shù)為10%～50%的葵花籽油包水型乳液的液滴尺寸分布。由于機器學習的發(fā)展和測量數(shù)據(jù)的提升，有學者嘗試用數(shù)據(jù)來尋找顆粒粒徑分布和輸入特征之間的關系，通過訓練機器學習模型來進行預測。Thompson等[11]利用基因表達編程和人工神經(jīng)網(wǎng)絡開發(fā)了一種新模型，用于預測活躍施工現(xiàn)場和露天采礦作業(yè)中裸露的土壤表面產(chǎn)生的雨水徑流中的總懸浮泥沙粒徑分布；Sch?fer等[12]使用卷積神經(jīng)網(wǎng)絡和人工圖像進行訓練，通過分割重疊顆粒來預測多相流中液滴的粒徑分布；Manee等[13]提出了一種深度學習傳感器，可對晶體生長過程中的顆粒進行在線監(jiān)測。

在本文中，我們提出了一種結合多輸出回歸算法和超聲衰減實驗的粒徑分布測量方法。首先根據(jù)超聲衰減實驗和 ECAH模型等先驗信息獲取實驗信號并提取特征，然后利用梯度提升決策樹算法預測單種粒徑，并組合全部粒徑預測模型構建多輸出回歸模型預測樣品的粒徑分布。最后將近似單峰分布、均勻分布和不規(guī)則分布的三種樣品粒徑分布預測結果與篩分法確定的粒徑分布進行對比，驗證多輸出回歸模型的準確性，可為粒徑分布測量提供一種新的參考方式。

1 基本理論

超聲法測量顆粒粒徑主要根據(jù)理論模型和優(yōu)化算法來獲得結果。理論模型介紹了超聲波在不同粒徑顆粒液體中的衰減情況，也成為了在機器學習算法中的特征基礎。接下來主要介紹 ECAH模型和超聲衰減實驗，了解粒徑測量的過程。

1.1 ECAH模型

ECAH模型說明了超聲波在水體中受到的粘滯力、熱損失和聲散射等因素的影響，可以描述超聲波具體的衰減過程。它忽略了顆粒之間的相互作用，適用于5%體積濃度以下、粒徑在1 000 um以下的懸移質(zhì)溶液環(huán)境，具體說明如式(1)[14-15]所示：

式中：α為超聲波在水中的衰減系數(shù)(單位：Np·m-1)；kc為連續(xù)介質(zhì)中的波數(shù)；φ為懸移質(zhì)溶液體積濃度，(單位：g·L-1)；R為懸移質(zhì)顆粒半徑(單位：μm)；An為懸移質(zhì)溶液中的超聲衰減系數(shù)，不同的下標n代表不同種衰減情況。

顆粒粒徑的求解過程是第一類 Fredholm 方程計算問題，可利用超聲衰減系數(shù)、濃度等參數(shù)之間的相關關系來求解，是典型的反問題，參數(shù)關系如式(2)[16]所示：

式中：q是顆粒粒度的頻度分數(shù)，是要求解的結果；K是消聲系數(shù)，是關于頻率和粒徑區(qū)間的函數(shù)。

通過對式(2)的離散化處理，可以轉化成矩陣形式來求解[16]：

式中： fi為不同測量頻率(單位：kHz)；Rj為不同粒徑區(qū)間(單位：Np·m-1)；α( fi, Rj)為在特定頻率和特定粒徑區(qū)間時的衰減系數(shù)(單位：Np·m-1)。

等式左邊是聲衰減過程，是超聲波在不同頻率和粒徑里的衰減過程，通過實驗獲得。等式右邊是理論部分，可以結合ECAH模型和各項物性參數(shù)來進行求解。式(1)～(3)描述了采用 ECAH模型求解顆粒粒徑分布的過程，是先驗信息。后續(xù)的特征部分就是根據(jù)這部分理論來選取與超聲衰減過程密切的相關參數(shù)作為數(shù)據(jù)集的輸入。

1.2 超聲衰減實驗

超聲衰減實驗采用兩個聚焦換能器進行聲電信號轉換，通過一發(fā)一收的工作方式，對電壓幅值進行測量，可探測到超聲波在懸移質(zhì)溶液中的衰減過程。聚焦換能器采用的材料為壓電材料PZT82；球面內(nèi)徑為 100 mm；投影直徑為 90 mm，換能器中心頻率為 750 kHz，單個聚焦換能器具體的結構如圖1所示。

圖1 超聲換能器結構Fig.1 Ultrasonic transducer structure

發(fā)射換能器接收信號發(fā)生器的電信號，將其轉化成聲信號發(fā)射，聲信號在穿過懸移質(zhì)溶液后到達接收換能器。在這個過程中超聲波會發(fā)生反射，從而多次到達接收換能器，采用示波器測量接收換能器提供的電壓信號，就可看出超聲波衰減過程，具體的測量裝置如圖2所示。

圖2 超聲測量裝置Fig.2 Ultrasonic measuring device

2 多輸出回歸算法

粒徑分布反演問題在機器學習中的表現(xiàn)是多輸入多輸出回歸問題，根據(jù)多個輸入找到對應的多個目標預測結果。本文采用組合梯度提升決策樹(Gradient Boosting Decision Tree, GBDT)算法來解決這一問題。首先利用 GBDT算法來找到單種粒徑的回歸預測結果，然后將所有GBDT算法組合起來，對整體粒徑分布進行預測。

2.1 梯度提升決策樹

梯度提升決策樹是一種集成學習算法，算法的基礎是分類回歸決策樹，利用損失函數(shù)的梯度方向決定后續(xù)回歸樹的優(yōu)化方向，通過貪心算法每一次找到下一步更優(yōu)的回歸樹，最后對所有回歸樹進行集成，輸出最后結果，就構成了梯度提升決策樹模型。算法訓練[17]的步驟如下所示：

算法每次對負梯度和輸入特征進行擬合，構建新的回歸樹來修正模型，其問題類型是多輸入單輸出的回歸預測模型，因此不適用于粒徑分布問題，只適用于單種粒徑分布的預測，需要將多個GBDT模型進行組合來解決問題。

2.2 GBDT回歸組合

對每種粒徑，單獨訓練一個GBDT模型，可以根據(jù)輸入特征和粒徑的數(shù)值進行回歸預測。為了預測樣品的粒徑分布，需要所有模型組合，每個輸出對應粒徑的結果，最后按照對整體粒徑的權重確定整體粒徑的頻度分布，多輸出回歸模型示意圖如圖3所示。

圖3 多輸出回歸模型示意圖Fig.3 Diagram of the multi-output regression model

3 特征選取

3.1 實驗信號

為了獲得實驗信號提取特征來制作數(shù)據(jù)集，搭建了如圖4所示的實驗測試系統(tǒng)。首先接收端聚焦換能器接收功率放大器放大的猝發(fā)波信號。超聲信號到達接收換能器時部分轉換成電壓幅值，部分反射回接收聚焦換能器。最后呈現(xiàn)在示波器上的是超聲波到達接收端聚焦換能器時的電壓幅值，最終將數(shù)據(jù)讀取并存儲在電腦中。實驗器材：發(fā)射和接收聚焦換能器、功率放大器、信號發(fā)生器、電磁攪拌器和示波器。

圖4 實驗系統(tǒng)Fig.4 Experimental system

實驗獲得的電壓信號經(jīng)去噪后的波形如圖5所示。同時考慮到換能器中心頻率的影響，提取其在頻域下的幅值，如圖6所示。

圖5 去噪后接收信號波形圖Fig.5 Waveforms of the received signals after denoising

圖6 去噪后接收信號頻譜圖Fig.6 Spectrums of the received signals after denoising

圖5中三個回波表示超聲波三次到達接收聚焦換能器時的電壓幅值，聲衰減系數(shù)就可通過這三個峰峰值計算而出，如式(7)[8]所示：

式中：αw為清水中的衰減系數(shù)(單位：Np·m-1)；V0是發(fā)射換能器電壓幅值(單位：V)，V1為接收換能器電壓幅值(單位：V)；L為兩個聚焦換能器之間的距離(單位：m)。

3.2 特征選取

多輸出回歸模型的特征來源于先驗信息，即顆粒的 ECAH模型和超聲衰減實驗的信號特征。接下來以各個參數(shù)與粒徑之間的關系來說明選擇它作為特征的理由。

根據(jù)ECAH模型反演來看，式(1)說明了粒徑與超聲衰減系數(shù)和體積濃度之間的關系，由圖5中的三次回波可計算出兩個衰減系數(shù)，因此將兩個聲衰減系數(shù)和體積濃度作為特征。式(3)說明了頻率可以影響聲衰減系數(shù)的計算，同時考慮到超聲換能器中心頻率對信號的影響，在選擇特征時加入頻率這一參數(shù)。

圖5說明了聲衰減變化的電壓衰減過程。圖6說明了圖5中三次回波在各頻率下的幅值，結果表明在中心頻率附近的傅里葉變換后幅值最大。三次回波的電壓峰峰值和三次回波的傅里葉變換最大幅值，是超聲衰減的關鍵因素，因此也選擇其作為特征。

綜上所述，選擇了懸移質(zhì)溶液體積濃度、兩次聲衰減系數(shù)、頻率、三個電壓峰峰值幅值和三次回波傅里葉變換后的最大幅值共 10個屬性作為數(shù)據(jù)集的特征，具體說明如表1所示。

表1 數(shù)據(jù)集輸入特征Table 1 Input features of data

為了避免特征數(shù)據(jù)量綱不同導致在回歸樹預測時帶來的偏差，需要對數(shù)據(jù)進行預處理。對特征為種類的頻率進行 onehot編碼，將 660、750、830 kHz三種頻率轉化為數(shù)值1、2、3；對其他連續(xù)數(shù)據(jù)進行標準化處理，將數(shù)據(jù)映射到[0,1]，轉換公式如式(8)[18]所示：

式中：x*是數(shù)據(jù)組中數(shù)據(jù)的轉換結果，x是數(shù)據(jù)組中數(shù)據(jù)初始值，xmin是數(shù)據(jù)組中數(shù)據(jù)最小值，xmax是數(shù)據(jù)組中數(shù)據(jù)最大值。

4 結果分析

4.1 數(shù)據(jù)集組成

本文的樣品通過篩分法進行配置，即先通過孔篩篩選出在各個粒徑區(qū)間范圍內(nèi)的顆粒粒徑，然后由精密天平稱重組合成近似各種分布類型的懸移質(zhì)樣品。以對照組的三種懸移質(zhì)樣本為例進行說明，三個樣本分別服從近似單峰分布、均勻分布和不規(guī)則分布，粒徑區(qū)間有 17種。單種粒徑區(qū)間質(zhì)量占樣本總質(zhì)量的百分比為該種粒徑區(qū)間的頻度，將 17種粒徑區(qū)間頻度進行組合，就構成了樣本的粒徑分布，具體的粒徑組成如表2所示。

表2 懸移質(zhì)顆粒樣品粒徑分布組成Table 2 Composition of the particle size distribution of suspended particles

本次實驗共配置了36組實驗樣品，其中33組用于制作訓練數(shù)據(jù)集，其余3組用于制作驗證數(shù)據(jù)集。對所有樣品進行三種頻率下的超聲衰減實驗，然后針對每一次實驗分別提取 10項特征作為數(shù)據(jù)集的輸入，并且將 17種粒徑分布的值作為數(shù)據(jù)集的輸出。

本文將訓練數(shù)據(jù)集用于訓練多輸出回歸模型，同時將驗證數(shù)據(jù)集輸入模型，即可得到粒徑分布預測結果。因為考慮到有三種頻率下的實驗，所以為了得到更加精確的結果，對三次預測結果進行均值化處理，作為真正的粒徑分布結果。

4.2 預測結果和誤差分析

中位徑D50代表一個樣品的累計粒度分布百分數(shù)達到50%時所對應的粒徑，中位徑誤差ε是常用的粒徑誤差計算公式，如式(9)[19]所示：

式中：D50,act為實際的中位徑(單位：μm)；D50,pre是模型預測確定的中位徑(單位：μm)。

首先利用訓練數(shù)據(jù)集訓練多輸出回歸模型，然后利用驗證數(shù)據(jù)集和模型預測三種樣品的粒徑分布結果。預測得到的單峰分布、均勻分布和不規(guī)則分布樣品篩分結果和多輸出回歸模型預測結果對比圖分別如圖7～9所示。

圖 7～9的橫坐標是粒徑區(qū)間的組合，柱狀圖分別為篩分法確定的該種粒徑頻率分布和模型預測的粒徑頻率分布，縱坐標為粒徑的頻率分布，即單種粒徑區(qū)間質(zhì)量占樣品質(zhì)量的百分比。所有粒徑的結果進行組合就構成了樣本粒徑分布的結果。接下來統(tǒng)計了每種樣品的單個粒徑的最大相對誤差、決定系數(shù)和中位徑誤差，如表3所示。

表3 不同樣本中位徑和測量誤差Table 3 Median diameters and measurement errors of different samples

圖7 單峰分布樣品兩種方法測量結果對比Fig.7 Comparison of measurement results of two methods for unimodal samples

三種樣品的單種顆粒相對誤差范圍都在±10%以內(nèi)，只有單峰分布的個別樣本在粒徑權重較小時出現(xiàn)了較大的偏差，為 12.74%。這表明，單種粒徑占總體粒徑分布的權重質(zhì)量較小時會導致相對較大的偏差。中位徑的誤差都較小，表明預測分布與實際分布較為一致，可以很好地反應懸移質(zhì)顆粒粒徑分布組成。

圖8 均勻分布樣品兩種方法測量結果對比Fig.8 Comparison of measurement results of two methods for uniformly distributed samples

圖9 不規(guī)則分布樣品兩種方法測量結果對比Fig.9 Comparison of measurement results of two methods for irregularly distributed samples

5 結 論

本文通過 ECAH模型和超聲衰減實驗等先驗信息，了解傳統(tǒng)粒徑分布測量的過程，并說明了顆粒粒徑與體積濃度、聲衰減系數(shù)之間存在的耦合關系，使用機器學習算法對粒徑分布進行預測來為粒徑分布測量提供一種新的方法。

粒徑分布預測作為一種多輸入多輸出回歸問題，采用組合梯度提升決策樹來進行預測。首先使用 GBDT對單種粒徑進行預測，然后組合所有粒徑預測過程，就構成了對粒徑分布結果的預測。然后選擇了與粒徑關系密切的體積濃度、二次聲衰減系數(shù)和超聲衰減實驗中最明顯的三次回波電壓峰峰值和傅里葉變換最大幅值共 10項參數(shù)作為數(shù)據(jù)集的特征。最后，采用近似單峰分布、均勻分布和不規(guī)則分布的三種樣品作為預測結果，三種樣品中單種粒徑相對誤差范圍在±10%以內(nèi)，中位徑誤差分別為0.07%、-0.1%和-2.2%，表明預測樣品結果與實際樣品結果較為一致，可以很好地測量樣品的粒徑分布。

本文利用的特征根據(jù)先驗信息進行選取，也可以進一步通過對超聲信號的分解來選取與粒徑關系密切的特征來進行訓練。同時，使用更大規(guī)模的數(shù)據(jù)集來進行訓練，可得到某個區(qū)域更準確的粒徑分布結果。本文的嘗試，可以為粒徑分布測量問題提供一種新的思路。