基于規(guī)劃求解的暴雨強(qiáng)度公式推求方法研究

李 莉， 冷 藝， 黃 俊， 張 賽， 王子奇， 鄒浩東， 洪毅怡暉， 符榮松

(1.重慶大學(xué) 三峽庫(kù)區(qū)生態(tài)環(huán)境教育部重點(diǎn)實(shí)驗(yàn)室, 重慶 400045； 2.中國(guó)交通建設(shè)股份有限公司， 北京 100088； 3.廣東海洋大學(xué) 海洋與氣象學(xué)院， 廣東 湛江 524088； 4.重慶工業(yè)職業(yè)技術(shù)學(xué)院， 重慶 401120)

1 研究背景

近5年來，我國(guó)300多個(gè)城市均發(fā)生了不同程度的內(nèi)澇災(zāi)害，而暴雨是形成洪澇災(zāi)害的主要原因[1]。 研究反映降雨規(guī)律的暴雨強(qiáng)度公式和參數(shù)對(duì)于城鎮(zhèn)防洪、海綿城市建設(shè)[2]有著重要的意義。目前我國(guó)暴雨強(qiáng)度公式的推求基本按照《室外排水設(shè)計(jì)規(guī)范》(GB50014-2006)推薦的基本步驟和方法進(jìn)行，另外還有不少文獻(xiàn)也對(duì)暴雨強(qiáng)度公式編制中的一些技術(shù)問題進(jìn)行了研究[2]，使得暴雨強(qiáng)度公式的編制逐步規(guī)范化。但目前對(duì)編制過程中基礎(chǔ)數(shù)據(jù)異常值的剔除、Cv、Cs值對(duì)擬合曲線的影響、有約束的規(guī)劃求解法在暴雨強(qiáng)度公式推求的應(yīng)用等問題的研究較少。

現(xiàn)有的皮爾遜-III分布適線研究[6]和公式推求方法[8-11]很少涉及異常值的剔除，隨著數(shù)據(jù)的更新、資料年限的增加，降雨資料中出現(xiàn)異常值的概率也隨之增加，這直接影響暴雨強(qiáng)度公式的準(zhǔn)確度。Cv、Cs值的大小不僅反映實(shí)測(cè)數(shù)據(jù)的分布規(guī)律，還影響暴雨強(qiáng)度-降雨歷時(shí)-重現(xiàn)期(q-t-P)的關(guān)系，只有在明確把握Cv、Cs的意義和影響的基礎(chǔ)上，才能根據(jù)樣本分布情況正確估計(jì)Cv、Cs的大小關(guān)系，得到準(zhǔn)確可靠的(q-t-P)表。此外，皮爾遜-III型曲線適線過程中推薦采用的目估適線法具有主觀性，工作量大，精度不夠高，通過約束條件采用規(guī)劃求解不僅能解決研究者發(fā)現(xiàn)參數(shù)出現(xiàn)負(fù)值的問題[12]，還能提高擬合速度。因此本文以重慶市酉陽(yáng)縣1993-2013年共21 a的年最大降雨資料為例，研究了異常值的剔除方法、Cv和Cs值對(duì)擬合曲線的影響規(guī)律以及目標(biāo)規(guī)劃求解的應(yīng)用，實(shí)現(xiàn)對(duì)公式的快速精確的推求。

2 降雨數(shù)據(jù)資料

本文收集了重慶市酉陽(yáng)縣1993-2013年共計(jì)21 a的逐分鐘降雨資料，按年最大值法取樣后將不同年份5～120 min共9個(gè)降雨歷時(shí)下的最大雨量統(tǒng)計(jì)如下(表1)。

表1 降雨數(shù)據(jù)資料

3 暴雨強(qiáng)度公式擬合與分析

3.1 異常值的剔除

我國(guó)制定的城市暴雨強(qiáng)度公式編制和設(shè)計(jì)暴雨雨型確定技術(shù)導(dǎo)則(以下簡(jiǎn)稱《導(dǎo)則》)[13]中提到了數(shù)據(jù)篩選，但是并未提出數(shù)據(jù)篩選的方法。異常值可能是由于雪融水、颶風(fēng)、人為等非正常因素造成的，將它們納入計(jì)算范圍會(huì)導(dǎo)致擬合曲線與普遍降雨規(guī)律偏離較遠(yuǎn)，因此需要根據(jù)氣象歷史資料剔除高異常值和低異常值。美國(guó)關(guān)于洪水控制的技術(shù)導(dǎo)則[14](以下簡(jiǎn)稱Bulletin 17B)中提出，選取樣本時(shí)會(huì)出現(xiàn)遇到特大值的情況，需要進(jìn)行專門分析確認(rèn)。一般特大值是指比相應(yīng)歷史資料序列的平均值(計(jì)入特大值后的均值)大2倍以上的稀遇暴雨值，合理取舍特大值有利于提升理論頻率曲線擬合的準(zhǔn)確度。通過分析研究，當(dāng)數(shù)據(jù)服從對(duì)數(shù)皮爾遜分布且偏態(tài)值在-3～3之間時(shí)，采用如下的方法剔除異常值較為合理：

(1)

(2)

小于qL的所有樣本數(shù)據(jù)應(yīng)被視為異常值而被剔除，而大于qH的數(shù)據(jù)只能在有氣象資料支撐的基礎(chǔ)上剔除最大值。當(dāng)樣本偏態(tài)系數(shù)大于0.4時(shí)，Bulletin 17B建議先剔除高異常值，偏態(tài)系數(shù)小于0.4時(shí)，先剔除低異常值，當(dāng)偏態(tài)系數(shù)介于-0.4～0.4時(shí)，對(duì)高、低異常值的剔除順序不作要求。

鑒于本文中酉陽(yáng)的樣本量?jī)H有21 a，異常值出現(xiàn)的概率較低，同時(shí)修訂公式時(shí)樣本量不應(yīng)小于20 a，對(duì)于每一個(gè)降雨歷時(shí)在本次擬合中僅需檢驗(yàn)是否存在一個(gè)異常值。同時(shí)，酉陽(yáng)城市防澇僅需達(dá)到20年一遇的等級(jí)，故采用0.1顯著水平下的KN值。在擬合酉陽(yáng)暴雨強(qiáng)度公式時(shí)，按照上述方法首先分析數(shù)據(jù)是否存在異常值，然后根據(jù)計(jì)算分別剔除了原始資料中45、60、90、120 min的最大暴雨強(qiáng)度值。

3.2 Cv、Cs值對(duì)擬合曲線的影響規(guī)律

皮爾遜-III型分布的概率密度函數(shù)為：

(3)

由于累積分布函數(shù)表示的是小于等于某一暴雨強(qiáng)度q的概率，而理論頻率是指大于等于某一暴雨強(qiáng)度q的概率,故按下列積分式計(jì)算理論頻率：

Ptheoretical=P(q≥qtheoetical)=

(4)

對(duì)于公式(4)的求解方法可以參考羅雅文[17]的方法，令t=β(q-b),公式(4)變?yōu)椋?/p>

Ptheoretical=P(t>ttheoretical)

(5)

通過以上分析，對(duì)Cv、Cs賦不同的數(shù)值計(jì)算皮爾遜-III型概率密度曲線及其積分后的累積概率曲線，曲線作圖詳見圖1所示。

由圖1(a)可知，當(dāng)Cs不變時(shí)，Cv越大，概率密度曲線矮且寬，累積分布曲線越來越緩，由此推導(dǎo)Cv越大，皮爾遜-III型曲線可擬合的范圍也越大，此時(shí)曲線能擬合出離散程度高的數(shù)據(jù)。當(dāng)實(shí)際降雨數(shù)據(jù)計(jì)算的概率集中在均值附近時(shí)，在進(jìn)行適線時(shí)可適當(dāng)減小Cv，而當(dāng)數(shù)據(jù)離散程度較大時(shí)，可以適當(dāng)增大Cv。由圖1(b)可知，當(dāng)Cv一定時(shí)，Cs越大，概率密度曲線越高且窄，左偏程度越高,累積分布曲線中部越向左偏，上段越陡，下段越平緩,由此可推導(dǎo)，當(dāng)數(shù)據(jù)正偏態(tài)程度越大，即小于平均值的數(shù)據(jù)占總數(shù)據(jù)的比例越大時(shí)，應(yīng)該用更大的Cs/Cv曲線去擬合數(shù)據(jù)。

綜上所述，在皮爾遜-III型曲線適線過程中，可通過計(jì)算樣本均值，分析數(shù)據(jù)的離散與偏態(tài)程度，結(jié)合圖1中數(shù)據(jù)離散程度與Cv、Cs變化的相關(guān)關(guān)系，來調(diào)整Cv/Cs比值，提高適線擬合速度。

3.3 基于規(guī)劃求解的分布曲線擬合及適線

為得到指定重現(xiàn)期及指定降雨歷時(shí)下的暴雨強(qiáng)度，需按照理論頻率分布擬合曲線，將分散的原始數(shù)據(jù)擬合得出內(nèi)在趨勢(shì)，為暴雨強(qiáng)度公式的外延計(jì)算創(chuàng)造條件。理論頻率分布曲線包括：皮爾遜-III型分布曲線(P-Ⅲ型曲線)、耿貝爾分布曲線(Gumbel曲線)和指數(shù)分布曲線。由于指數(shù)分布曲線適用于非年最大值法取樣的數(shù)據(jù)[18]，耿貝爾曲線是皮爾遜-III型曲線Cs=1.140的一個(gè)特例，其彈性不足，對(duì)樣本適應(yīng)性較差[19]，因此本研究采用皮爾遜-III型曲線來擬合暴雨強(qiáng)度-降雨歷時(shí)-重現(xiàn)期(q-t-P)的關(guān)系。

將通過規(guī)劃求解優(yōu)化擬合的皮爾遜-III型分布曲線與經(jīng)驗(yàn)頻率點(diǎn)共同體現(xiàn)在海森頻率格紙[20]上，得到如圖2所示的暴雨強(qiáng)度-降雨歷時(shí)-重現(xiàn)期(q-t-P)關(guān)系圖。

根據(jù)以上方法按照皮爾遜-III型分布擬合曲線，擬合的絕對(duì)均方誤差為0.0557 mm/min，相對(duì)均方誤差為6.27%。同時(shí)對(duì)未剔除異常值時(shí)誤差進(jìn)行了對(duì)比分析，發(fā)現(xiàn)未剔除異常值時(shí)理論頻率曲線的絕對(duì)均方誤差為0.0704 mm/min，相對(duì)均方誤差為7.03%；剔除異常值后絕對(duì)均方誤差值減小0.0147 mm/min，相對(duì)均方誤差值減小0.76%。

3.4 暴雨強(qiáng)度公式推求

通過計(jì)算qtheoretical或者由圖2查詢得到一組包含5個(gè)重現(xiàn)期(酉陽(yáng)縣管網(wǎng)按20年一遇防澇等級(jí)設(shè)計(jì)，故只需20、10、5、3、2年的重現(xiàn)期)、9個(gè)降雨歷時(shí)(5、10、15、20、30、45、60、90、120 min)的數(shù)據(jù)，分別采用最小二乘法和高斯牛頓法推導(dǎo)暴雨強(qiáng)度公式。

依據(jù)《室外排水設(shè)計(jì)規(guī)范》(GB50014-2006)[21]，暴雨強(qiáng)度公式定義為：

(6)

式中：q為暴雨強(qiáng)度，L/(s·hm2)；P為重現(xiàn)期，a；t為降雨歷時(shí)，min；A1、b、C、n為與地方暴雨特性有關(guān)且需求解的參數(shù)，其中A1為雨力參數(shù)，即重現(xiàn)期為1 a時(shí)的1 min設(shè)計(jì)雨量，mm ；C為雨力變動(dòng)參數(shù)；b為降雨歷時(shí)修正參數(shù)，即對(duì)暴雨強(qiáng)度公式兩側(cè)求對(duì)數(shù)后能使曲線轉(zhuǎn)為直線所加的時(shí)間參數(shù)，min；n為暴雨衰減指數(shù)，與重現(xiàn)期有關(guān)。

3.4.1 最小二乘法推求結(jié)果 為求得暴雨強(qiáng)度公式，對(duì)于公式(6)兩側(cè)取對(duì)數(shù)得：

lnq=ln 167A1+ln(1+ClgP)-nln(t+b)

(7)

令y=lnq，b0=ln 167A1，x1=ln(1+ClgP)，b2=-n,x2=ln (t+b)，公式(7)變?yōu)椋?/p>

y=β0+x1+β2x2

(8)

(9)

(10)

應(yīng)用最小二乘法求解暴雨強(qiáng)度公式曲線與皮爾遜-III頻率分布設(shè)計(jì)值的比較如圖2所示。最小二乘法求解暴雨強(qiáng)度公式的絕對(duì)均方誤差為0.0654 mm/min，相對(duì)均方誤差為3.96%。最小二乘法是以殘差平方和最小為目標(biāo)，對(duì)估量值求偏導(dǎo)從而求解參數(shù)，這種方法對(duì)異常值特別敏感，在計(jì)算中收斂速度不是很高，此外對(duì)于非線性函數(shù)需要進(jìn)行線性轉(zhuǎn)化，故最小二乘法在暴雨強(qiáng)度公式的擬合上具有一定局限性。

圖1 不同Cv、Cs對(duì)概率密度函數(shù)和累積分布函數(shù)的影響

圖2 暴雨強(qiáng)度-降雨歷時(shí)-重現(xiàn)期(q-t-P)關(guān)系的經(jīng)驗(yàn)頻率點(diǎn)和理論頻率曲線

(11)

高斯牛頓法采用Matlab程序求解，得暴雨強(qiáng)度公式為：

(12)

應(yīng)用高斯牛頓法求解暴雨強(qiáng)度公式曲線與皮爾遜-III頻率分布設(shè)計(jì)值的比較如圖4所示。高斯牛頓法求解暴雨強(qiáng)度公式的絕對(duì)均方誤差為0.0502 mm/min，相對(duì)均方誤差為5.38%。

圖3最小二乘法求解暴雨強(qiáng)度公式推求值 圖4高斯牛頓法求解暴雨強(qiáng)度公式推求值

與P-III擬合值的比較 與P-III擬合值的比較

4 結(jié) 論

針對(duì)暴雨強(qiáng)度公式推求中存在的問題，以實(shí)際暴雨強(qiáng)度數(shù)據(jù)結(jié)合理論知識(shí)對(duì)推求過程的幾個(gè)問題進(jìn)行了研究，得到以下幾個(gè)結(jié)論：

(1)通過研究分析提出了異常值的剔除方法，結(jié)合分析Cv、Cs的意義總結(jié)其對(duì)適線的影響規(guī)律。當(dāng)數(shù)據(jù)服從對(duì)數(shù)皮爾遜分布且偏態(tài)值在-3～3之間時(shí)，采用對(duì)數(shù)形式的高低限值方法剔除異常值，合理取舍異常值有利于提升理論頻率曲線擬合的準(zhǔn)確度；當(dāng)數(shù)據(jù)集中分布在均值左右時(shí)適當(dāng)減小Cv，當(dāng)數(shù)據(jù)正偏程度較大時(shí)，適當(dāng)增大Cs/Cv，為方便準(zhǔn)確地?cái)M合出暴雨強(qiáng)度-降雨歷時(shí)-重現(xiàn)期(q-t-P)曲線提供了支撐。

(2)通過約束條件，設(shè)置目標(biāo)函數(shù)采用Excel規(guī)劃求解Cv、Cs，與目估適線法相比，明顯提高了擬合速度和精度。

(3)基于規(guī)劃求解的最小二乘法推導(dǎo)的暴雨強(qiáng)度公式絕對(duì)均方誤差為0.0654 mm/min，相對(duì)均方誤差為3.96%；高斯牛頓法推導(dǎo)的暴雨強(qiáng)度公式絕對(duì)均方誤差為0.0502 mm/min，相對(duì)均方誤差為5.38%，通過高斯牛頓法推求的酉陽(yáng)暴雨強(qiáng)度公式較合理。高斯牛頓法可采用最小二乘法的計(jì)算值作為初始值進(jìn)行迭代計(jì)算，相比將而言高斯牛頓法具有更快的收斂速度和更高的擬合精度。