內(nèi)蒙古河套灌區(qū)土壤水鹽空間變異性的ik法閾值選擇

內(nèi)蒙古河套灌區(qū)土壤水鹽空間變異性的ik法閾值選擇

指導思想的不規(guī)則性（縮寫ik）法是條庫爾提出的一種非參數(shù)估計方法。史密斯和霍夫最初用于土壤和生態(tài)科學。IK法與其他插值估計方法比較有明顯優(yōu)點:它不依賴空間現(xiàn)象的平穩(wěn)性;也不要求區(qū)域化變量服從某種分布假設,由于對隨機變量經(jīng)過特有的指示變換后,使隨機函數(shù)的貢獻都一樣,可自然消除特異值(Outlier)的影響;可對非取樣點處的不確定性進行條件概率分布函數(shù)(ConditionalCumulativeDistributionFunction簡寫CCDF)估計,尤其是當人們感興趣的不是某些點的具體估值,而是大于或小于某閾值(CutoffValue)或某范圍內(nèi)物質(zhì)品位在空間的分布概率時,可簡捷獲得在一定風險條件下未知量分布概率的估值,符合地學科學家用人工智能思維處理空間問題的思路。IK法自創(chuàng)立以來,在國外廣泛用于地質(zhì)、土壤、水文(Neuman)等領域,成為一種直接估計局部區(qū)域品位分布的數(shù)學地質(zhì)新方法。20世紀90年代開始在我國地學及土壤學科應用,也有一些報導:如地質(zhì)礦產(chǎn)方面侯景儒、土壤水文方面李保國等,水庫工程地質(zhì)方面吳蓉等,農(nóng)業(yè)水土工程方面劉全明等、徐英等均有一定成效,但多為局部估計的初步試探,而對IK法的運用規(guī)則及區(qū)域性IK法閾值及其函數(shù)關系的系統(tǒng)研究尚少見。本文以黃河大型灌區(qū)(河套)土壤水鹽監(jiān)測實驗的水土資源空間變異為案例,重點分析IK方程三項主要結構因子:指示閾值、指示變異函數(shù)和條件概率分布函數(shù)的變化趨勢、相互關系和取值規(guī)則,為IK法的應用提供參考。1ik平均估計的確定設在某實驗區(qū)D上取樣并測定其品位觀測值,若該區(qū)域上的邊界門限(即閾值)為Z,則D內(nèi)每一個樣品點x∈D上定義一個Z的階梯函數(shù),即指示函數(shù):i(x;Ζ)={1當x點上的觀測值Ζ(x)≤Ζ0當x點上的觀測值Ζ(x)>Ζ(1)D內(nèi)任一區(qū)域A∈D,低于閾值Z的觀測值Z(x)所占A區(qū)的比例表示為:?(A;Ζ)=1A∫Ai(x;Ζ)dx∈(2)式中?(A;Z)是關于Z(x)和Z的二元函數(shù),即小于閾值Z的全部i(x;Z),(x∈A)的平均值。指示函數(shù)i(x;Z)在給定閾值Z的條件下服從二項分布,期望值是:E{Ι(x;Ζ)}=1?Ρrob{Ζ(x)≤Ζ}+0?Ρrob{Ζ(x)>Ζ}=Ρrob{Ζ(x)≤Ζ}=F(Ζ)(3)式中F(Z)稱觀測值Z(x)的分布函數(shù)在閾值Z處的數(shù)值。當I(x+h;Z)和I(x;Z)是被矢量h分隔的兩個隨機變量,即兩個指示變換的數(shù)據(jù)時,則指示變異函數(shù)(IndicatorVariogrom)定義為:γΙ(h;Ζ)=12E{[Ι(x+h;Ζ)-Ι(x;Ζ)]2}(4)用指示變異函數(shù)可描述IK的空間變異程度,并可獲得條件概率分布函數(shù)(CCDF)的估計值:F*[z|(n)]=prob*[Ζ≤z|(n)]=1nn∑α=1i(Ζα;Ζ)?∈[0?1](5)這個估計能被接受的前提是:各樣本觀測值Zα=zα是獨立的,它們每一個對Z的作用是相同的,否則應采用不等權估計:F*[z|(n)]=prob*[Ζ≤z|(n)]=n∑α=1λα(z)?i(za；z)α=1,?,n(6)可以看出指示值i(zα;z)是用Zα-zα單獨估計Z時的ccdf:F*[z|(n)=zα]=prob*[Ζ≤z|Ζ=Ζα]=i(zα;z)α=1,?,n(7)故式(6)可以寫成:F*[z|(n)]=prob*[Ζ≤z|(n)]=n∑α=1λα(z)?F*(z|Ζ=Ζα)(8)上式說明:用多個觀測值對CCDF做出的估計等于用單個觀測值進行估計時所得CCDF的加權平均。指示函數(shù)值的方差S2(Z),即指示函數(shù)的基臺值C0,可從γ*I(h;Z)的圖形求得,故分布函數(shù)F*(z|(n))可從解以下二次方程得到F2(z|(n))-F(z|(n))-C0=0(9)同普通Kriging(OK)法一樣,再運用線性無偏最優(yōu)插值和估計方差最小的OK要求,結合指示函數(shù)的有關條件建立類似普通克立格(OK)方程的IK方程組式(10),可獲得待估域A的IK平均估計。{n∑β=1λβˉγi(xα,xβ;Ζ)+μ=ˉγi(xα,A;Ζ)n∑α=1λα=1(α=1,2,?,n)(10)用式(10)可求出這些權系數(shù)λα及未知點的i*(x;Z)估值,即Z(x)≤Z在該點出現(xiàn)的概率和整個CCDF的估值,可不必剔除特異值而獲得具有一定風險條件下的未知Z(x)的估計量及空間分布。2案例分析2.1土樣及測定方法研究區(qū)域設在黃河河套灌區(qū)內(nèi)的沙壕渠試驗區(qū)。試驗區(qū)采樣面積為4km2,用手持GPS接收機定位。以間隔h=200m(采樣點數(shù)為121個,中等田間尺度)的網(wǎng)狀進行采樣,采樣分兩層取土,深度為0～20、20～40cm。土樣的測定分析項目為含水率θ和電導率EC、pH值。其中EC的測定采用水土比5∶1浸提,電導法測定水溶性鹽總量。采樣點布置見圖1。2.2數(shù)據(jù)的比較及預處理IK法關鍵是閾值的合理確定,目前還沒有從理論上研究水鹽的閾值合理選擇問題。閾值的選擇主要是根據(jù)研究目標的專業(yè)門限要求、估計誤差和允許風險大小進行選擇,在缺乏經(jīng)驗時,有人建議用中位數(shù)作為初步估計參考值,經(jīng)作者初步實踐認為:中位數(shù)并非最佳閾值,估計結果具有平滑性,且接近OK法,也可能會產(chǎn)生IK法與OK法無明顯區(qū)別的誤導,為此應進行不同方案的比較及探索研究。表1水分數(shù)據(jù)的偏態(tài)系數(shù)接近于零,表明數(shù)據(jù)分布接近于正態(tài)分布的左右對稱,變異系數(shù)小反映水分變異性不大;表2鹽分數(shù)據(jù)的偏態(tài)系數(shù)均大于零,表明數(shù)據(jù)分布呈右偏態(tài),變異系數(shù)大反映鹽分變異性較大;水鹽數(shù)據(jù)的峰度系數(shù)均大于零,表明數(shù)據(jù)分布比正態(tài)分布高聳而狹窄。按照均值的3倍標準差鄰域判斷,發(fā)現(xiàn)鹽分數(shù)據(jù)有明顯特異值存在,0～20cm土層水分有特異值存在,適宜用IK法估計。對于大、中尺度樣本應根據(jù)不同作物、不同土壤和氣候及地理條件,統(tǒng)一選定合理閾值,如水分閾值可選用田間持水量、凋萎系數(shù)、適宜含水量等;鹽分閾值可采用FAO推薦的作物全生育期平均土壤濃度臨界值、作物能容忍的土壤最大鹽分值等。本文綜合考慮數(shù)據(jù)的分位數(shù)特點,水分初步采用20%、22%、25%、27%4個閾值;鹽分初步采用0.2、0.4、0.6、0.8、1.0、1.2mS/cm6個閾值進行參數(shù)的分析與比較,旨在探索選擇合理水鹽閾值的原則。下文僅列出0～40cm土層土壤的水鹽計算成果。3閾值與函數(shù)關系的研究3.1水鹽閾值對鹽脅迫下自適應的影響水鹽指示條件概率的計算采用ARCGIS9.2軟件的地質(zhì)統(tǒng)計學模塊完成。表3中水分預測平均條件概率隨著閾值的增加由99.2%減少至18.2%,反映了超過大閾值的概率逐漸減少;RMSE結果表明隨著閾值的增加,均方根誤差由9.5%增加至22.4%、51.1%,后又減小至39.0%,有先增加、后減小的變化趨勢;表3中鹽分預測平均條件概率隨著閾值的增加由82.0%減少至1.8%,反映了超過大閾值的概率逐漸減少;RMSE結果表明隨著閾值的增加,均方根誤差由38.4%增加至52.2%,后減小至9.4%,也有相似的先增加、后減小的趨勢。從圖2可以直觀地看出水鹽閾值與指示概率均值的關系:水分預測平均概率隨著閾值的增加而減少,減小的速度不同,超過22%閾值對應的指示概率均值減小速度加大。鹽分預測平均概率隨著閾值的增加而減少,減小的速度在0.6mS/cm閾值處變慢。水分預測概率均方根誤差隨著閾值的增加,先增加后減少,在25%～26%(接近中位數(shù)25.5%)區(qū)間概率預測誤差最大,接近50%,這將會降低風險評價與IK法估計的可靠性,在20%閾值處概率預測誤差最小,可見單從提高概率預測精度考慮,水分指示閾值應該盡量選擇低于中位數(shù)的偏小值。鹽分預測概率均方根誤差隨著閾值的增加,先增加后減少,在0.4mS/cm(接近中位數(shù)0.38mS/cm)處最大超過了50%,也將會降低風險評價與IK法估計的精度,可見單從提高概率預測精度考慮,鹽分指示閾值應該盡量選擇超過中位數(shù)的偏大值。圖3—6為不同閾值下土壤水鹽指示概率預測圖。分析后發(fā)現(xiàn)高閾值的大概率區(qū)均包含在低閾值的大概率區(qū)內(nèi),其分布位置隨閾值的增加而逐漸減少和集中,可用于不同閾值目標下的概率與風險評價。分析水分閾值25%、27%的指示概率預測圖,含水率超過閾值25%、27%的概率在試驗區(qū)中部及左上角較大,含水率超過27%的概率最大值為50%,有4個明顯大概率區(qū)均包含在超過25%的大概率區(qū)域內(nèi)。分析鹽分閾值0.4、0.6mS/cm相應的指示概率預測圖,可以了解鹽分積聚的可能性,如電導率超過閾值0.6mS/cm的概率最大為36%,分布在試驗區(qū)北部兩側及西南角,且均包含在超過閾值0.4mS/cm的大概率區(qū)域內(nèi)。3.2中性關系變異函數(shù)性變異性函數(shù)的檢測由表3可見:當水分閾值由小逐漸增大時,指示變異函數(shù)值有一種由小到大、再到小的趨勢,轉折點在25%處,可見接近中位數(shù)閾值的指示變異函數(shù)值最大,其指示變異函數(shù)結構性最差,所反映的指示概率的空間變化性最大,由此得到的指示概率預測誤差最大,這與上文RMSE在閾值接近中位數(shù)處最大的結果相一致。鹽分指示閾值與指示變異函數(shù)的關系同水分,指示變異函數(shù)值轉折點在0.4mS/cm處。3.3閾值對結構性分析結果的影響從表5分析水分指示變異函數(shù)的結構性,閾值22%的指示變異函數(shù)結構性表現(xiàn)為純隨機,因為其C0/(C0+C)=1;結構性較強對應的閾值為20%、25%、27%,隨著閾值增大,結構性有由強到弱,再由弱到強的變化趨勢,為非降函數(shù),其他兩層指示變異函數(shù)的結構性分析具有相同的結論。從表6分析鹽分指示變異函數(shù)的結構性,無純隨機現(xiàn)象,隨著閾值增大,結構性有由強到弱,再由弱到強的變化趨勢,0.4～0.6mS/cm區(qū)間內(nèi)閾值對應的指示變異函數(shù)結構性最差,其指標值接近于1。4與其他方法插值對比本文選用含水率25%、電導率0.6mS/cm分別作為水鹽的指示閾值進行IK法區(qū)域水鹽估計,首先借助于IK方程組求解?(A;Z)的線性估計量?*(A;Z);然后計算區(qū)域A平均值的估計值:Z**(A;Ζe)=n∑e=1[?*(A;Ze)×Z*(A;Ze)],最后與OK法、傳統(tǒng)反距離加權平均法(簡寫IDW法)的插值結果做了對比。圖7—8為土壤水鹽IK法插值圖,可以反映試驗區(qū)的水鹽分布條件。水分分布比較連續(xù),有3處明顯濕潤區(qū);鹽分在西南、東北及試區(qū)中部(對角線方向)有明顯條帶狀積聚。從表7可以分析3種方法的水分預測均值與實測值基本相同,在121個采樣點上的預測平均誤差為1.6%,說明了IK法同樣具有精確的插值功能。IK法的估計方差、變異系數(shù)最大,但卻小于采樣點對應參數(shù),說明IK法對水分普通克立格插值平滑效應有一定程度的消弱。從表8可以分析3種方法的鹽分預測均值與實測值基本相同,在121個采樣點上的預測平均誤差為0.230mS/cm,說明了IK法同樣具有精確的插值功能。IK法的估計方差、變異系數(shù)大于OK法,但卻小于采樣點對應參數(shù),說明IK法對鹽分普通克立格插值平滑效應有一定程度的消弱。5加指示概率估計的優(yōu)缺點(1)IK法的難度在于閾值的合理選擇,不僅涉及研究目標和專業(yè)門限要求,也與研究允許風險概率及變異程度有關,目前尚無理論分析方法,多數(shù)情況仍然需要試算探求。(2)經(jīng)過模型與參數(shù)分析,發(fā)現(xiàn)水鹽指示Kriging平均概率與指示閾值、指示變異函數(shù)的關系為:隨著閾值的增加指示Kriging平均概率減少;當閾值由小逐漸增大時,指示變異函數(shù)值有一種由小到大,再到小的趨勢。(3)作者認為選擇中位數(shù)作閾值并非最佳估計,接近中位數(shù)閾值的指示變異函數(shù)值最大,其指示變異函數(shù)結構性最差,所反映的指示概率的空間變化性最大,由此得到的指示概率預測誤差最大。本