本發(fā)明涉及水文站網(wǎng)技術(shù)，具體涉及一種基于克里金和信息熵理論相結(jié)合的水文站網(wǎng)優(yōu)化模型。

背景技術(shù)：

水文站是一個在河流或者流域上設(shè)立的主要用于觀測以及搜集河流湖泊和水庫等水體相關(guān)水文和氣象資料的基層水文機構(gòu)，在很大程度上滿足了水文預(yù)報、水文情報、水資源評價工作和水科學(xué)研究的基本需求。因此規(guī)劃合理的水文站網(wǎng)能夠充分反映水文時空變異特征，使之能收集準(zhǔn)確詳盡的水文信息，這顯然是有必要探究較為客觀地理論方法支撐水文站網(wǎng)的合理規(guī)劃。一個有效水文站網(wǎng)體系一方面能反映流域時空分布特性并能夠預(yù)測未知的點到面的相關(guān)數(shù)據(jù)的，另一方面在經(jīng)濟合理的前提下的站點組合能夠反映最優(yōu)的信息量。前人在站網(wǎng)規(guī)劃的研究方法上大體可以歸納為以下幾種：

1、數(shù)理統(tǒng)計方法，該方法應(yīng)用最早，但要求研究者對水資源系統(tǒng)結(jié)構(gòu)有充足的認(rèn)識；同時，受到數(shù)理統(tǒng)計的原理的限制，統(tǒng)計分析技術(shù)的選擇和樣本數(shù)量都將對數(shù)據(jù)分析的結(jié)論產(chǎn)生較大影響；并且該方法只能通過估計精度與樣本數(shù)量之間的關(guān)系來確定站點的數(shù)量，達不到對站點的空間布局優(yōu)化的目的。

2、克里格插值法，需要在圖上對水文變量估算誤差的改進作出主觀評價，若估計誤差的標(biāo)準(zhǔn)差大于給定的標(biāo)準(zhǔn)差的范圍，則空間中需要增設(shè)站點，反之就應(yīng)該減少站點，實際應(yīng)用中單一采用該方法主觀性比較大精度不夠高。

3、信息熵方法，在以往基于信息熵進行水文站網(wǎng)優(yōu)化研究中，常以信息傳遞量為依據(jù)，選擇某一指標(biāo)進行站網(wǎng)評價優(yōu)化，而單一信息熵指標(biāo)往往不能全面反映站網(wǎng)所包含的信息總量、信息冗余程度等決定站網(wǎng)效用的重要因素。站對間的傳遞熵的估計很大程度上受限于聯(lián)合概率密度函數(shù)的合理估計。

技術(shù)實現(xiàn)要素：

發(fā)明目的：本發(fā)明的目的在于解決現(xiàn)有技術(shù)中存在的不足，提供一種基于克里金和信息熵理論相結(jié)合的水文站網(wǎng)優(yōu)化模型。

技術(shù)方案：本發(fā)明一種基于克里金法和信息熵理論耦合的水文站網(wǎng)優(yōu)化模型，其特征在于：依次包括如下步驟：

(1)變異函數(shù)的優(yōu)選：基于最小累計殘差的準(zhǔn)則(smallestresidualsumofsquares,sserr)對已有的站網(wǎng)空間分布關(guān)系進行變異函數(shù)的模擬；從潛在的4個變異函數(shù)模型即指數(shù)(exponential)、球形(spherical)，高斯(gaussian)和matern變異函數(shù)中依據(jù)殘差最小原則選擇最優(yōu)的變異函數(shù)模型，其中克里金方差定義為：

其中，h0i為已知站點和未知站點間的距離，γ(h0i)是距離h0i下的變異值，μx為拉格朗日乘子，n為空間分布劃分的樣點個數(shù)，wi為權(quán)重值，σx(y0)為克里金標(biāo)準(zhǔn)差；

(2)決定增設(shè)站點區(qū)域和具體位置：首先，得到克里金標(biāo)準(zhǔn)差分布圖(ksemap)，依據(jù)圖中標(biāo)準(zhǔn)差最大的原則準(zhǔn)確定位出需要增加的站點；

(3)增設(shè)站點雨量的插值：得到增加站點的日雨量值，并依據(jù)交叉檢驗數(shù)據(jù)統(tǒng)計量評價變異函數(shù)模型的效果；基于均方根誤差(rootmeansquarepredictionerror(rmse)),均標(biāo)準(zhǔn)差(meansquarednormalizederror,mnse)和相關(guān)系數(shù)(correlationcoefficient)：

其中，為位置yi處的雨量估計值而x(yi)為觀測值，σ²(yi)為相關(guān)的克里金方差；

(4)基于ni-kse準(zhǔn)則優(yōu)選出最優(yōu)的站點組合：首先，最小克里金標(biāo)準(zhǔn)差目標(biāo)函數(shù)即可由第一步中的克里金方差開方得出，而最大信息量ni的目標(biāo)函數(shù)通過以下公式得出：

其中，ni(x1，x1，...，xk)代表從新增站點后得到的新的站網(wǎng)中選出的k個站點組成的子站網(wǎng)的信息總量，x1，x2，...，xk則是各個雨量站，xj是沒有被選入備選集合的站點，μ1是信息權(quán)重指標(biāo)，h(x1，x2，...，xk)為聯(lián)合信息熵值，t(xi，xj)為信息傳遞量，c(x1，x1，...，xk)為總的信息關(guān)聯(lián)量(totalcorrelation)；

接著通過信息量目標(biāo)函數(shù)權(quán)重值μ1確定ni-kse準(zhǔn)則的計算值(trni-kse)：

mintrni-kse＝0.6rankdescend(ni)+0.4rankascend(kse)

其中，rankdescend(ni)表示信息總量按降序下得到序數(shù)，rank(kseascend)為kse在升序下得到的序數(shù)，所以ni-kse準(zhǔn)則就是綜合兩個目標(biāo)函數(shù)得到單目標(biāo)函數(shù)下最小值時，該站點組合為最優(yōu)的站點組合；

(5)得出和原始站網(wǎng)同等數(shù)目下的最優(yōu)站點組合后，將最優(yōu)的組合和原始站點組合的雨量估計值和觀測值的誤差進行對比分析，以此檢驗克里金和信息熵理論結(jié)合的站網(wǎng)優(yōu)化模型的效果。

進一步的，所述步驟(1)中變異函數(shù)的優(yōu)選步驟具體為：

(1.1)將每個站點的觀測數(shù)據(jù)作為隨機變量的觀測值xi(i＝1,2,…,n)，求得邊緣熵值h(xi)；模擬處樣本變異函數(shù)模型：對原站網(wǎng)每天的數(shù)據(jù)空間分布進行經(jīng)驗?zāi)M，并畫出適合的樣本分布曲線；

(1.2)基于原始數(shù)據(jù)得出四種變異函數(shù)模型下的理論分布曲線；

(1.3)計算四種理論模型和樣本分布模型的最小累計殘差值，優(yōu)選出其中最小的殘差值對應(yīng)的理論變異函數(shù)模型。

進一步的，所述步驟(4)中基于ni-kse準(zhǔn)則的優(yōu)選步驟為：

(4.1)初始化優(yōu)選站點集合s0并設(shè)置為空集，潛在的站網(wǎng)集合為f,并將f包含了原始站點和新增站點；

(4.2)由函數(shù)max{h(xi)}，得出最大邊緣熵條件下對應(yīng)的站點i，并將此站點放到集合s中設(shè)定為s1，站點集f刪除站點i，更新站點集合s和f；

(4.3)依據(jù)最大信息量函數(shù)(maxni(x1，x2，...，xk))從f中得出新的優(yōu)選站點，實時更新站點集合s和f，并實時計算出各個站點組合下的kse和ni值；

(4.4)依次按照以上步驟可以將原始集合f的的站點按照重要程度進行排序；

(4.5)將排序后的站點組合按照ni-kse準(zhǔn)則下的評價指標(biāo)mintrni-kse重新排序得出最優(yōu)的站點組合。

有益效果：本發(fā)明將最大化信息總量和最小化估計誤差值作為目標(biāo)函數(shù)，并提出了ni-kse準(zhǔn)則下的綜合評價指標(biāo)值(trni-kse)得出最優(yōu)的站點組合。其中信息總量這一評價因素兼顧到了最大信息聯(lián)合熵、最小信息冗余量和最大信息傳遞能力三個子目標(biāo)?？死锝饦?biāo)準(zhǔn)差(kse)只需要得出合理的變異函數(shù)模型就可以得到對應(yīng)的數(shù)值。與現(xiàn)有技術(shù)下相比，本發(fā)明具有以下優(yōu)點：

(1)克服了傳統(tǒng)的單一方法的局限性。一方面，單一的克里金插值法僅僅從誤差最小的方面評選出增設(shè)和精簡站點的區(qū)域，單一采用該方法主觀性比較大精度不夠高，另外前人采用的克里金法較多采用特定的變異函數(shù)模型，本專利考慮了4個潛在的變異函數(shù)模型顯得較為客觀全面；另一方面，單一的信息熵方法僅僅從信息聯(lián)合熵、最小化互信息的角度上優(yōu)選最優(yōu)的站點組合，沒有考慮到站點間的空間變異性，也無法增設(shè)站點只能在原有站點基礎(chǔ)上進行精簡，站網(wǎng)的預(yù)測誤差和估計精度不能得到保證。而本專利系統(tǒng)性的考慮到信息量最大化和估計誤差最小化兩大目標(biāo)，綜合并最大化地發(fā)揮了了克里金和信息熵理論的優(yōu)勢，是優(yōu)化模型更加的系統(tǒng)全面。

(2)綜合性地提出三個站網(wǎng)評價指標(biāo)maxni，minkse，mintrni-kse，對站網(wǎng)的站點合理分布提供較為客觀依據(jù)。同時采用了交叉檢驗統(tǒng)計數(shù)據(jù)評價變異函數(shù)模型的效果。

綜上所述，本發(fā)明將克里金和信息熵理論結(jié)合對站網(wǎng)空間分布特性和系統(tǒng)的信息總量的不確定性進行了定量地描述，既可實現(xiàn)對站網(wǎng)信息的定量分析，又能合理布置站點的空間分布，具有合理性和有效性。

附圖說明

圖1為本發(fā)明的整理流程圖；

圖2為實施例中流域站點示意圖；

圖3為實施例中變異函數(shù)模擬圖；

圖4為實施例中克里金誤差分布圖；

圖5為實施例中最優(yōu)的站點分布圖；

其中，圖3(a)為算例1的變異函數(shù)模擬圖，圖3(b)為算例2的變異函數(shù)模擬圖；圖3(c)為算例3的變異函數(shù)模擬圖；圖3(d)為算例4的變異函數(shù)模擬圖。

具體實施方式

下面對本發(fā)明技術(shù)方案進行詳細說明，但是本發(fā)明的保護范圍不局限于所述實施例。

為便于理解本發(fā)明，做以下說明：

四種變異函數(shù)模型歸納如下:

●指數(shù)函數(shù)模型：

●球形函數(shù)模型：

●高斯函數(shù)模型：

●matern函數(shù)模型：

γ(h)＝(s-n)[1-(1-k)exp(1r|h|)-kexp(-rh²)]+n，

其中s，n，r，k是變異函數(shù)的參數(shù)。

假定(x1，x2，...，xn)代表離散化的隨機變量組合(n個站點組成的原始站網(wǎng))，其聯(lián)合概率密度函數(shù)為p(x1，x2，...，xn)，xi代表站點i處的水文時間序列值，邊緣概率密度函數(shù)分別為p(x1)，p(x2)，...，p(xn)，那么在水文站網(wǎng)優(yōu)化過程中，邊緣熵可以定義為

多變量的聯(lián)合信息熵h(x1，x2，...，xn)定義為：

傳遞量指標(biāo)(或者互信息)可以表示為：

為了刻畫系統(tǒng)的信息冗余量，引入了信息冗余指標(biāo)(totalcorrelation,tc),其可以定義為：

其中，h(xi)為站點i的雨量時間序列xi的邊緣熵。

如圖1所示，本發(fā)明的基于克里金和信息熵理論結(jié)合的水文站網(wǎng)優(yōu)化模型，具體包含以下步驟：

1.變異函數(shù)的優(yōu)選：基于最小累計殘差的準(zhǔn)則(smallestresidualsumofsquares,sserr)對已有的站網(wǎng)空間分布關(guān)系進行變異函數(shù)的模擬。從潛在的4個變異函數(shù)模型即指數(shù)(exponential)，球形(spherical)，高斯(gaussian)和matern變異函數(shù)中依據(jù)殘差最小原則選擇最優(yōu)的變異函數(shù)模型。具體的變異函數(shù)的優(yōu)選步驟如下：

將每個站點的觀測數(shù)據(jù)作為隨機變量的觀測值xi(i＝1,2,…,n)，求得邊緣熵值h(xi)。

第1步，模擬處樣本變異函數(shù)模型：對原站網(wǎng)每天的數(shù)據(jù)空間分布進行經(jīng)驗?zāi)M，并畫出適合的樣本分布曲線；

第2步，基于原始數(shù)據(jù)得出四種變異函數(shù)模型下的理論分布曲線；

第3步，計算四種理論模型和樣本分布模型的最小累計殘差值，優(yōu)選出其中最小的殘差值對應(yīng)的理論變異函數(shù)模型。

2.決定增設(shè)站點區(qū)域和具體位置：首先，在得出適合的理論變異函數(shù)后得到克里金標(biāo)準(zhǔn)差分布圖(ksemap)，依據(jù)圖中標(biāo)準(zhǔn)差最大的原則準(zhǔn)確定位出需要增加的站點。

盡管，這一步簡單但是確實較為關(guān)鍵的一步，因為這一步將直接影響到接下來的幾步研究對象的具體位置。顯而易見地是，kse分布圖很好地詮釋了任意站點的空間變異性，一般來說站點密度較大的區(qū)域，kse值就會相對較小，而站點稀疏的地區(qū)誤差值就會較大需要增加站點滿足估計精度的要求。

3.增設(shè)站點雨量的插值：得到增加站點的日雨量值，并依據(jù)交叉檢驗數(shù)據(jù)統(tǒng)計量評價變異函數(shù)模型的效果?；诰礁`差(rootmeansquarepredictionerror(rmse)),均標(biāo)準(zhǔn)差(meansquarednormalizederror,mnse)和相關(guān)系數(shù)(correlationcoefficient)：

其中，為位置yi處的雨量估計值而x(yi)為觀測值，σ²(yi)為相關(guān)的克里金方差。

4.基于ni-kse準(zhǔn)則優(yōu)選出最優(yōu)的站點組合：首先，最小克里金標(biāo)準(zhǔn)差目標(biāo)函數(shù)即可由第一步中的克里金方差開方得出，而最大信息量(ni)的目標(biāo)函數(shù)需要以下公式得出：

其中，ni(x1，x1，...，xk)代表從新增站點后得到的新的站網(wǎng)中選出的k個站點組成的子站網(wǎng)的信息總量，x1，x2，...，xk則是各個雨量站下的水文時間序列，xj是沒有被選入備選集合的站點j處的水文時間序列，μ1是信息權(quán)重指標(biāo)，h(x1，x2，...，xk)為聯(lián)合信息熵值，t(xi，xj)為信息傳遞量，tc(x1，x1，...，xk)為總的信息冗余量(totalcorrelation)。

一旦確定了信息量目標(biāo)函數(shù)權(quán)重值μ1，就能得到ni-kse準(zhǔn)則的計算值(trni-kse),

mintrni-kse＝0.6rankdescend(ni)+0.4rankascend(kse)，

其中，rankdescend(ni)表示信息總量按降序下得到序數(shù)，rank(kseascend)為kse在升序下得到的序數(shù)。所以ni-kse準(zhǔn)則就是綜合兩個目標(biāo)函數(shù)得到單目標(biāo)函數(shù)下最小值時，該站點組合為最優(yōu)的站點組合。具體的ni-kse的優(yōu)選步驟如下:

(1)初始化優(yōu)選站點集合s0并設(shè)置為空集，潛在的站網(wǎng)集合為f,并將f包含了原始站點和新增站點。

(2)由函數(shù)max{h(xi)}，得出最大邊緣熵條件下對應(yīng)的站點i，并將此站點放到集合s中設(shè)定為s1。站點集f刪除站點i。更新站點集合s和f。

(3)依據(jù)最大信息量函數(shù)(maxni(x1，x2，...，xk))從f中得出新的優(yōu)選站點。實時更新站點集合s和f。并實時計算出各個站點組合下的kse和ni值。

(4)依次按照以上步驟可以將原始集合f的的站點按照重要程度進行排序。

(5)將排序后的站點組合按照ni-kse準(zhǔn)則下的評價指標(biāo)mintrni-kse重新排序得出最優(yōu)的站點組合。

5、雨量估計值和觀測值的誤差進行對比：得出和原始站網(wǎng)同等數(shù)目下的最優(yōu)站點組合后，將最優(yōu)的組合和原始站點組合的雨量估計值和觀測值的誤差進行對比分析，以此檢驗克里金和信息熵理論結(jié)合的站網(wǎng)優(yōu)化模型的效果。

實施例1：本實施以上海市水文站網(wǎng)優(yōu)化作為實際應(yīng)用

以上海16個水文站組成的站網(wǎng)為例，以2012年1月1日到12月31日的日雨量序列為樣本，用基于克里金和信息熵理論結(jié)合的水文站網(wǎng)優(yōu)化模型對該站網(wǎng)進行評價和優(yōu)化。

(1)流域概況

本文的數(shù)據(jù)資料來源于上海市，2012年1月1日-2012年12月31日逐日雨量數(shù)據(jù)。上海市全市面積6340.5平方公里，位于長江的入?？诤吞饔虻臇|側(cè)，位于長三角經(jīng)濟發(fā)展區(qū)域的核心。河流湖泊總面積達500平方公里，其中河流覆蓋率達到9％-10％(見圖2、表1)。圖2中有16個已有的站點，經(jīng)過克里金插值分析新增了5個站點。站點基本信息量和坐標(biāo)設(shè)置見表1。

表1上海市雨量站一覽表

(2)模型運行

首先對上海水文站網(wǎng)內(nèi)的16個站點進行編號(1～16)，并選取了四個算例(scenario)用于展示。(1)scenario1:meandailyrainfalldata；(2)scenario2:rainfalldataon2012/1/21；(3)scenario3:rainfalldataon2012/4/23；(4)scenario4:rainfalldataon2012/11/22。變異函數(shù)模擬結(jié)果見圖3。kse圖見圖4。

表2不同權(quán)重指標(biāo)下的站點排序結(jié)果

表3基于ni-kse準(zhǔn)則下的站點優(yōu)化結(jié)果

(3)站網(wǎng)評價

圖3中可得不同的空間分布結(jié)構(gòu)適合不同的變異函數(shù)模型，為此選擇了4個潛在的理論模型非常有必要。

圖4中可以得出需要在17,18,19,20,21點(表1中)的坐標(biāo)位置增加5個站點，因為這些未知站點的標(biāo)準(zhǔn)差值最大。接著需要利用克里金插值法得出這些站點的時間序列值進而進行接下來的ni-kse準(zhǔn)則分析得到最優(yōu)站點的組合。

(4)站網(wǎng)優(yōu)化

依據(jù)表2可得在權(quán)重指標(biāo)μ1為0.8和0.9時，排序結(jié)果趨于穩(wěn)定，可知權(quán)重指標(biāo)為0.8或0.9為宜。為了強調(diào)信息傳遞能力的重要性，選擇了0.9作為后續(xù)的權(quán)重系數(shù)。表3是基于表2在站的個數(shù)為9個基礎(chǔ)上重新進行ni-kse準(zhǔn)則的進一步優(yōu)化的。利用最大化信息量和最小化克里金標(biāo)準(zhǔn)差之后，能夠得到更為合理客觀的結(jié)果，在站點個數(shù)為10時，最優(yōu)的站點組合是s10＝(10,19,20,17,7,21,3,9,6，11)。依照相同的原理和方法得到了站點數(shù)目為16的情況下最優(yōu)的站點組合為s16＝(1,3,6,7,9,10,11,13,17,19,20,21,5,16,2,15)。其分布圖為圖5所示，由圖中可知經(jīng)過ni-kse準(zhǔn)則優(yōu)選后需要新增4個站點(17,19,20,21)，精簡去掉4個舊站點(4,8,12,14)。

(5)模型效果測試

為了能夠通過類似后驗的方法得到模型算出的最優(yōu)組合確實能夠達到比原始的站網(wǎng)達到更好的效果。由表4，可以看出最優(yōu)的站點組合只得出了0.96％的誤差百分率而原始的站網(wǎng)有2.24％的誤差。在信息指標(biāo)方面，新的站網(wǎng)組合也比舊站網(wǎng)提供出較大的聯(lián)合熵和信息傳遞能力，較小的信息冗余量值。為此，基于克里金和信息熵的水文模型優(yōu)選出的站點組合確實能提供較為出色的估計能力和信息量。

表4站點組合效果對比

^a聯(lián)合熵；^b互信息；^c多維冗余信息量(totalcorrelation,tc)；^d站網(wǎng)總的信息量。通過上述實施例可以看出，本發(fā)明在考慮到信息量最大化的同時也能考慮到站點間的空間分布特性，并最終能夠提供比較好的模擬效果，既可實現(xiàn)對站網(wǎng)信息的定量分析，又能合理布置站點的空間分布，具有合理性和有效性，并最終為水文站網(wǎng)合理規(guī)劃和站點增刪提供技術(shù)和決策支持。