依據(jù)主成分和協(xié)整性的大壩變形奇異診斷*

楊光，李姝昱，孫錦

（1.華北水利水電大學(xué)水利學(xué)院 鄭州，450046）

（2.河海大學(xué)水文水資源與水利工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室 南京，210098）

（3.黃河水利科學(xué)研究院 鄭州，450003）

（4.華北水利水電大學(xué)測(cè)繪與地理信息學(xué)院 鄭州，450046）

引言

變形是大壩服役性態(tài)變化的綜合反映，是衡量結(jié)構(gòu)安全與否的重要標(biāo)志。科學(xué)地分析變形原位監(jiān)測(cè)信息，不僅是現(xiàn)行規(guī)范的要求，而且是監(jiān)控大壩安全的有效手段［1-2］。受多種因素的干擾，原位監(jiān)測(cè)資料中常包含奇異成分，表現(xiàn)為測(cè)值在某時(shí)刻或某時(shí)段的異常突跳，即孤立型或斑點(diǎn)型奇異成分。若未有效診斷出監(jiān)測(cè)資料中的奇異成分，極可能影響大壩安全監(jiān)控結(jié)論的客觀性，造成虛假報(bào)警或?qū)⑽ｋU(xiǎn)狀況遺漏。尤其對(duì)于運(yùn)行多年的老壩，監(jiān)測(cè)儀器的性能和穩(wěn)定性較差，此問(wèn)題尤為突出。因此，提出科學(xué)的大壩變形原位監(jiān)測(cè)數(shù)據(jù)奇異成分診斷方法，具有重要的理論意義和應(yīng)用價(jià)值。

大壩變形原位監(jiān)測(cè)數(shù)據(jù)奇異成分診斷包括辨識(shí)和估計(jì)2個(gè)環(huán)節(jié)。常規(guī)的辨識(shí)方法有［3-4］：①過(guò)程線法，通過(guò)繪制變形過(guò)程線，直觀辨識(shí)奇異成分；②假設(shè)檢驗(yàn)法，如拉依達(dá)檢驗(yàn)、狄克松檢驗(yàn)。傳統(tǒng)的估計(jì)方法包括［3-4］：①忽略法，若無(wú)法還原或難度較大，可直接忽略；②似然估計(jì)法，如臨近插值、三次Hermite插值、統(tǒng)計(jì)模型及自回歸模型?？傮w來(lái)看，常規(guī)診斷方法大多以變形時(shí)間序列為分析對(duì)象，無(wú)法有效地捕捉多測(cè)點(diǎn)間隱含的變形關(guān)系，且主觀性較強(qiáng)，易誤診。

相較于1維時(shí)間序列，多維數(shù)據(jù)序列具有信息量豐富、自由度高及穩(wěn)定性強(qiáng)的優(yōu)勢(shì)，包含了隱匿、新穎及有潛在價(jià)值的信息［5］。PCA［6-9］和CA［10-12］是兩種多維數(shù)據(jù)特征挖掘方法。PCA可將多維數(shù)據(jù)序列投影到不同方向上，得到主要規(guī)律和無(wú)法解釋成分，其中無(wú)法解釋成分包含了奇異數(shù)據(jù)、噪聲等。對(duì)于大壩上臨近的監(jiān)測(cè)點(diǎn)，變形監(jiān)測(cè)序列存在關(guān)聯(lián)性，若未出現(xiàn)超標(biāo)準(zhǔn)洪水、極端天氣等非常規(guī)狀況，某監(jiān)測(cè)點(diǎn)變形數(shù)據(jù)異常，其他臨近點(diǎn)的數(shù)據(jù)未見(jiàn)異常，則可判定該監(jiān)測(cè)點(diǎn)的變形數(shù)據(jù)奇異。CA刻畫(huà)了多維數(shù)據(jù)序列的長(zhǎng)期關(guān)系，若未發(fā)生非常規(guī)狀況，多測(cè)點(diǎn)變形監(jiān)測(cè)序列的協(xié)整關(guān)系一般不會(huì)明顯改變，據(jù)此可實(shí)現(xiàn)對(duì)奇異成分的估計(jì)。

綜上所述，筆者利用多維數(shù)據(jù)的優(yōu)勢(shì)，基于PCA，研究奇異成分辨識(shí)準(zhǔn)則；依據(jù)CA，探究奇異成分似然估計(jì)模型；在此基礎(chǔ)上，結(jié)合實(shí)際工程，檢驗(yàn)PCA-CA方法的有效性。

1 診斷方法

1.1 奇異成分辨識(shí)

依據(jù)2維數(shù)據(jù)格式，建立大壩變形原位監(jiān)測(cè)數(shù)據(jù)集合，即

其中：m為監(jiān)測(cè)點(diǎn)總數(shù)；n為監(jiān)測(cè)時(shí)長(zhǎng)；xi，xj分別為測(cè)點(diǎn)i和j的監(jiān)測(cè)序列；xki，xkj分別為測(cè)點(diǎn)i和j在k時(shí)刻的監(jiān)測(cè)值。

PCA辨識(shí)步驟如下。

1）分析環(huán)境量監(jiān)測(cè)資料，若未發(fā)生非常規(guī)狀況，則執(zhí)行后續(xù)步驟。

2）對(duì)矩陣X進(jìn)行標(biāo)準(zhǔn)化，得到矩陣,為標(biāo)準(zhǔn)化的xi，為 標(biāo) 準(zhǔn)化的xj。相關(guān)系數(shù)矩陣R為

其 中：rij為和的相關(guān)系數(shù)，計(jì)算公式為rij=(i，j=1，2，…，m)；和為標(biāo)準(zhǔn)化 的xki和xkj。

3）計(jì)算相關(guān)系數(shù)矩陣R的特征向量矩陣P和特征值矩陣λ，計(jì)算公式為

其中：0為m階零矩陣；I為m階單位矩陣。

式（3）中的P和λ分別表示為

其 中：λi為矩陣R的第i個(gè)特征值，滿(mǎn)足λ1≥λ2≥…≥λm；pi為λi對(duì)應(yīng)的特征向量。

4）建立n×m矩陣T，記為(t1，t2，…，tj，…，tm)

5）計(jì)算第j項(xiàng)tj(j=1，2，…，m)對(duì)原始序列X的解釋能力ej

由式（7）可知：①前j項(xiàng)對(duì)X的累積解釋能力為，所有m項(xiàng)的累積解釋能力為1；②e1≥e2≥…≥em，即t1→tm對(duì)X的解釋能力逐漸變?nèi)?；③e1越大，X的主要規(guī)律越明顯。因此，為保證PCA辨識(shí)精度，要求e1盡可能大。

6）建立SPE統(tǒng)計(jì)量。由于P為正交矩陣，滿(mǎn)足

則式（6）可變換為

k時(shí)刻SPE統(tǒng)計(jì)量記為SPEk，計(jì)算公式為

7）依據(jù)假設(shè)檢驗(yàn)，在顯著性水平α下，建立SPE控制限［13-14］，記為SPEα，即

假設(shè)檢驗(yàn)中α常取為0.05或0.01，即接受原假設(shè)時(shí)，正確的概率為95%或99%。

8）圖1為SPE統(tǒng)計(jì)與控制限的關(guān)系。若SPEk未超過(guò)SPEα，如圖1（a）所示，則不存在奇異成分；若SPEk在某時(shí)段超過(guò)SPEα，如圖1（b）所示，則為斑點(diǎn)型奇異；若SPEk在某時(shí)刻超過(guò)SPEα，如圖1（c）所示，則為孤立型奇異。若奇異，則執(zhí)行步驟9。

圖1 SPE統(tǒng)計(jì)量與控制限的關(guān)系Fig.1 Relationships between SPE and its control limit

9）測(cè)點(diǎn)j對(duì)k時(shí)刻SPEk超過(guò)SPEα的貢獻(xiàn)度用CSPEkj統(tǒng)計(jì)量表征，即

CSPEkj越大，則測(cè)點(diǎn)j對(duì)k時(shí)刻SPEk超限 的貢獻(xiàn)度越大，而CSPEkj最大處，即為存在奇異成分的測(cè)點(diǎn)。

1.2 奇異成分似然估計(jì)

假設(shè)m個(gè)變形監(jiān)測(cè)序列中有f個(gè)存在奇異成分，記為x′1…x′h…x′f，剩余p=m-f個(gè)不存在奇異成分，記為。CA似然估計(jì)步驟如下。

1）建立多測(cè)點(diǎn)變形原位監(jiān)測(cè)序列的協(xié)整模型，即

其中：cqg為的系數(shù)；aq為的最高次數(shù)，依據(jù)變形散點(diǎn)關(guān)系確定；ε為余量序列。

2）應(yīng)用逐步回歸算法，確定式（14）的系數(shù)cqg，得到余量序列ε。

3）利用ADF檢驗(yàn)［15-16］，測(cè)試ε的平穩(wěn)性。若ε平穩(wěn)，則協(xié)整，執(zhí)行步驟4；若非平穩(wěn)，則不協(xié)整，此時(shí)算法不適用。

4）利用式（14），估計(jì)x′h的 奇 異成分，修正矩陣X。

5）利用PCA準(zhǔn)則，辨識(shí)矩陣X中是否仍然存在奇異成分，若存在，重復(fù)執(zhí)行步驟1～4，直至無(wú)奇異成分。

2 工程實(shí)例分析

2.1 工程實(shí)例1

陳村拱壩最大壩高為76.3 m，設(shè)計(jì)汛限水位為117.0 m，汛后正常高水位為119.0 m，500年校核水位為123.8 m，5 000年校核水位為126.1 m，保壩水位為127.7 m。18正下（1號(hào)測(cè)點(diǎn)），26正下（2號(hào)測(cè)點(diǎn)）和29倒（3號(hào)測(cè)點(diǎn)）測(cè)點(diǎn)的位置如圖2所示，徑向變形過(guò)程線和置信區(qū)間法計(jì)算結(jié)果如圖3所示，環(huán)境量過(guò)程線如圖4所示。

圖2 29#，26#和18#壩段垂線監(jiān)測(cè)點(diǎn)Fig.2 Vertical monitoring points of the 29#,26# and 18#dam sections

依據(jù)水壓-溫變-時(shí)效（hydraulic-seasonal-time，簡(jiǎn)稱(chēng)HST）模式［17］，建立變形監(jiān)控模型，HST模型精度如表1所示。由表可知，3個(gè)模型的復(fù)相關(guān)系數(shù)R均大于0.95，且標(biāo)準(zhǔn)差S較小，擬合精度較高。應(yīng)用置信區(qū)間法［17］，評(píng)判變形變化是否正常，2S為正常與基本正常的臨界值，3S為基本正常與異常的臨界值。由圖3（b）～（d）可知：①1號(hào)測(cè)點(diǎn)未超出3S，在42天超出了2S，未報(bào)警；②2號(hào)測(cè)點(diǎn)未超出3S，在23天超出了2S，未報(bào)警；③3號(hào)測(cè)點(diǎn)在2天超出了3S，且在148天超出了2S，發(fā)生報(bào)警。

表1 HST模型精度Tab.1 Accuracy of the HST models

圖3 徑向變形過(guò)程線和置信區(qū)間法計(jì)算結(jié)果Fig.3 Time series of radial deformation and calculation results of confidence interval approach

由圖4可知，該時(shí)段未出現(xiàn)非常規(guī)水位和溫度狀況，因此受奇異成分影響，監(jiān)控模型虛假報(bào)警。利用本研究方法，診斷3個(gè)測(cè)點(diǎn)的奇異監(jiān)測(cè)數(shù)據(jù)。經(jīng)計(jì)算，e1為95.23%，且余量序列協(xié)整，診斷出的奇異監(jiān)測(cè)數(shù)據(jù)有斑點(diǎn)型和孤立型，如：1號(hào)測(cè)點(diǎn)的變形監(jiān)測(cè)值在2012-07-11為孤立型奇異；3號(hào)測(cè)點(diǎn)的測(cè)值在2011-05-25至2011-07-19期間為斑點(diǎn)型奇異。

圖4 環(huán)境量過(guò)程線Fig.4 Time series of environment factors

重新建立監(jiān)控模型，結(jié)果顯示：①1號(hào)測(cè)點(diǎn)在3天超出了2S，未超過(guò)3S；②2號(hào)測(cè)點(diǎn)在5天超出了2S，未超 過(guò)3S；③3號(hào) 測(cè)點(diǎn) 在14天超 出 了2S，未超過(guò)3S。

2.2 工程實(shí)例2

斑點(diǎn)型奇異可視為孤立型奇異的集合，因此孤立型奇異辨識(shí)和斑點(diǎn)型奇異估計(jì)具有分析的代表性。結(jié)合錦屏一級(jí)拱壩3個(gè)垂線點(diǎn)PL11-3，PL11-4和PL13-3的徑向變形監(jiān)測(cè)資料，通過(guò)人為構(gòu)造奇異成分［18］，檢驗(yàn)PCA-CA方法的性能。測(cè)點(diǎn)布置如圖5所示，圖6為徑向變形原位監(jiān)測(cè)過(guò)程線。

圖5 11#，13#壩段垂線監(jiān)測(cè)點(diǎn)布置Fig.5 Vertical monitoring points of the 11# and 13# dam sections

圖6 PL11-3，PL11-4和PL13-3的徑向變形監(jiān)測(cè)序列Fig.6 Time series of radial deformation of PL11-3,PL11-4 and PL13-3

2.2.1 PCA辨識(shí)性能分析

圖7為PCA方法的診斷結(jié)果。依據(jù)PCA辨識(shí)步驟1～5，得到了t1～t3效應(yīng)量，如圖7（a）～（c）所示，測(cè)點(diǎn)變形時(shí)間序列的相關(guān)系數(shù)矩陣如表2所示，t1～t3的特征值、解釋能力以及累積解釋能力如表3所示。依據(jù)步驟6，可得到SPEk過(guò)程線，在執(zhí)行步驟7時(shí)，α取為0.01，計(jì)算得到控制限SPE0.01=0.11，由圖7（d）可知，SPEk＜SPE0.01，因此3個(gè)測(cè)點(diǎn)的變形原位監(jiān)測(cè)數(shù)據(jù)不包含奇異成分。

圖7 PCA方法的診斷結(jié)果Fig.7 Diagnosis results of the PCA method

表2 相關(guān)系數(shù)矩陣Tab.2 Correlation coefficient matrix

表3 特征值、解釋能力及累積解釋能力Tab.3 Eigenvalues,explanatory capabilities and accumulated explanatory capabilities

PL11-4測(cè)點(diǎn)在2014-10-03的徑向變形監(jiān)測(cè)值為39.4 mm，將其構(gòu)造為不同程度的孤立型奇異成分，記作0#～6#情況，分別為39.6，39.9，40.4，40.9，41.4，42.4和43.4 mm。圖8為1#～5#情況下SPEk與SPE0.01的關(guān)系圖，圖9為1#～6#情況的CSPE統(tǒng)計(jì)量。由圖8可知，在2#～5#情況下，SPEk超過(guò)了SPE0.01，且 從 圖9可 以 看 出，在2#～6#情 況 下，PL11-4測(cè)點(diǎn)的CSPE統(tǒng)計(jì)量明顯高于PL11-3和PL13-3，因此2#～6#可判定為孤立型奇異成分。表4比較了PCA準(zhǔn)則、拉依達(dá)準(zhǔn)則、狄克松準(zhǔn)則和t準(zhǔn)則的性能，可以看出，4種方法分別可辨識(shí)出相對(duì)誤差為3.81%，7.61%，7.61%和5.08%的孤立型奇異成分。

表4 4種方法性能對(duì)比Tab.4 Performance comparison of four methods

圖8 1#～5#情況下SPEk與SPE0.01的關(guān)系Fig.8 Relations between SPEk and SPE0.01 under 1# case～5# case

圖9 1#～6#情況的CSPE統(tǒng)計(jì)量Fig.9 CSPE statistics of 1# case～6# case

2.2.2 CA估計(jì)性能分析

假設(shè)PL11-4測(cè)點(diǎn)在2016-07-15至2016-09-15期間的變形監(jiān)測(cè)數(shù)據(jù)為斑點(diǎn)型奇異，利用2014-01-01至2016-07-14的原位監(jiān)測(cè)資料建立協(xié)整模型，對(duì)奇異成分進(jìn)行估計(jì)。圖10為3個(gè)測(cè)點(diǎn)徑向變形散點(diǎn)關(guān)系。由圖可知，式（16）的最高次數(shù)aq均取1，利用逐步回歸法，可得

圖10 PL11-3，PL11-4和PL13-3測(cè)點(diǎn)徑向變形散點(diǎn)關(guān)系Fig.10 Scatter relations among radial deformation of PL11-3,PL11-4 and PL13-3

圖11為CA方法的估計(jì)結(jié)果。PL11-4的監(jiān)測(cè)值、擬合值和余量序列如圖11（a）所示，復(fù)相關(guān)系數(shù)為0.997 1。ADF檢驗(yàn)結(jié)果顯示：t統(tǒng)計(jì)量的值為-2.724 4，1%和5%兩個(gè)顯著性水平的臨界值分別為-2.400和-2.150，因此t統(tǒng)計(jì)量小于臨界值，即原假設(shè)均被拒絕，余量序列ε平穩(wěn)。圖11（b）為自回歸模型、統(tǒng)計(jì)模型和CA模型對(duì)奇異成分的估計(jì)結(jié)果，3種模型的復(fù)相關(guān)系數(shù)分別為0.743 2，0.871 5和0.994 5。

圖11 CA方法的估計(jì)結(jié)果Fig.11 Estimation results of the CA method

3 結(jié)論

1）為消除奇異成分對(duì)大壩變形安全監(jiān)控結(jié)論客觀性的影響，運(yùn)用PCA理論，建立了SPE和CSPE，借助假設(shè)檢驗(yàn)，提出了PCA辨識(shí)準(zhǔn)則，在此基礎(chǔ)上，基于CA原理，應(yīng)用ADF檢驗(yàn)和逐步回歸法，提出了CA似然估計(jì)模型。利用陳村拱壩和錦屏一級(jí)拱壩徑向變形原位監(jiān)測(cè)資料，檢驗(yàn)了PCA-CA方法的有效性。

2）拉依達(dá)準(zhǔn)則、狄克松準(zhǔn)則和t準(zhǔn)則以1維時(shí)間序列為分析對(duì)象，辨識(shí)性能對(duì)樣本分布型式的依賴(lài)性較大，僅體現(xiàn)了概率意義。PCA準(zhǔn)則針對(duì)多維數(shù)據(jù)序列進(jìn)行分析，不僅包含概率含義，而且考慮了多測(cè)點(diǎn)變形的主成分關(guān)系，辨識(shí)性能最佳。

3）自回歸模型和統(tǒng)計(jì)模型著重分析1維時(shí)間序列，其中：自回歸模型僅考察了前期監(jiān)測(cè)數(shù)據(jù)的時(shí)間波動(dòng)規(guī)律，外延性較差；統(tǒng)計(jì)模型刻畫(huà)了水壓、溫變和時(shí)效因素，外延性?xún)?yōu)于自回歸模型，但屬于半經(jīng)驗(yàn)?zāi)Ｐ?，估?jì)精度仍不足；CA模型以多維數(shù)據(jù)序列為分析對(duì)象，既刻畫(huà)了變形隨時(shí)間的波動(dòng)特征，亦表征了多測(cè)點(diǎn)變形間的協(xié)整關(guān)系，同時(shí)也反映出外部作用的綜合效應(yīng)，外延性較好，取得了最佳的估計(jì)精度。

4）PCA-CA方法要求：①t1效應(yīng)量對(duì)原始變形序列X的解釋程度e1盡可能大，通過(guò)工程實(shí)例1可知，當(dāng)e1=95.23%時(shí)，PCA準(zhǔn)則可有效辨識(shí)出奇異數(shù)據(jù)；②余量序列平穩(wěn)。對(duì)于我國(guó)已修筑的部分老壩，布設(shè)的變形監(jiān)測(cè)儀器較少，不易滿(mǎn)足上述條件，此時(shí)，仍需采用常規(guī)方法。