劉秀林 張行南 方園皓 黃晴
摘要:隨著水文自動測報技術的發(fā)展,越來越多的自記雨量站投入使用,而研究自記雨量站的數據質量對水利工程的運行具有重要意義。采用格拉布斯準則和K-medoids聚類方法對金沙江下游流域收集的2008~2015年遙測雨量站逐時降雨數據系列進行數據質量分析研究。結果顯示,該方法在年累積雨量異常值的確定以及逐時降雨數據奇異點的尋找方面效果顯著,能快速判別異常值和奇異點,且有統計學理論基礎。所提出的方法為顆粒度越來越細的水文數據質量分析提供了一定參考依據。
關鍵詞:遙測雨量;逐時降雨數據;格拉布斯準則;K-medoids聚類;奇異點;金沙江下游流域
中圖法分類號:P33
文獻標志碼:A
DOI:10.16232/j.cnki.1001-4179.2019.03.023
隨著我國社會經濟水平的不斷發(fā)展,水利基礎工程設施建設的腳步也越來越快,截至2017年,長江流域的雨量站數量已達9959個,其中包括大量為水情自動測報系統服務的遙測雨量站。遙測雨量站的建設和應用加密了站網密度,滿足了流域實時洪水預報的工程需求,在一定程度上為面雨量計算方法、流域洪水預報的研究提供了大量的數據支持。隨著長江,上游水庫工程的進一步規(guī)劃建設,為滿足水庫的施工期洪水預報對預見期和預報精度以及水庫建成后樞紐信息實時監(jiān)測等計算要求,需進一步加密現有的雨量站網。由此可見,現階段我國的雨量站數目仍將增加,同時,雨量數據的時序顆粒度也越來越細,目前的遙測雨量站數據儲存皆是逐時雨量,一個雨量站每年約8760個值存儲于數據庫中。
雨量數據經遙測雨量站實時傳輸入數據庫,大量的數據未經質量分析,給工程項目和科學研究帶來了極大的不便。單站逐時降雨量時序數據不同于水位或者流量的時間系列,具有連續(xù)性的特點,在分析系列時可通過平均變化量(或者平均變化率)來確定異常點,逐時降雨量數據存在大量的零值,是離散的、隨機的周期性過程。序列多會出現兩類問題:模式異常和點異常。模式異常指在一條序列上與其他模式存在顯著差異的、具有異常行為的模式;點異常指在某段時序區(qū)間內與其他序列點存在顯著差異的、具有異常特征的序列點。這兩類問題在單站逐時降雨量序列中具體表現形式分別為年累積雨量異常和區(qū)間奇異點。
本文以金沙江下游流域自建雨量站(指金沙江下游梯級水電站水情自動測報系統三期建設的站點,以及部分可直接或間接獲取數據的信息共享站點)為例,提出面對大量數據時,降雨量數據質量分析的方法,以期為高效分析降雨量數據質量提供一定的參考依據。
1 研究區(qū)域與數據
1.1 研究區(qū)域
金沙江流域地處青藏高原、云貴高原和四川盆地西部邊緣,呈東西短、南北長的狹長形狀。金沙江下段自攀枝花站起,至四川省宜賓岷江口止,地理位置位于東經100°~105°、北緯24°~29°之間,河長782.4km,落差729m。該河段水量大、落差集中,是長江流域水能資源最豐富的河段,同時也是雨量站點較為密集的區(qū)域。
金沙江下游流域氣候為亞熱帶季風濕潤氣候,多年平均降水量為893mm,水量豐沛穩(wěn)定,暴雨多為兩次以上的連陰雨天氣形成,一次暴雨過程的歷時約為3~6d,最大1d降雨量高達100~200mm。本次研究收集到金沙江下游流域95個自記雨量站點逐時雨量數據,站點分布情況見圖1。
1.2 雨量站點信息
自記雨量站點信息及數據均為三峽集團公司梯調通信中心提供,因篇幅限制,各站點經緯度信息和高程信息均未列出,具體位置可參見圖1。在95個雨量站中,大沙店、尼格、向家壩(專)和龍街(三)站點資料系列為2012~2015年,底壩和細沙站資料系列為2013~2015年,其余均為2008~2015年,數據總量為349651條。
2逐時降雨數據質量分析方法
在一定區(qū)域范圍內,鄰近雨量站之間的降雨量有一定的相關性,故在橫向比較(單站各年份比較)的基礎上,針對相應情況,還可進行縱向比較(鄰近站點相應年份比較),以確定數據質量。WMO《水文實踐指南》第1卷中提出:在溫帶、內陸熱帶山區(qū),最合適的做法是按大約500m高差來規(guī)劃高度帶。針對金沙江下游流域雨量站分布位置高程差較大的特點,將與研究站點高程差500m作為搜索范圍,尋找5個鄰近站點進行對比分析。
2.1 累積雨量異常判別
累積雨量異常表現為相較于其他年份該年總量偏大或者偏小,即該年的逐時降雨量數據可能呈現系統性的偏大或者偏小,屬于時間序列異常研究中的模式異常問題。這類問題需要特別警惕,應認真核實數據的正確性和可信度,因為它會影響水文研究中的水量平衡計算。對于該問題,考慮到本文單站年累積雨量統計值個數在3~8之間,故采用改進的格拉布斯準則來確定異常年份。
格拉布斯準則適用于測量次數較少的情況(3≤n<100),可一次性求出多個異常值。改進的格拉布斯準則是將原準則公式中的平均值用中位數代替,可有效消除同側異常值的屏蔽效應,是更為穩(wěn)健的處理方法。其判別方法如下。
先將樣本從小至大排序為新的系列X=(x,x2,.x,),統計臨界系數G(a,n)的值G,(查臨界值表獲得),然后計算G,G,:
公式
式中,a為顯著性水平,n為測量次數,Xψ為樣本中位數,σ為標準誤差。
若G≥G,且G>Go,則x應予以剔除;若G,≥G且G>Go,則x。應予以剔除;若G 根據以上統計學方法,金沙江下游流域雨量站累積雨量異常分析總體思路為:先對單個站點計算各年的累積雨量,挑出異常年份;再用各年的累積雨量與鄰近站點的同年累積雨量做格拉布斯準則分析。由于累積雨量在年際間相差可能較大,而鄰近雨量站的雨量總體反映了豐水年、平水年或者枯水年,可作為較好的參考,故設定兩者輸出相同的異常年份則判定為需核查的數據。 2.2 奇異點識別 奇異點在本次研究中是指在某段時序區(qū)間內的極大值遠遠大于該段內的其他點,可能存在單點數據有誤的隨機誤差。奇異點的存在并不是錯誤值,而是需要進一步核實是不是正確數據的極大值,以避免后續(xù)研究(如水文分析中的次洪參數率定等)中帶來不必要的影響。 針對逐時降雨量數據存在大量的零值,且各個測量值之間離散相互獨立的特點,采用時間序列異常研究中的基于特征空間的方法來識別奇異點:首先對時間序列進行分段,即分為非汛期和汛期,從分段中提取特征,然后在特征空間中應用無序數據集合中的異常點檢測方法一聚類法來尋找奇異點,分段的思想考慮到了非汛期中的奇異點在汛期時段中顯得平庸的特點,能夠有效檢測出各段奇異點。 聚類法是將數據集根據相似度劃分成若干組的統計學方法,不同組中相似度低,相似度可用距離進行度量。K-medoids聚類法是系統聚類法中最為常用的一種,因其算法簡單、收斂速度快、中心點明確以及局部搜索能力強的優(yōu)點被應用到很多方面。K-medoids算法步驟如下: (1) 針對數據集{y,,y2,.,yn},適當選擇k個樣本作為初始聚類中心2,2,g}; (2) 對每個樣本y;找到離它最近的聚類中心z。,并將其分配到z。所標明的類u; (3) 更新每個類的中心: 公式 (5)如果D值收斂,則返回(z,z,*.,zn,U),并終止算法,否則轉至步驟(2)。 K-medoids算法聚類的顯著缺點是需提前指定分類數目,采用最優(yōu)聚類數的評價指標Silhouette來確定分類數目。該指標反映了聚類結構的類內緊密性和類間分離性,既可用于估計最佳聚類數,也可用于評價聚類質量,Silhouette指標值在范圍內變動,指標值越大表示聚類質量越好,最大值對應的類數為最佳聚類數。Silhouette指標值的計算公式為 公式 式中,a(i)是樣本i與類內所有其他樣本的平均距離,b(i)為樣本i到其他每個類中樣本平均距離的最小值。 初步分析金沙江下游流域雨量站逐時降雨量(也稱降雨強度,mm/h)數據,可以發(fā)現大量的零值數據(無降雨),故初步處理應將零值去掉,形成無零值的數據系列文件;對新的數據系列用K-medoids算法進行分類。實踐中僅將分類數k取為2或3,采用Silhouette指標值進行比較確定最優(yōu)分類。小時降雨強度r的等級劃分標準為:r<2.5mm(小雨),2.5mm≤r<8mm(中雨),8.0mm≤r<15mm(中雨),15mm≤r(暴雨)。由此可看出,逐時降雨量的特點是存在大量的較小值較多的中間值和很少的較大值,奇異點存在于較大值的聚類中10]。本次研究對較大值所在的聚類需要做進一步的分析。資料顯示,最大1h降雨強度極少超過100mm,故在該聚類中直接將100mm以上的數據標記為奇異點。另外,經人工逐一在excel中繪制降雨量柱狀圖發(fā)現,一年中各個分段出現奇異點的概率不大,故將該聚類中的閾值設定為3,若該區(qū)段中較大值的個數不大于3個,那么皆標記為奇異點,若大于3個且小于100個,則采用格拉布斯準則來尋找異常數據(奇異點),找到異常數據則輸出奇異點,尋找不到則輸出無奇異點,認為存在100個及以上的較大值在同一聚類中,屬于正常降雨數據。 根據金沙江下游流域降雨特點,將降雨量系列分為3個區(qū)段,1~3月和11~12月降雨較少,分為第1區(qū)段;4,5月和10月降雨量中等,為第2區(qū)段;6~9月降雨量較多,為第3區(qū)段。根據數據特點,各區(qū)段正常降雨強度閾值設定為:≤15mm,≤25mm和≤35mm。由于相鄰區(qū)域降雨存在相似性,即暴雨出現時間的近同步性,且存在每年汛期偏枯偏澇等不同的特點,故單站奇異點的分析采用同年鄰近雨量站降雨對比尋找異常點,而不是該站各年份之間的比較。 2.3 實現方式 根據統計學方法搜尋金沙江下游流域雨量站逐時降雨量存在的累積量異常和奇異點問題,不僅將人工目估的定性方法轉為定量,使得問題尋找有據可循,而且易于編程實現,節(jié)省了大量的人力時間。研究選用Matlab進行編程,利用內置的utmzone和mfwdtran函數將各雨量站的地理坐標轉換為地圖坐標,以便雨量站之間距離的測算,來搜尋鄰近雨量站;利用內置的K-medoids和Silhouette函數實現聚類分析,完成累積量異常和奇異點問題年份和問題點的提取。應用格拉布斯準則搜尋異常累積雨量時,考慮到累積雨量值分布密集,顯著性水平a取0.05即可,分析奇異點時,由于數據差異較大,離散程度高,顯著性水平a取0.005。 3 結果及討論 3.1 計算結果 利用本文提出的檢測累積雨量異常和奇異點的方法,對金沙江下游95個雨量站點數據進行逐一核查,檢測到2015年龍山村站、2012年大沙店站、2012年八家村水庫站、2011年地索(二)站、2014年封過站、2011年后布列托站、2013年細沙站以及2012年龍街(三)站累積雨量異常。經核查原始數據,其中,大沙店站、八家村水庫站、地索(二)站、細沙站和龍街(三)站均是由于該年建設站點投入使用時間較遲,數據系列大多從9月開始記錄使得年累積雨量異常;龍山村站和淌塘站均是年內累積量異常偏小。這里僅詳細分析龍山村站,在累積雨量表1中,先進行各年份比較,查找到2015年累積雨量(79mm)和2011年累積雨量(1031.5mm)異常,再計算同一年份鄰近站點累積雨量,龍山村站2009,2015年的累積雨量異常,根據兩者輸出相同的異常年份為異常值,判定2015年數據需核查矯正;封過站2014年累積雨量達2555.6mm,逐時降雨數據系統性偏大,金沙江下游流域年雨量一般不超過2000mm,故該站點數據需核查校正。 由于奇異點數量較多,這里僅列出表現異常的第1區(qū)段(1月至3月和11月至12月)奇異點數據表和數量最多的第3區(qū)段(6~9月)奇異點分布直方圖。從表2中可以看出,檢測到的奇異點,100%的降雨強度大于15mm(暴雨等級),30.3%的降雨強度大于50mm,18.2%大于90mm??梢哉J為,這些奇異點均屬異常數據,需詳細核實。圖2為第3區(qū)段奇異點降雨量級直方圖,從95個雨量站共734個年逐時雨量數據中找出164個異常點,可以清晰地看到奇異點逐時降雨量范圍在30~50mm居多,有107個,占比65.2%。實踐中可核查該區(qū)域降雨強度在該段內是否屬于正?,F象,若屬正常,那么采用該方法時正常降雨閾值設定可放寬到50mm;地索(二)站2011年5月22日11時降雨數據為235.9mm,考慮到直方圖間距未加入圖中,該數據需校正。 3.2 討論 (1) 格拉布斯準則。剔除粗大誤差的統計學方法有很多。雖然經證明,格拉布斯準則適用于測量次數n大于3小于100的情況,但是異常值判別標準G(a,n)與顯著性水平a的選取同樣相關。在累積雨量異常判別時,為嚴格起見,取a=0.05。對于該例而言,敏感程度較高,從表2中可看到尋找到了較多的異常值,而該異常值在縱向對比時又為正常值,所以綜合考量,實驗判定橫向縱向對比均為異常數據的值為可疑值。 在奇異點判別時,顯著性水平a=0.005,這主要是考慮了逐時降雨數據離差系數大的特點。但從圖2中可看到挑選出的異常值65.2%聚集在30~50mm之間。這種現象可能預示著正常數據被當作可疑數據被挑選出來,表明較小的顯著性水平依舊無法正確判別,而必須加入閾值設定環(huán)節(jié)才能使判別更加高效。 (2) 聚類分析。在逐時雨量異常數據判別時,聚類分析契合了逐時降雨數據量大且類間差別大,類內差距小的特點,可以高效地將較大值分成一類,相當于剪枝的思想,而將實驗研究對象快速聚焦到奇異點的分析上。但從表2可以看到,區(qū)段劃分具有極強的主觀性,3月中下旬在第1區(qū)段中表現異常的值有一部分小于25mm,若是將該時段劃分到第2區(qū)段,可能屬于正常數據而被淹沒。在實踐中,還需根據研究流域實際降雨特征來劃分區(qū)段以及選取閾值。 (3) 問題核實。用年鑒資料來核查數據,以2013年向家壩(專)站為例,本文方法搜尋到奇異點2013年2月25日21:00降雨36mm,當日降雨36mm,而年鑒當日顯示未降雨;9月3日13:00降雨75.5mm,當日總降雨99.5mm,而年鑒資料當日僅降雨11mm??梢?,為嚴謹使用數據,前期降雨數據質量分析工作非常重要。 本文僅研究極大值出現異常的情況,而未對降雨量小的值進行研究,一方面是方法限制,無法獲得正確的參照值;另一方面是因為儀器測量本身存在一定誤差,較小的降雨量產生異常一般在誤差允許范圍內。 查找到的問題數據除了與年鑒進行比較外,也可與其他可靠的氣象產品比較,將錯誤數據矯正,矯正的方法可直接采用正確數據,無準確參考值的情況下可采用鄰近站點降雨插值校正等。 4 結語 (1) 本文詳細討論了金沙江下游流域逐時降雨數據質量分析方法,結合格拉布斯準則對累積雨量異常年份進行確定,利用K-medoids聚類方法挑選可疑極大值數據,確定奇異點。年累積雨量挑選方法科學易行且準確率高,奇異點分析方法由于數據離散程度較大,需結合閾值選取以避免將正常值作為異常值的錯誤。 (2) 基于特征空間的方法來識別奇異點,關鍵在于分段,需綜合考慮流域降雨特點,將降雨特性相近的時段合并,差異大的時段分開,可提高奇異點提取準確性。 (3) 遙測雨量站逐時降雨數據產生的誤差分系統誤差和粗大誤差兩類,即累積雨量異常和奇異點兩類問題。本文提出的方法是對初步處理大量該類數據的一次探索,以期為越來越精細化的水文數據質量分析提供一定參考依據。 參考文獻: [1]胡海洪.遙測雨量站在青海省中小河流監(jiān)測站網中的應用[J].中國高新技術企業(yè),2013(21):59-61. [2]吳騫,吳紹春.基于離群分析的水位異常識別研究[J].硅谷,2010(24):45. [3]詹艷艷,徐榮聰.時間序列異常模式的K-均距異常因子檢測[J].計算機工程與應用,2009,45(9):141-145. [4]張睿,周建中,肖舸,等.金沙江下游梯級和三峽梯級水電站群聯合調度補償效益分析[J].電網技術,2013,37(10):2738-2744. [5]高俊剛,吳雪,張鐿鋰,等.基于等級層次分析法的金沙江下游地區(qū)生態(tài)功能分區(qū)[J].生態(tài)學報,2016,36(1):134-137. [6]世界氣象組織(WMO)著,趙珂經等譯.水文實踐指南[M].北京:水利電力出版社,1987., [7]熊艷艷,吳先球.粗大誤差四種判別準則的比較和應用[J].大學物理實驗,2010,23(1):66-68. [8]夏寧霞,蘇一丹,覃希.一種高效的K-medoids聚類算法[J].計算機應用研究,2010,27(12):4517-4519. [9]周世兵,徐振源,唐旭清.新的K-均值算法最佳聚類數確定方法[J].計算機工程與應用,2010,46(16):27-31. [10]楊茂林,盧炎生.基于剪枝的海量數據離群點挖掘[J].計算機科學,2012,39(10):152-156. [1]RoySS,RouaultM.Spatialpattermsofseasonalscaletrendsinex- tremehourlyprecipitationinSouthAfrica[J].AppliedGeography,2013(39):151-157. [12]Gwo-FongLin,Ming-JuiChang,Chian-FuWang.ANovel SpatiotemporalStatisticalDownscalingMethodforHourlyRainfall[J].WaterResourManage,2017(31):3465-3489. [13]周國良,高唯清,黃昌興.2016年我國極端暴雨事件淺析[J].中國防汛抗旱,2017,27(1):75-78,87. [14]趙超,包為民,瞿思敏,等.遙測系統降雨觀測粗差修正研究[J].人民長江,2003,34(2):4-5,55. [15]李朋軍.遙測與虹吸雨量計降水數據對比分析[J].水科學與工程技術,2016(2):51-53. [16]楊旭,劉志武,李波.多源降水數據在長江上游流城比較研究[J].長江流域資源與環(huán)境,2016,25(1):131-139. 引用本文:劉秀林,張行南,方園皓,黃晴.金沙江下游遙測雨量站數據質量研究[J].人民長江,2019,50(3):131-135. Study on data quality of hourly rainfall of telemetry rainfall stations in lower reaches of Jinsha River LIU Xiulin',ZHANG Xingnan' 1,2,3,FANG Yuanhao',HUANG Qing* (1. College of Hydrology and Water Resources,Hohai University,Nanjing 210098,China;2. National Cooperative Innovation Center for Water Safety & Hydro-Science,Hohai University,Nanjing 210098,China;3. National Engineering Research Center of Water Resources Efficient Utilization and Engineering Safety,Hohai University,Nanjing 210098,China;4. Geographic and Oceano-graphic Sciences College,Nanjing University,Nanjing 210098,China) Abstract:Along with the development of automatic measurement and forecast technology,more and more automatic rainfallstations are put into operation. Analyzing the quality of data measured by automatic rainfall station is of great significance to the operation of water project. The analysis on the hourly rainfall series from telemetry rainfall stations in the lower reaches of the Jinsha River from 2008 to 2015 was carried out by using the Grubbs criterion and k-medoids clustering algorithm. Based on the statistical theory,the results showed that the methods have a remarkable effect in determining the abnormal value of annual accumulated rainfall and the singular points of hourly rainfall data. The method provides some references for hydrological data qualityanalysis that is changed by increasingly fine particles. Key words:telemetry rainfall data;hourly rainfall series;Grubbs criterion;k-medoids clustering;singular points;lower reaches of Jinsha River