烏魯瓦提大壩監(jiān)測數(shù)據(jù)去噪統(tǒng)計模型研究

毛建剛，文 俊，朱明遠(yuǎn)

(1.新疆水利水電科學(xué)研究院，烏魯木齊 830049；2.新疆農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院，烏魯木齊 830052)

0 前 言

大壩安全監(jiān)測數(shù)據(jù)是判斷大壩安全、掌握大壩運行規(guī)律的依據(jù)。通過監(jiān)測數(shù)據(jù)構(gòu)建大壩安全預(yù)測模型，對實際工程運行監(jiān)管有著重要的指導(dǎo)意義[1]。近年來，統(tǒng)計回歸、遺傳算法、神經(jīng)網(wǎng)絡(luò)、時間序列[2]、統(tǒng)計回歸分析、小波分析[3]和卡爾曼濾波[4]等技術(shù)方法廣泛應(yīng)用于大壩變形監(jiān)測分析的工程實踐中。但是單一模型很難解決實際工程問題，多種方法結(jié)合可提高大壩監(jiān)測數(shù)據(jù)分析的精度和可靠性。

高斯白噪聲是數(shù)據(jù)序列中功率譜密度在整個頻域內(nèi)均勻分布的噪聲，監(jiān)測數(shù)據(jù)采集過程中受多因素的影響，導(dǎo)致采集的數(shù)據(jù)序列中摻雜多種噪聲，常用的數(shù)字信號噪聲處理技術(shù)有小波分析、卡爾曼濾波、時間序列、經(jīng)驗?zāi)Ｊ椒纸?LMD)等。傳統(tǒng)回歸模型僅考慮過程噪聲，缺乏高斯白噪聲的考慮，這降低了回歸模型的精度，直接影響分析結(jié)果的可靠性，小波分析被譽為“數(shù)學(xué)顯微鏡”，它通過對時間序列的多時間尺度分析，能夠有效識別、提取數(shù)據(jù)序列的主要頻率成分以及局部信息，但這一過程也容易導(dǎo)致數(shù)據(jù)失真[5]。卡爾曼濾波是去除噪聲還原真實數(shù)據(jù)的一種數(shù)據(jù)處理技術(shù)，是在最小均方差原理下對數(shù)據(jù)序列進(jìn)行最優(yōu)估計的一種數(shù)據(jù)處理技術(shù)[6]。因此用小波去噪和卡爾曼濾波進(jìn)行結(jié)合，將去噪后的數(shù)據(jù)導(dǎo)出后做統(tǒng)計回歸分析，提高預(yù)測的精度和可靠性。

隨著壩工理論的發(fā)展，單一方法和途徑不再適合于復(fù)雜的監(jiān)測數(shù)據(jù)分析中，而多種理論和方法的有機(jī)結(jié)合將更有利于提高數(shù)據(jù)分析的精度和掌握大壩變形規(guī)律[7]。本文針對烏魯瓦提大壩自動化升級改造后監(jiān)測數(shù)據(jù)受噪聲污染情況復(fù)雜、采用傳統(tǒng)統(tǒng)計模型分析的結(jié)果不能滿足實際需要等問題，利用原型觀測資料，以小波分析和卡爾曼濾波濾除原始含噪信號中的高斯白噪聲和粗差，結(jié)合前人經(jīng)驗選取影響因子建立統(tǒng)計回歸模型，以此提高監(jiān)測數(shù)據(jù)分析處理的精確度和可靠性。

1 工程概況

1.1 基本概況

烏魯瓦提水利樞紐工程位于新疆和田地區(qū)喀拉喀什河中游段，是喀拉喀什河的控制性工程，壩址距和田市71 km，水庫總庫容3.47億m3。其中，主壩最大壩高133.0 m，壩頂高程1 965.80 m，壩頂長365.0 m，壩頂寬8.9 m，具有灌溉、發(fā)電、生態(tài)、防洪等綜合效益。工程建成投運后，徹底解決了和田地區(qū)的近13.33 hm2(200萬畝)的農(nóng)業(yè)灌溉用水問題，提高了下游河道的防洪能力，減輕了洪水災(zāi)害,工程經(jīng)濟(jì)效益、社會效益、環(huán)境效益十分顯著。

1.2 數(shù)據(jù)選取

烏魯瓦提大壩安全監(jiān)測數(shù)據(jù)采集系統(tǒng)是基于Windows NT網(wǎng)絡(luò)環(huán)境下，采用南瑞DAMS-IV型智能分布式數(shù)據(jù)采集系統(tǒng)進(jìn)行數(shù)據(jù)采集，直接使用去除采粗差后的采集數(shù)據(jù)進(jìn)行計算分析，結(jié)果精確度較低。因此本文選取河床最大斷面0+190.00 m斷面(如圖1所示)H2-2和H2-3兩測點共235組水平位移監(jiān)測數(shù)據(jù)(后10期數(shù)據(jù)驗證)為樣本，甄選出精度較高的去噪模型并建立統(tǒng)計預(yù)測模型，以還原真實數(shù)據(jù)，降低誤差，提高數(shù)據(jù)分析的精度。

2 監(jiān)測數(shù)據(jù)的去噪處理

2.1 小波分析

小波分析是一種良好的時頻分析工具，它通過準(zhǔn)確的分析、診斷、編碼壓縮和量化、快速傳遞或存儲、精確重構(gòu)原始信號，最終實現(xiàn)恢復(fù)真實數(shù)據(jù)信號的目的。對于任意函數(shù)或者信號f(x)，其小波變換為[8]：

式中：Wf(a,b)為小波系數(shù);a反映函數(shù)的尺度;b是沿時間軸或位置軸平移的位置。

為了應(yīng)用方便，通常將連續(xù)小波及其變換離散化，此時可表示為：

圖1 0+190.00 m斷面水平位移測點布置 單位：m

小波變換表示為離散小波變換：

小波去噪通常是將大壩的原始監(jiān)測數(shù)據(jù)看作是一串信號，由有用信號和噪聲共同組成[9]，通常表達(dá)為：

f(t)=s(t)+n(t)

(4)

式中:f(t)為原始信號，s(t)是有用信號，n(t)是噪聲且符合N(0，σ2)。

實際應(yīng)用中，小波變換具有帶通濾波的功能，可以將信號劃分為不同的頻帶，這里若設(shè)原始信號的頻率為f，在尺度參數(shù) j = 1，2，…，J下，應(yīng)用小波分解，其對應(yīng)的頻帶數(shù)為2J，相應(yīng)的頻率范圍為2-J(i-1)f～2-Jif,其中i = 1，2，…，2J，表示分解信號的頻帶序列，分解后可以看出從低頻到高頻的信號信息，各頻帶互不重疊。

根據(jù)前人研究成果和經(jīng)驗[10]，選用dbN小波時，最優(yōu)分解尺度N選取越大，被濾掉的噪聲越多，越有利于噪聲分離，但是在重構(gòu)時信號失真也越大，所以必須選擇一個合適的分解尺度。實測數(shù)據(jù)一般信噪比不能預(yù)知，文鴻雁[3]提出逐漸增大尺度，根據(jù)均方根誤差值的變化是否趨于穩(wěn)定來確定最大尺度，一般情況下，當(dāng)r接近于1時，噪聲基本去除。本文選取db4小波對原始數(shù)據(jù)進(jìn)行1～4層分解得到均方根誤差RMSE和r值,具體見表1。因此選取db4小波對原始數(shù)據(jù)進(jìn)行1層分解。

表1 不同尺度下均方根誤差

具體分解步驟如下：

(1)監(jiān)測數(shù)據(jù)的小波分解。應(yīng)用Matlab小波工具箱wavede函數(shù)，對信號進(jìn)行分解，調(diào)用格式如下：

[c,l]=wavedec(x,N,'wname',)

(5)

式中：c表示各層分量，包括近似系數(shù)和細(xì)節(jié)系數(shù)，l表示各層分量長度，x表示原始信號，N分解的層數(shù)，wname小波基名稱。

(2)小波閾值量化。對第1層到第N層的每一層高頻系數(shù)，選擇一個閾值進(jìn)行量化處理。

(3)重構(gòu)。重構(gòu)則是分解的逆過程，對低頻系數(shù)、高頻系數(shù)分別重構(gòu)處理。

H2-2(a)、H2-3(b)測點原始信號、分解后低頻信號和高頻信號對比圖分解結(jié)果如圖2。

圖2 兩測點原始信號、分解后低頻信號和高頻信號對比

由圖2分解的信號可以得出，兩測點信號服從高斯分布，功率譜密度均勻分布，是含高斯白噪聲的監(jiān)測序列。小波分解與重構(gòu)去噪法能夠?qū)υ肼曔M(jìn)行去除，提取出信號的高頻和低頻部分，但是由于它在分解的過程中只對低頻信號再分解，對高頻信號不再實施分解，使得它的頻率分辨率隨頻率升高而降低，丟失部分有用信息。該方法雖然能去除噪聲，但對于有用信號和噪聲的頻帶相互重疊的情況，效果就不甚理想。

2.2 卡爾曼濾波

小波去噪適用于處理已知噪聲的頻率范圍且信號和噪聲的頻帶相互分離時的數(shù)據(jù)序列，但實際監(jiān)測數(shù)據(jù)中廣泛存在高斯白噪聲，其效果并不理想?？柭鼮V波是去除噪聲還原真實數(shù)據(jù)的一種數(shù)據(jù)處理技術(shù)，能夠從一系列存在噪聲的數(shù)據(jù)中進(jìn)行數(shù)值最優(yōu)估計。因此在小波去噪后導(dǎo)出數(shù)據(jù)并銜接卡爾曼濾波原理，借助Matlab平臺編輯濾波程序，進(jìn)行小波-卡爾曼濾處理，濾除數(shù)據(jù)序列中存在的高斯白噪聲，以期獲得較高的去噪效果[11]。多數(shù)情況下，主要運用的是離散卡爾曼濾波，其數(shù)學(xué)模型[12]：

Xk=Fk/k-1Xk-1+Gk-1Wk-1

(6)

Lk=HkXk+Vk

(7)

上式為卡爾曼濾波的函數(shù)模型，由最小二乘原理，可推導(dǎo)出卡爾曼濾波的遞推公式：

狀態(tài)向量的一步預(yù)測：

狀態(tài)向量一步預(yù)測值方差矩陣：

狀態(tài)向量估計值：

狀態(tài)向量估計方差矩陣：

Pk=(1-JkHk)Pk/k-1

(11)

其中，Jk為濾波增益矩陣，具體形式為：

式中:J是濾波增益矩陣;Hk為觀測矩陣;Xk為系統(tǒng)在tk時刻的狀態(tài)n×1維向量;Fk/k-1為狀態(tài)轉(zhuǎn)移矩陣;Lk為觀測值向量。

經(jīng)過小波去噪后可以看出2個測點含有高斯白噪聲，為了驗證單一模型和組合模型對含高斯白噪聲信號處理的精度，擬采用2種方案進(jìn)行濾波處理，方案如下：

(1) 方案1用原始數(shù)據(jù)直接卡爾曼濾波處理;

(2) 方案2用小波去噪后的數(shù)據(jù)結(jié)合卡爾曼濾波處理。

濾波效果如圖3。

圖3 兩測點2種濾波方案對比

圖3濾波效果表示，經(jīng)過小波預(yù)處理后再結(jié)合卡爾曼濾波對監(jiān)測信號高斯白噪聲的處理效果更好，信號更加光滑且不失真。

2.3 精度評定

為了更加科學(xué)地驗證2種去噪方案的精度，對兩測點水平位移監(jiān)測數(shù)據(jù)進(jìn)行濾波后模型精度評定。模型精度評定公式如下[13]:

(1) 均方誤差(MSE)

(2) 平方和誤差(SSE)

(3) 平均相對誤差(MRE)

表2兩個測點的濾波方案中，方案2的均方誤差(MSE)、平方和誤差(SSE)和平均相對誤差(MRE)均較方案1小，表明當(dāng)監(jiān)測信號含有高斯白噪聲時小波分析與卡爾曼濾波組合模型明顯優(yōu)于單一的去噪方法，去噪精度更高。

表2 精度評定

3 統(tǒng)計回歸模型的建立

3.1 回歸分析

統(tǒng)計回歸分析在相關(guān)變量中將一個變量視為因變量，其他一個或多個變量視為自變量，建立線性或非線性數(shù)學(xué)模型的一種統(tǒng)計分析方法。這里設(shè)因變量為Y，影響因變量的n個自變量分別為x1,x2,…xn，假設(shè)n個自變量里面每一個自變量對因變量Y的影響都是線性的變化，得數(shù)學(xué)表達(dá)式為[14-15]：

Y=β0+β1x1+β2x2+…+βnxn+ε

(16)

式中:Y是因變量;β0,β1,…,βn為回歸參數(shù);ε為隨機(jī)變量或誤差;εi～N(0,σ2)，i=1，2，...，n。

3.2 因子選取

大壩位移影響因子中最優(yōu)影響因子的選取是構(gòu)建大壩安全統(tǒng)計回歸模型的關(guān)鍵。模型越復(fù)雜、包含參數(shù)越多，擬合和預(yù)測數(shù)據(jù)的功能越強[16]。烏魯瓦提大壩位移變形影響因素包括壩型、壩基、填筑材料、揚壓力、時效、水位變化及溫度變化等。面板堆石壩的填筑堆石料內(nèi)部變形受到外界溫度變化的影響，變形和力學(xué)特性會發(fā)生改變，溫度循環(huán)造成堆石料一定程度上的損傷，導(dǎo)致其顆粒強度降低，位移量增大[17]。諸多因素中有些難以量化，本文根據(jù)已有的面板堆石壩安全監(jiān)測研究成果，選取水位、溫度、時效等作為影響大壩位移變形的因子，表達(dá)為:

y(t)=yH(t)+yT(t)+yθ(t)

(17)

式中:y(t)為監(jiān)測效應(yīng)量y在時刻t的統(tǒng)計估計值; yH(t)、yT(t)、yθ(t)分別為y(t)的水壓分量、溫度分量和時效分量。

3.3 模型的建立

本文選取上下游水位差H的1、2、3次方為水壓因子(H1、H2、H3)。在溫度因子選擇上，引起堆石壩內(nèi)部溫度場變化的主要原因是氣溫和水溫，而水溫的變化也是由壩址的氣溫變化引起的，只是在時間上有滯后，因此溫度分量選取可用變形監(jiān)測日的平均氣溫來作為周期因子，i取當(dāng)天或前3、7、15、30、45、60d的平均值[18-19]；時效因子選取觀測當(dāng)天至觀測起測日的累計天數(shù)除以100(θ=t/100)的θ和lnθ值，綜上模型可表達(dá)為：

[b1i(siniwt-siniwt0)+(cosiwt-cosiwt0)]

式中:α0，αi，bi，ci為回歸系數(shù)，w=2π/365;H為上下游水位差;Ti為觀測當(dāng)天溫度及前i天的平均溫度;t從某天起算的時間;c1θ是時效分量隨時間呈線性變化;c1lnθ是時效分量隨時間呈非線性變化;m為多項式階數(shù)。

根據(jù)上述分析，擬采用2種方案建立統(tǒng)計回歸模型，方案如下：

(1) 方案1為用原始數(shù)據(jù)直接建立的統(tǒng)計模型；

(2) 方案2為用小波-卡爾曼濾波后數(shù)據(jù)建立的統(tǒng)計模型。

借助SPSS軟件平臺，運用多元回歸分析得出兩種方案模型的匯總表，分別見表3和表4。表3和表4顯示，經(jīng)小波-卡爾曼去噪后建立的統(tǒng)計回歸模型的各參數(shù)值提高，標(biāo)準(zhǔn)估計的誤差降低，模型的回歸精度提升。

表3 方案1模型匯總

表4 方案2模型匯總

3.4 監(jiān)測數(shù)據(jù)的預(yù)測

為了進(jìn)一步驗證去噪預(yù)處理后模型的實用性和優(yōu)越性，對2種方案的數(shù)據(jù)進(jìn)行對比分析。用后10期實測值為樣本與預(yù)測值對比，繪制監(jiān)測點的預(yù)測效果如圖4所示。

計算絕對誤差平均(MAE)和均方根誤差(RMSE)，作為預(yù)測模型評價的標(biāo)準(zhǔn)。具體預(yù)測值和2種方案的MAE、RMSE值見表5。

圖4 兩測點2種方案的預(yù)測值對比

表5 2種方案的實測值和預(yù)測值對比

圖4顯示，相比方案1，方案2各時段的預(yù)測結(jié)果更接近實測值；表5顯示，方案1兩測點的絕對誤差平均(MAE)和均方根誤差(RMSE)分別為0.33、0.08和0.35、0.09，方案2兩測點的絕對誤差平均(MAE)和均方根誤差(RMSE)分別為0.08、0.03和0.09、0.03，兩測點的MAE和RMSE都有所降低，但方案2降低得較為明顯，上述分析結(jié)果表明方案2預(yù)測精度更高，用小波-卡爾曼濾波后數(shù)據(jù)建立的統(tǒng)計模型的預(yù)測結(jié)果更符合工程實際。

4 結(jié) 論

本文以工程實際為研究對象，針對監(jiān)測數(shù)據(jù)中含高斯白噪聲導(dǎo)致傳統(tǒng)統(tǒng)計模型預(yù)測精度低的問題，探討小波分析和卡爾曼濾波在監(jiān)測數(shù)據(jù)噪聲處理中的適用性，得出以下結(jié)論：

(1) 小波分析能夠剔除原始監(jiān)測信號中的粗差和噪聲，但易丟失部分有用信息；卡爾曼濾波能夠有效濾除監(jiān)測信號中的高斯白噪聲，在小波分解的基礎(chǔ)上引入卡爾曼濾波處理，其結(jié)果較單一濾波方法好。

(2) 經(jīng)小波-卡爾曼濾波方法數(shù)據(jù)預(yù)處理后建立的統(tǒng)計模型，回歸參數(shù)明顯增大，標(biāo)準(zhǔn)估計誤差明顯降低，提高了模型的回歸精度；預(yù)處理后數(shù)據(jù)建模的預(yù)測結(jié)果與實測結(jié)果較為接近，MAE和RMSE都有所降低，預(yù)測精度提高明顯，該模型對工程實際具有一定的實用性，有參考價值。