基于面板數(shù)據(jù)分析方法的混凝土拱壩變形數(shù)據(jù)時空聚類模型

胡添翼

(上海市政工程設(shè)計研究總院(集團(tuán))有限公司 水利水運(yùn)設(shè)計研究院,上海 200092)

1 研究背景

拱壩作為最常見的混凝土壩類型之一，是一種高次超靜定空間殼體結(jié)構(gòu)[1]，變形狀態(tài)非常復(fù)雜，其安全狀況非常值得關(guān)注。為了解拱壩安全狀況，科研人員常對拱壩進(jìn)行應(yīng)力和變形監(jiān)測，其中以變形監(jiān)測數(shù)據(jù)最為直觀可靠[2]。隨著大壩安全監(jiān)控理論發(fā)展，關(guān)于拱壩變形的分析和預(yù)測陸續(xù)涌現(xiàn)出許多成果[3-5]，如統(tǒng)計模型、混合模型等。目前針對拱壩的安全監(jiān)測已經(jīng)實現(xiàn)了對整個壩區(qū)的監(jiān)測，在拱壩壩體及周圍布設(shè)有大量變形監(jiān)測點，如小灣拱壩安全監(jiān)測系統(tǒng)已成為世界上最大的安全監(jiān)測系統(tǒng)之一[6]。大壩海量的變形監(jiān)測數(shù)據(jù)具有廣闊的數(shù)據(jù)挖掘空間。

面板數(shù)據(jù)，指在時間序列上取多個截面，在這些截面上同時選取樣本觀測值所構(gòu)成的樣板數(shù)據(jù)；同時考慮樣本時空性質(zhì)的面板數(shù)據(jù)被稱為時空面板數(shù)據(jù)。在多測點監(jiān)測狀況下，高拱壩測點變形數(shù)據(jù)同時包含多測點截面數(shù)據(jù)和單測點時間序列，具有時間、空間2種維度，因此可以被看作是一種時空面板數(shù)據(jù)[7-9]。時空面板數(shù)據(jù)聚類特性的挖掘和時空劃分是時空數(shù)據(jù)挖掘的重要課題之一，其主要內(nèi)容為根據(jù)時空數(shù)據(jù)的一些特征，對空間中的對象進(jìn)行聚類或劃分，使具有相似特征的對象能夠被劃分到同一個類別中去。近年來，面板數(shù)據(jù)分析方法已經(jīng)在經(jīng)濟(jì)、農(nóng)業(yè)、工程等領(lǐng)域獲得發(fā)展，在大壩安全監(jiān)控領(lǐng)域，傳統(tǒng)的聚類分析已經(jīng)得到了應(yīng)用，如施玉群等[10]采用聚類分析方法描述多效應(yīng)量之間的關(guān)聯(lián)關(guān)系，建立了相應(yīng)的融合診斷模型；陳繼光[11]基于灰色關(guān)聯(lián)度計算和模糊聚類分析原理，對變形數(shù)據(jù)進(jìn)行了動態(tài)分析；陳毅[12]建立了基于基因表達(dá)式編程的大壩位移強(qiáng)度聚類分析模型，對大壩的整體變形規(guī)律進(jìn)行分析。但這些研究一般僅僅考慮變形的時間性或者空間性，沒有同時從時、空2個方面進(jìn)行考量，不是典型的時空聚類分析。

本文基于時空數(shù)據(jù)挖掘領(lǐng)域的聚類方法，分析混凝土拱壩變形序列的變化過程，提取了變形序列的相似性特征；提出了混凝土壩變形數(shù)據(jù)不同時間截面、不同測點變形序列的3種相似性指標(biāo)(絕對距離、增量距離、增速距離)，并計算了相應(yīng)的綜合距離指標(biāo)，由此定量分析變形時間截面和變形空間測點的相似程度；利用Ward聯(lián)結(jié)聚類方法，對混凝土壩變形時段及其相應(yīng)變形區(qū)域進(jìn)行劃分；最終建立基于面板數(shù)據(jù)分析方法的高混凝土壩變形測點聚類分析模型，對混凝土拱壩變形狀況進(jìn)行綜合分析。結(jié)合工程實例，對時空聚類模型的運(yùn)用效果進(jìn)行驗證。

2 混凝土拱壩變形相似性特征

2.1 相似性特征選取

假設(shè)混凝土壩上3個不同位置的測點分別為A、B、C，其變形數(shù)值隨時間的變化過程如圖1所示。

圖1 不同變形數(shù)據(jù)的變化趨勢Fig.1 Trends of deformation data at different monitoring points

圖1中，如果單純考慮每個測點的變形過程和變化趨勢，A點和C點是最為接近的，應(yīng)該被歸為一類；如果單純考慮變形的初始值和最終值(以及相應(yīng)的平均值大小)，A點和B點雖然變化趨勢不完全相同，但初始最終的變形值相同，應(yīng)該被歸為一類；如果從更大的范圍來看，3個測點的變形值都有增大的趨勢，也有可能被聚類為同一類。事實上，合理聚類方法應(yīng)該綜合考慮幾種不同的變形特征，包括各個測點時段內(nèi)變形值的變化數(shù)值、變化速率和變化趨勢等各種變形特征需要被綜合考慮。

綜上，假設(shè)測點i在t時刻的變形為δit，本文著重考慮測點變形下面幾個方面的變化特征。

(1)測點i在t時刻的變形值大小，記為xit，防止變形值差異過大的時間截面/測點被聚到一類,即

xit=δit。

(1)

(2)測點i在t時刻的數(shù)值變化，記為yit，防止變化值差異過大的時間截面/測點被聚到一類,即

yit=xit-xi,t-1=δit-δi,t-1。

(2)

(3)測點i在t時刻的數(shù)值變化幅度，記為zit，防止變化趨勢差異過大的時間截面/測點被聚到一類,即

(3)

2.2 相似性特征的標(biāo)準(zhǔn)化

容易發(fā)現(xiàn)前文提取的3種變形特征量單位和數(shù)量級不一致，如果直接相加求所謂的“綜合聚類指標(biāo)”顯然是不盡合理的，因此計算綜合相異性指標(biāo)之前需要對xit、yit、zit進(jìn)行標(biāo)準(zhǔn)化操作，本文采用Z-score方法對基礎(chǔ)指標(biāo)進(jìn)行標(biāo)準(zhǔn)化操作。假設(shè)共有N個監(jiān)測點、T個時間截面的變形監(jiān)測數(shù)據(jù)。以變形值xit為例，計算出所有xit的算術(shù)平均值(數(shù)學(xué)期望)，記為μx；計算所有xit的標(biāo)準(zhǔn)差，記為σx，則Z-score標(biāo)準(zhǔn)化公式為

(4)

標(biāo)準(zhǔn)化后的變量值在0上下波動，平均數(shù)為0，標(biāo)準(zhǔn)差為1。經(jīng)過標(biāo)準(zhǔn)化操作之后，數(shù)據(jù)的取值范圍得到了進(jìn)一步限制，從而避免了量綱和取值范圍所帶來的影響。

3 混凝土拱壩變形時空相似性指標(biāo)

本文運(yùn)用聚類分析方法，一方面是為了對時間截面的測點變形狀態(tài)進(jìn)行聚類,從而實現(xiàn)拱壩變形時段劃分；另一方面是為了在時段劃分的基礎(chǔ)上，根據(jù)測點變形序列對空間測點進(jìn)行聚類從而實現(xiàn)變形區(qū)域劃分。根據(jù)測量方式的不同，衡量變量的尺度一般分為間隔尺度、名義尺度和有序尺度[7]。考慮拱壩變形監(jiān)測通常采用每周一測或每日一測等監(jiān)測間隔，本文對于變形面板數(shù)據(jù)聚類分析的探討僅僅基于間隔尺度。

假設(shè)用δit表示t時刻測點i的變形指標(biāo)，用dij表示測點i和測點j變形之間的相似程度大小。連續(xù)變量常用的相似性指標(biāo)有歐氏距離、歐氏距離的平方(簡稱SED)、絕對值距離和相關(guān)系數(shù)等。其中，SED因為不涉及到開方等計算，在計算上相對歐氏距離更加簡便，本文所有的相異性指標(biāo)均基于SED。以測點i和測點j在t時刻的變形值為例，其SED為

dijt(SED)=(δit-δjt)2。

(5)

式中δit和δjt分別表示測點i和j在t時刻的變形值。

3.1 時間序列劃分指標(biāo)

3.1.1 基本指標(biāo)

為得到更加合理的相似性指標(biāo)，本文從變形值、變形值的波動大小及其波動幅度3方面出發(fā)，定義了3種基本相似性指標(biāo)，即不同測點的橫截面絕對距離、橫截面增量距離和橫截面增速距離，詳情如下。

(6)

式中：xnp=δnp；xnq=δnq。

(7)

式中：ynp=xnp-xn,p-1；ynq=xnq-xn,q-1。

(8)

3.1.2 綜合指標(biāo)

(9)

式中α1、α2、α3分別表示3種基本指標(biāo)的權(quán)重系數(shù)，α1+α2+α3=1。

3.2 空間測點聚類指標(biāo)

3.2.1 基本指標(biāo)

參考針對時間序列劃分的基本相似性指標(biāo)，這里定義3種基本相似性指標(biāo)，即測點聚類的測點絕對距離、測點增量距離和測點增速距離，詳情如下。

(10)

式中：xet=δet；xft=δft。

(11)

式中：yet=xet-xe,t-1；yft=xft-xf,t-1。

(12)

3.2.2 綜合指標(biāo)

(13)

式中β1、β2、β3分別表示3種基本指標(biāo)的權(quán)重系數(shù)，β1+β2+β3=1。

式(9)和式(13)中的權(quán)重系數(shù)可以根據(jù)實際情況進(jìn)行主觀或客觀賦值，本文采用客觀賦權(quán)熵權(quán)法來計算式(9)和式(13)中的權(quán)重系數(shù)[13-14]。

4 混凝土拱壩時空聚類模型

常見的聚類分析方法分為劃分聚類和層次聚類。對混凝土拱壩變形數(shù)據(jù)進(jìn)行聚類分析，因為事先不能確定聚類的數(shù)量，所以劃分聚類方法并不適用于壩體測點聚類，層次聚類方法會相對更加合適一些。在實際應(yīng)用中，Ward聯(lián)結(jié)法的分類效果相對較好，應(yīng)用也比較廣泛，本文選用該方法進(jìn)行聚類。聚類內(nèi)容分為2部分：一部分是將時間序列離散化，根據(jù)每個時間截面大壩測點的變形情況對時間截面進(jìn)行聚類，從而對變形時間序列進(jìn)行有效劃分；另一部分是將空間變形區(qū)域離散化，根據(jù)每一個劃分時段測點的變形序列對測點進(jìn)行聚類，從而對大壩的變形區(qū)域進(jìn)行劃分。

4.1 時間序列劃分方法

(14)

(15)

當(dāng)k確定時，要選擇使W達(dá)到極小的分類。

基于Ward方法思想，對于包含T個時期、N個測點的大壩變形數(shù)據(jù)，可得到面板數(shù)據(jù)分類Gl內(nèi)部Tl個不同時間截面的離差平方和為

(16)

若分為k類，則全局總離差平方和為

(17)

4.2 空間測點聚類方法

(18)

(19)

當(dāng)v確定時，要選擇使W達(dá)到極小的分類。

基于Ward方法思想，對于包含T個時期，N個測點的大壩變形數(shù)據(jù)，可得到面板數(shù)據(jù)分類Hu內(nèi)部Nl個測點離差平方和為

(20)

若分為v類，則全局總離差平方和為

(21)

4.3 聚類數(shù)目的確定

聚類數(shù)目需要綜合考慮譜系聚類樹狀圖和測點變形狀況進(jìn)行最后的確定，實際操作中，往往通過閾值法確定層次聚類的數(shù)目。假設(shè)聚類過程共進(jìn)行h次合并，第g次與最后一次分區(qū)的區(qū)域距離之比為

Sg=Wg/Wh-1。

(22)

式中：Wg為第g次分區(qū)全局總離差平方和；Wh-1為最后一次分區(qū)全局總離差平方和。

如果Sg與Sg+1相差較小，g與Sg-1相差較大，相應(yīng)的Wg可以作為分區(qū)的閾值，這樣就可以確定聚類的數(shù)目。

4.4 模型運(yùn)作流程

本模型基于Ward層次聚類分析方法，主要分為時段劃分和空間聚類2個部分。首先是變形時段劃分，根據(jù)每個時間截面所有測點的變形相似性特征，對混凝土拱壩變形時間截面進(jìn)行聚類分析，從而實現(xiàn)對變形序列的有效時段劃分。在時段劃分之后，基于每一個時段內(nèi)測點變形序列的變化情況，再對混凝土拱壩的空間測點進(jìn)行聚類分析，從而可以在空間上對混凝土拱壩的變形區(qū)域進(jìn)行劃分。模型運(yùn)作流程如圖2所示。

圖2 模型運(yùn)作流程Fig.2 Process of model operation

通過時空聚類分析模型，可以識別混凝土拱壩變形的相似時段和相似性區(qū)域，從而了解混凝土拱壩變形狀態(tài)。

5 工程實例

為驗證前文提出的時空聚類模型，以國內(nèi)某300 m級混凝土拱壩變形情況為例進(jìn)行分析。

5.1 工程概況

該拱壩壩體共設(shè)置了34個垂線測點監(jiān)測大壩安全狀況，測點具體布設(shè)情況如圖3所示。

圖3 某高拱壩變形測點布置示意圖Fig.3 Layout of deformation monitoring points for an arch dam

以該拱壩34個測點從2014年5月1日開始的時間截面為300 d的徑向變形數(shù)據(jù)和相關(guān)衍生指標(biāo)共計30 600(3×300×34)個數(shù)據(jù)為基礎(chǔ)，進(jìn)行聚類分析。根據(jù)相關(guān)文獻(xiàn)，庫水位是拱壩位移變化的最主要因素，本文主要考慮水位變動條件。由于水位變動條件下拱壩變形大小受壩前水深影響最為顯著[15-16]，本文主要考察壩前水深變化過程中變形數(shù)據(jù)的聚類情況，壩前水深和所有測點的平均變形隨時間變化過程線如圖4所示。

圖4 水位變化及變形過程線Fig.4 Process lines of water level fluctuation and deformation

5.2 變形時段劃分

首先對拱壩的變形時段進(jìn)行劃分。對不同時間截面的變形狀態(tài)進(jìn)行聚類分析，得到聚類譜系樹狀圖，如圖5所示(因時間截面太多，圖中僅取10類以內(nèi)的聚類結(jié)果)。

圖5 拱壩變形序列劃分譜系聚類樹狀圖Fig.5 Tree plot of deformation sequence clustering for the arch dam

從聚類譜系樹狀圖可以得知，變形時段可劃分為3個主要時段：第1時段為2014年5月1日—2014年9月9日，此段時間內(nèi)上游水深較深，大壩的整體變形比較平穩(wěn)；第2時間段為2014年9月10日—2014年12月28日，此段時間內(nèi)水深下降，大壩變形隨水深下降不斷減??；第3時間段為2014年12月29日—2015年2月24日，此段時間內(nèi)水深較淺，大壩變形再次趨于穩(wěn)定。模型劃分的3個時段如圖6所示。

圖6 拱壩變形時間序列劃分結(jié)果Fig.6 Results of time series division of arch dam deformation

從圖6可以發(fā)現(xiàn)，本模型劃分出來的3個變形時段和大壩的整體變形情況以及上游水深的變化情況高度吻合，充分體現(xiàn)出每一個時段自身的特點，說明本模型對變形時段的劃分是比較合理的。

5.3 空間測點聚類

根據(jù)前面聚類劃分出來的3個變形時間段，分別在每一段時間段內(nèi)根據(jù)測點的變形特征對測點進(jìn)行聚類分析。

5.3.1 高水位期

針對高水位運(yùn)行的第1時段，空間測點聚類結(jié)果如圖7所示。

圖7 拱壩高水位變形聚類分區(qū)Fig.7 Arch dam’s deformation zone under high water level

從測點聚類的情況來看，主要出現(xiàn)了4個比較明顯的變形相似區(qū)域，其中Ⅰ區(qū)變形最大，Ⅱ區(qū)次之，Ⅲ區(qū)變形略小于Ⅱ區(qū)，Ⅳ區(qū)變形最小。變形最大區(qū)域集中在拱冠梁附近，呈長條狀分布；其它聚集區(qū)域以壩軸線為中心對稱分布。

5.3.2 水位消落期

針對水位不斷消落的第2時段，空間測點聚類結(jié)果如圖8所示。

圖8 拱壩水位消落期變形聚類分區(qū)Fig.8 Arch dam’s deformation zone during water level drawdown

從測點聚類的結(jié)果來看，拱壩被分為3個變形相似區(qū)域，其中Ⅰ區(qū)變形最大，Ⅱ區(qū)次之，Ⅲ區(qū)變形最小。變形最大區(qū)域為大壩的中部；變形區(qū)域分布的對稱性較好，但是Ⅱ區(qū)局部測點似乎出現(xiàn)了變形過大的問題。由此可以發(fā)現(xiàn)，隨著上游水深的下降，拱壩變形最大區(qū)域已經(jīng)開始出現(xiàn)下移，說明拱壩的應(yīng)力狀態(tài)已經(jīng)出現(xiàn)變化。

5.3.3 低水位期

針對水位消落后低水位運(yùn)行的第3時段，空間測點聚類結(jié)果如圖9所示。

圖9 拱壩低水位變形聚類分區(qū)Fig.9 Arch dam’s deformation zone under low water level

從測點聚類的結(jié)果來看，拱壩被分為3個變形相似區(qū)域，其中Ⅰ區(qū)變形最大，Ⅱ區(qū)次之，Ⅲ區(qū)變形最小。相較于上一個時段，變形最大的區(qū)域已經(jīng)轉(zhuǎn)移至拱壩的中下部，且該區(qū)域的面積有所減??；變形區(qū)域分布仍然呈高度對稱，且之前變形過大的測點的變形已經(jīng)有所減小。

6 結(jié) 論

本文根據(jù)3種變形數(shù)據(jù)的主要變化特征，從時段劃分和測點聚類角度，構(gòu)建了基于歐氏距離平方的3種基礎(chǔ)相似性指標(biāo)，在此基礎(chǔ)上，提出了對應(yīng)的綜合相異性指標(biāo)及相應(yīng)的計算方法；基于Ward聯(lián)結(jié)聚類法，分別提出了針對混凝土拱壩變形時段劃分和變形空間測點聚類的聚類方法和步驟，建立了混凝土拱壩變形時空聚類模型；通過對某混凝土拱壩變形的時空聚類分析，驗證了該模型的有效性。

根據(jù)工程實例的計算結(jié)果，可以看出：

(1)本文提出的混凝土拱壩變形特征以及由此提出的聚類相似性指標(biāo)是比較合理的，能夠綜合反映大壩的變形狀態(tài)；所采用的層次聚類方法計算簡單快捷，可以對大壩的變形序列進(jìn)行時段劃分，對每個變形時段的空間測點進(jìn)行聚類，從而了解大壩在不同時段的變形狀態(tài)。

(2)不同時段特別是不同水位條件下的聚類分析會有比較明顯的差異。從變形區(qū)域的分布情況來看，高水位和低水位等水位變化較小的時段區(qū)域分布的對稱性良好，水位下降等水位變化時段區(qū)域分布的對稱性相對差。因此，需要重點關(guān)注水位變化時期大壩的安全狀況。