水文頻率計(jì)算中窗寬和核函數(shù)對(duì)密度函數(shù)估計(jì)影響分析

董 潔，刁艷芳，譚秀翠，徐 妍

(山東農(nóng)業(yè)大學(xué)水利土木工程學(xué)院，山東 泰安 271018)

水利建設(shè)規(guī)劃設(shè)計(jì)時(shí)，要確定一個(gè)合理的設(shè)計(jì)標(biāo)準(zhǔn)或工程規(guī)模。一般采用水文頻率計(jì)算方法計(jì)算設(shè)計(jì)值，目前主要有參數(shù)和非參數(shù)數(shù)理統(tǒng)計(jì)方法[1]。水文頻率計(jì)算需要對(duì)水文變量的總體密度進(jìn)行估計(jì)。如果密度函數(shù)結(jié)構(gòu)已知而只有其中某些參數(shù)未知，此時(shí)的密度估計(jì)就是傳統(tǒng)的參數(shù)估計(jì)問題。如果密度函數(shù)未知(或最多只知道連續(xù)、可微等條件)，僅從即有的樣本出發(fā)得出密度函數(shù)的表達(dá)式，這就是非參數(shù)密度估計(jì)。非參數(shù)密度估計(jì)始于直方圖法，后來發(fā)展為最近鄰法、核估計(jì)法等，其中理論發(fā)展最完善的是密度的核估計(jì)法[2]。

1 核估計(jì)

1.1 概 述

設(shè){x1,…,xn}為離散的隨機(jī)樣本，單變量核密度估計(jì)為：

(1)

核估計(jì)既與樣本有關(guān)，又與核及窗寬的選取有關(guān)。在給定樣本以后，一個(gè)核估計(jì)的好壞，取決于核及窗寬的選取是否得當(dāng)。窗寬和核的選擇直接影響密度函數(shù)的估計(jì)精度[3]。

(2)

1.2 最優(yōu)窗寬的確定

1.2.1 理論界定

(3)

從而可以得到這種意義下的最優(yōu)窗寬表達(dá)式：

(4)

由此可以看出：

(1)最優(yōu)窗寬應(yīng)隨樣本的增大而不斷減小且速度為o(n)。

(2)f″(x)反映密度函數(shù)的震動(dòng)速率，劇烈震動(dòng)的密度函數(shù)應(yīng)對(duì)應(yīng)較小的最優(yōu)窗寬。

(3)表達(dá)式中含有未知量f″(x)，因此無法得到具體的窗寬數(shù)值。

1.2.2 窗寬的選擇

窗寬的選擇一般分為固定和變窗寬[4,5]。

(1)固定窗寬。就是在每一個(gè)擬合點(diǎn)取等窗寬，缺點(diǎn)是所估計(jì)量不能充分利用變量X的設(shè)計(jì)密度所提供的信息，且對(duì)復(fù)雜曲線的擬合效果欠佳。

(2)變窗寬。有局部變窗寬和全局變窗寬2類。局部變窗寬h(x0)隨位置x0的變化而變化，全局變窗寬h(xj)隨數(shù)據(jù)點(diǎn)xj的變化而變化。變窗寬的引入可以反映不同點(diǎn)的光滑程度，降低擬合曲線在峰頂區(qū)域的偏差以及尾部區(qū)域的方差，提高擬合曲線的靈活性，適用于空間非齊次曲線的擬合。例如交叉證實(shí)法Cross-Validation(CV)。

前面提到的窗寬選擇需要對(duì)擬合的密度函數(shù)有一定的假設(shè)，而CV法是一種數(shù)據(jù)本源(data based)方法，不需要對(duì)擬合密度函數(shù)假設(shè)，而是從現(xiàn)有的數(shù)據(jù)直接得到合理的窗寬。由樣本{X1,X2,…,Xn}作缺值估計(jì)：

(5)

但是，當(dāng)核不是密度函數(shù)時(shí)，估計(jì)量已經(jīng)不是密度函數(shù)，進(jìn)而不能用極大似然的思想求得，我們可以由以下最小平方差的思想LSCV(Least Square CV)求之，算出積分平方差I(lǐng)SE(Integrated Square Error):

(6)

好的密度估計(jì)函數(shù)應(yīng)對(duì)應(yīng)較小的ISE，或:

(7)

(8)

可以得到“最優(yōu)”窗寬。但是在實(shí)踐中常會(huì)出現(xiàn)不夠光滑的現(xiàn)象，而且這種窗寬的計(jì)算量太大，占用的時(shí)間太長(zhǎng)，因而，下面給出簡(jiǎn)便可行的方法。

1.2.3 確定最優(yōu)窗寬的具體方法

表1 擬合不同密度函數(shù)時(shí)窗寬的計(jì)算結(jié)果Tab.1 The window width calculation results of fittingdifferent density function

1.3 核的選擇

(1)一般核函數(shù)屬于對(duì)稱的密度函數(shù)族P，即q(·)滿足如下條件：

(9)

從減小積分均方誤差(L2)的角度來看，Silverman[6]指出P族中不同核對(duì)減小積分均方誤差沒有明顯差別，因此一般可根據(jù)其他需要(如計(jì)算方便)選定合適的核。后面的實(shí)例中就是考慮到計(jì)算方便以及水文的特點(diǎn)，我們選用了指數(shù)函數(shù)作為核。

(2)核函數(shù)為高階函數(shù)族Hs，即其中q(·)滿足如下條件：

(12)

引入這種函數(shù)的道理是基于以下命題。

命題：設(shè)核q(·)∈Hs具有s階導(dǎo)數(shù)，則積分均方誤差為：

(13)

從命題可以看出這種核的優(yōu)勢(shì)在于隨著階s的增大：

(14)

隨之減小，進(jìn)而積分均方差減小，不足是由它做成的核不是非負(fù)函數(shù)，進(jìn)而不是密度函數(shù)。因此，在水文計(jì)算中，我們通常不選此類核。

2 模擬研究

2.1 窗寬對(duì)模擬結(jié)果的影響

為了考證不同窗寬對(duì)密度擬合的影響，以下模擬以對(duì)數(shù)正態(tài)密度曲線為所要估計(jì)的曲線，以指數(shù)核(2)為核函數(shù)，分別取4種依次遞增窗寬(1，2，3，4)進(jìn)行了模擬。從圖1～圖4中可以看出，隨著窗寬的增大，模擬曲線也越光滑，所以窗寬的選擇對(duì)估計(jì)結(jié)果有較大的影響。

2.2 核對(duì)估計(jì)誤差的影響

前面已經(jīng)指出，不同核的估計(jì)量的影響在實(shí)踐中不大，為驗(yàn)證它，我們用計(jì)算機(jī)隨機(jī)生成500個(gè)服從對(duì)數(shù)正態(tài)分布的隨機(jī)數(shù)，分別以不同的核作為密度估計(jì)量，為了便于比較，取同樣的窗寬，并做其圖形，考察這些模擬圖形對(duì)理論對(duì)數(shù)正態(tài)密度函數(shù)的擬合情況。在模擬中分別取以下核。

(1)均勻核：

圖1 窗寬1擬合曲線Fig.1 Window width 1 fitting curve

圖2 窗寬2擬合曲線Fig.2 Window width 2 fitting curve

圖3 窗寬3擬合曲線Fig.3 Window width3 fitting curve

圖4 窗寬4擬合曲線Fig.4 Window width 4 fitting curve

(2)指數(shù)核:

(3)cauchy核：

(4)EV1核：

q(x)=exp [-x2-exp(-x2)]

圖5～圖8中光滑曲線是對(duì)數(shù)正態(tài)密度曲線，其余曲線是分別以上述核為函數(shù)的模擬曲線。

從模擬的結(jié)果可以看出，均勻核模擬曲線不夠光滑，指數(shù)核最好，其他核比較光滑，差別也不太多。

圖5 均勻核擬合曲線Fig.5 Uniform kernel fitting curve

圖6 指數(shù)核擬合曲線Fig.6 Exponential kernel fitting curve

圖7 cauchy核擬合曲線Fig.7 Cauchy kernel fitting curve

圖8 EV1核擬合曲線Fig.8 EV1 nuclear fitting curve

3 水文頻率計(jì)算模型實(shí)例

表2 五龍口水文站徑流量頻率計(jì)算Tab.2 runoff frequency calculation of Wulongkou station

經(jīng)過用參數(shù)適線法(LP)、線性距法(LM)和非參數(shù)核估計(jì)方法(NP)分析論證，頻率小于1%時(shí)，非參數(shù)核估計(jì)法計(jì)算的流量設(shè)計(jì)值要大于參數(shù)法，頻率大于1%時(shí)，非參數(shù)核估計(jì)法的流量設(shè)計(jì)值要小于參數(shù)法。3種方法結(jié)合使用，可以對(duì)比分析，以確定較合理的流量設(shè)計(jì)值。

4 結(jié) 語

(1)分析了核估計(jì)方法中的核和窗寬的選擇問題，并針對(duì)不同的核、窗寬對(duì)估計(jì)結(jié)果產(chǎn)生的影響進(jìn)行了模擬研究。結(jié)果表明，指數(shù)核的擬合較好；窗寬的選取是非常重要的，它決定著擬合的計(jì)精度。最優(yōu)窗寬應(yīng)當(dāng)使統(tǒng)計(jì)量的積分均方差(MISE)最小，并推出了最優(yōu)窗寬的近似計(jì)算公式。

(2)根據(jù)核函數(shù)和最優(yōu)窗寬的討論，用非參數(shù)核估計(jì)和參數(shù)適線法、線性距法對(duì)五龍口水文站的流量進(jìn)行了頻率計(jì)算，通過3種方法對(duì)比分析，可以確定合理的流量設(shè)計(jì)值，為解決水文頻率計(jì)算提供了有效方法。

□

[1] 黃振平，陳元芳.水文統(tǒng)計(jì)學(xué)[M].北京：中國(guó)水利水電出版社，2011.

[2] 陳希儒，柴根象.非參數(shù)統(tǒng)計(jì)教程[M].上海：華東師范大學(xué)出版社，1993.

[3] C G Lambert. Efficient on-line nonparametric kernel density estimation[J]. Algorithmica, 2004,25(1).

[4] Adamowski k. Regional analysis of analysis of annual maximum and partial duration flood data by nonparametric and l-moment methods[J].Journal of Hydrology，2000，(229)：219-231.

[5] Goel N K. A derived flood frequency distribution for correlated rainfall intensity and duration[J].Journal of Hydrology 2000,(228):56-67.

[6] Silverman B W. Choosing window width when estimating a density[J]. Biometrika，1978,(65):1-11.

[7] 閆寶偉，劉彥成. 基于兩變量核密度估計(jì)的枯水頻率分析[J]. 人民珠江,2017,38(4):15-20.

[8] 王雪妮，周 晶.一種新的洪水頻率分析方法研究[J].水利學(xué)報(bào)，2016,47(6):798-802.