基于分形插值的城市降雨時間序列精細化估算

孫欣躍

（遼寧省水資源管理集團有限責任公司，沈陽110003）

分形插值是基于分形基礎理論的具體應用研究，能很好地描述生活中廣泛存在的對象，往往具有復雜性和粗糙性。除此之外，分形插值函數(shù)靈活性強和穩(wěn)定性高的特點使其廣泛應用于水文地質、 動畫仿真、圖像處理等眾多研究領域［1］。

1 研究區(qū)概況

濟寧市山東省西南部的一個地級市， 濟寧市氣候四季分明，氣溫宜人。受熱帶海洋氣團或熱帶海洋氣團影響， 天氣炎熱多雨， 年平均降水量597～820mm；冬季多為偏北風，受極地大陸氣團影響，多為晴冷天氣；春秋兩季為大氣環(huán)流調整期，春季易干旱多風，變化迅速溫暖；秋天涼爽，有時下雨。濟寧市的水資源較為充沛，年降水量大，對于防洪減災工作難度較大。

2 計算方法和數(shù)據(jù)

2.1 計算原理

分形插值是對具有分形特征的曲線進行的快速插值方法，其根據(jù)分形曲線的自相似性，基于迭代函數(shù)系統(tǒng)（IFS）實現(xiàn)。 首先要構造一個迭代函數(shù)系，以此來表示采用分形的插值，令吸引子A為插值函數(shù)的圖形。對于{（xi，yi），i=0，1，···，N}∈R2，分形插值函數(shù)是 插 值 于{（xi，yi），i=0，1，···，N}的 一 個 連 續(xù) 函數(shù)f［x0，xn］→R，且其圖形就是雙曲迭代函數(shù)系統(tǒng)的吸引子A。

創(chuàng)建一個吸引子為插值函數(shù)f（x）的迭代函數(shù)系{R2：Wi，i=1，2···，N}，設函數(shù)系中的每個函數(shù)Wi是仿射變換，表示為［2］：

Wi∶R2→R2，i=1，2，…，N

2.2 降雨時間序列的提取

氣象數(shù)據(jù)由中國氣象信息中心官網(wǎng)下載， 數(shù)據(jù)類型為某氣象站點氣候資料的日值數(shù)據(jù)。 以數(shù)值形式的氣象資料為每日平均資料， 時降雨資料可基于ArcMap于時刻降雨資料雨量圖中提取， 逐小時雨量分布如圖1。

圖1 逐小時雨量分布

每幅氣象圖代表某一小時的氣象數(shù)據(jù)， 在某一具體氣象站點多幅圖的數(shù)據(jù)的提取就構成了該點某時間區(qū)間的序列。 本文采用的數(shù)據(jù)是小時時間間隔的降雨數(shù)據(jù)。 下載2018年6月30日12：00到2018年7月5日12：00時刻雨量圖， 選取或創(chuàng)建指定氣象站點（濟寧為54326站點）的Shape 文件，提取該點（或多點）的時間序列數(shù)據(jù)。 其中點Shape數(shù)據(jù)為設定的氣象數(shù)據(jù)提取位置。通過Spitial Analyst 中的提取分析工具，提取圖像素數(shù)值，建立換算表，取圖例中的中間值，得到的換算如表1。

表1 降雨數(shù)據(jù)提取對應

最終得到濟寧市（氣象站點58326）的120h連續(xù)降雨值，轉換為降雨強度單位，降雨時間序列如圖2。

圖2 120h降雨時間序列

3 研究結果

3.1 降雨數(shù)據(jù)分形特性分析

分形分析采取盒子覆蓋法計算曲線的盒維數(shù)，通過繪制單位長度和覆蓋所需的總正方形數(shù)量的雙對數(shù)曲線得到［3］。選取若干不同尺度的r，一般情況下選取離散點時間序列樣點周期的整數(shù)倍。 用邊長為r的正方形對曲線進行覆蓋，計算不同尺度r覆蓋時間序列曲線所需單位正方形的總數(shù)N（r）。繪制ln N（r）-ln（1/r）的雙對數(shù)散點圖，采取最小二乘法估計擬合出散點斜率。 回歸直線與散點的擬合就反應了降雨時間序列的自相似特性［4］。

根據(jù)分形維數(shù)的盒維數(shù)計算法， 采用MATLAB軟件引入這些公式，改變小正方形的邊長r，給定的單位長度r為兩個時間單位（r=2），單位長度為r=2i/N（i=1，2，…，5），統(tǒng)計正方形數(shù)目N（r）。

表2 r覆蓋統(tǒng)計

將表中的ln（1/ri）和ln（N（ri））數(shù)據(jù)值，添加趨勢擬合線，可得到趨勢擬合的回歸直線。得到的回歸方程為：y=1.012x+6.069。

圖3 時間序列曲線盒覆蓋雙對數(shù)

根據(jù)分形盒維數(shù)的定義與原理， 得知該擬合直線的斜率D即為時間序列曲線的盒維數(shù)估計數(shù)值。所以維數(shù)為1.012。 雙對數(shù)曲線的擬合相關系數(shù)R2為0.998。 從擬合直線的直觀表象來看， 對數(shù)點分布在直線的附近，擬合相關率高達99.8%，說明此次降雨過程的時間序列曲線具有很強的局部與整體的相似性，定量反映了該時間序列曲線具有分形特性。

3.2 降雨時間序列的分形插值計算

由于對缺少實測降雨數(shù)據(jù)， 中國氣象數(shù)據(jù)網(wǎng)站提供的地面逐小時降雨資料， 所以可依據(jù)時降雨時間序列的分形特性， 采用分形插值函數(shù)計算提取得到120個小時每兩個時刻之間的降雨量。

對時刻進行排序，即構造120個插值點為p0（x0，y0）=（0，0.10）、p1（x1，y1）=（1，0.10）、p2（x2，y2）=（2，0.10）···p119（x119，y119）=（119，0.10）。 上述120個點是實現(xiàn)仿射分形插值函數(shù)構造的所需插值點，根據(jù)垂直比例因子非確定的取值方法，可初步得到各個仿射系數(shù)［5］。

以上計算過程可基于MATLAB軟件程序文件實現(xiàn)，對已知初始點進行迭代，根據(jù)分形插值函數(shù)的計算規(guī)律，需要進行一次迭代過程，即在原兩個數(shù)據(jù)點之間，插入118個數(shù)據(jù)點，可將原來以小時為單位的降雨時間序列細化到30s，得到插值后曲線，如圖4。

圖4 分形插值后降雨時間序列

插值得到細化時間序列將插值分形曲線圖與實測曲線圖進行對比可看出，曲線走勢相同，變化規(guī)律基本一致。 表3是插值數(shù)據(jù)和實測數(shù)據(jù)的平均值、方差和分形維數(shù)的比較。

表3 原始數(shù)據(jù)與插值數(shù)據(jù)比較

式中 m為插值數(shù)據(jù)；x為實測數(shù)據(jù)；r為相對誤差。

3.3 分形插值與傳統(tǒng)插值方法結果對比

與傳統(tǒng)插值方法進行對比， 可看出分形插值在插值效果上具有明顯優(yōu)越性， 為120個數(shù)據(jù)線性插值的結果。 線性插值沒有考慮已知連點間的具體變化過程，將兩點直接用直線進行連接，而分形插值基于部分與整體的相似性， 將兩點間變化的過程具體化，反映了局部特征。拉格朗日插值也屬于一種常用的傳統(tǒng)插值方法，又稱為多項式插值。該插值方法的缺點在于插值節(jié)點個數(shù)的增加的同時， 不能確保兩個插值節(jié)點之間插值函數(shù)一定能夠很好地與被插值函數(shù)吻合。從理論角度來看，冪次數(shù)多項式需要的數(shù)據(jù)多，擬合結果比低冪次準確，但從舍入誤差看，高次插值計算過程復雜，而且計算量較大，計算過程中會產(chǎn)生較大程度的誤差累積。擬合過程中出現(xiàn)了“龍格現(xiàn)象”， 使得靠近區(qū)間端點處的擬合數(shù)值很大，曲線的變化速度較高，出現(xiàn)了較大誤差。而分形插值法的誤差值較小，所得結果較為精準如圖5。

圖5 線性插值后降雨時間序列

4 結語

本章從分形理論與分形的基本性質出發(fā)， 對濟寧市2018年6月30日12：00 至2018年7月5日12：00連續(xù)120h的降雨時間序列進行了分形性分析， 得到該降雨時間序列曲線的分形盒維數(shù)為1.012，雙對數(shù)曲線擬合度為99.8%，說明了該曲線具有很好的分形特性。進而應用分形插值理論，對該分形特性曲線運用迭代系統(tǒng)進行分形插值， 根據(jù)最終所需時間序列的精確度， 確定進行一次迭代， 實現(xiàn)每兩個小時之間118 個插值點的生成， 最終將以小時為單位的時間序列細化處理為30s為單位的新時間序列。 將源數(shù)據(jù)與插值后數(shù)據(jù)進行誤差統(tǒng)計，平均值、均方差和分形維數(shù)的誤差分別為0.007，0.066，0.18，誤差均在可接受范圍內， 可滿足當前中小型城市的城市雨洪建模的降雨數(shù)據(jù)需求。