萊州灣南岸土壤鹽分空間變異性指示克立格閾值研究

（河南工程學院，河南鄭州市，491191） 畢延鳳

萊州灣地區(qū)有濃度最高、分布最廣、儲量最大的濱海相地下鹵水，是我國重要的海鹽生產基地和海水養(yǎng)殖基地，還是山東省糧食高產區(qū)和全國蔬菜市場的重要基地。這些經濟農業(yè)活動擴大了咸水入侵范圍，2007年7月開始，在國家908專項支持下，對萊州灣沿岸海水入侵嚴重地區(qū)進行了為期兩年的咸水入侵災害監(jiān)測，結果顯示萊州灣地區(qū)咸水入侵范圍已達2 500km2。咸水入侵導致大量土地發(fā)生鹽漬化，植物群落和生態(tài)系統(tǒng)相應發(fā)生改變。

土壤鹽分高值區(qū)以特異值形式存在，數(shù)據(jù)往往不呈正態(tài)分布，導致變異函數(shù)不穩(wěn)健[1]，傳統(tǒng)的高斯克里格估值存在困難。Journel于1983年創(chuàng)立的指示克里格(IK)，不要求數(shù)據(jù)服從正態(tài)分布且不受特異值的影響，它對連續(xù)數(shù)據(jù)通過閾值變換轉變成0和1的二進制數(shù)據(jù)，并以指示變異函數(shù)為基本結構分析工具，預測精度超過了普通克里格法[2]。IK法自創(chuàng)立以來廣泛應用于地質[3]、土壤[4、5]、水文[6、7]等領域的風險預測研究中，但多為局部估計的初步探討，對IK法閾值的選擇及與函數(shù)關系的系統(tǒng)研究尚少見[8、9]。

以萊州灣南岸濰坊北部平原的土壤鹽分空間變異為案例，以ArcGIS9.3、GS+7.0為操作平臺，重點分析了不同分位數(shù)閾值下鹽分空間結構、IK插值精度及指示閾值與條件概率分布函數(shù)的相互關系和取值規(guī)則，為IK法的應用提供參考。

1 研究區(qū)概況

研究區(qū)位于渤海萊州灣南岸，濰坊北部地區(qū)，面積1 800km2。地勢南高北低，地形平坦廣闊，屬濱海堆積平原地貌。自西向東依次有彌河、白浪河、濰河和膠萊河由南向北注入渤海，海岸線為東南－西北走向，呈弧形曲線狀。本區(qū)屬干旱半干旱氣候，年均降雨量652.8mm，年均蒸發(fā)量1 774mm，蒸發(fā)量遠大于降水量。

自南至北分布著棕壤、褐土、潮土、砂姜黑土和鹽土5大土類，下分15個亞類、34個土屬、110個土種。其中鹽土類面積約149×104畝，占可利用土壤面積的7.43%，集中分布在市域北部濱海地帶，即壽光、寒亭、昌邑3縣、區(qū)的北部。鹽土土類下分濱海潮鹽土亞類及2個土屬、8個土種。

2 材料與方法

自西向東沿垂直海岸線方向共布設了13條監(jiān)測剖面，每條剖面上設置監(jiān)測點5～7個不等，均由海水入侵區(qū)過渡至非海水入侵區(qū)，共計80個，各相鄰剖面及相鄰點之間的直線距離約2～3km（圖1）。取樣點大部分選在農田，亦有部分采樣點選在鹽荒地、林地及草地。土樣的采集利用手動土壤采樣器，每個監(jiān)測點取地表以下1m土樣，分6層，分別為0～10cm、10～20cm、20～40cm、40～60cm、60～80cm及80～100cm，每份樣品重量約1kg，編號，帶回實驗室，風干、研磨、過篩，進行土壤全鹽含量的測定[10]。

指示克里格法的原理詳見文獻[5、7]，首先確定評價指標的閾值，在此選取20、30、40、50、60、70、80分位數(shù)的鹽分閾值，并為評價指標確定指示函數(shù)，將采樣點數(shù)據(jù)進行二態(tài)指示變換（指示變換值1或0），然后將指示變換值在GS+7.0中進行變異函數(shù)模擬計算，最后將得到的最佳變異函數(shù)模型參數(shù)輸入ArcGIS9.3進行IK插值，得到土壤鹽分滿足一定閾值下的概率空間分布圖。

經典統(tǒng)計采用SPSS，變異函數(shù)的模型擬合采用GS+7.0，概率空間分布圖的繪制以及平均概率、預測概率均方根誤差（RMSE）的計算采用ArcGIS9.3。

3 結果與分析

3.1 土壤鹽分含量的描述性統(tǒng)計

從表1中K-S_P值可知，20～40cm和80～100cm土壤鹽分經自然對數(shù)轉換后K-S檢驗的p值＞0.05（雙截尾），服從對數(shù)正態(tài)分布，其余土層以及0～100cm加權平均含鹽量（表1最后一行）的分布均呈明顯的偏態(tài)分布。土壤剖面中，上層0～10cm，10～20cm和底層80～100cm土層土壤鹽分含量的最大值和最小值差別最大。上層0～10cm和10～20cm鹽分均值為中度鹽漬化，往下層鹽分有增加的趨勢，呈現(xiàn)重度鹽漬化。標準差也呈現(xiàn)往下層增加的趨勢。從各層變異系數(shù)來看，范圍為2.767～3.699，均為強變異，20～40cm變異系數(shù)最大，其他各層總體呈現(xiàn)往下層遞增的趨勢。

在秋季持續(xù)干旱的地表蒸發(fā)作用下，表層0～10cm全鹽量均值大于10～20cm，而總體剖面呈底聚類型，這與萊州灣南岸廣大平原區(qū)除緊鄰萊州灣的沿海狹窄條帶，由于鹵水開發(fā)逐年降低的地下水位，干旱少雨的氣候和灌溉水的緊缺密切相關。同時也體現(xiàn)了鹽漬化程度正趨于和緩。

3.2 不同閾值下土壤鹽分空間結構的變化

用0～100cm加權平均含鹽量為研究對象，由于數(shù)據(jù)不呈正態(tài)分布，擬對其進行指示變異函數(shù)的研究。

利用表2中的鹽分閾值，分別對土壤鹽分實測值進行二態(tài)指示變換，在GS+7.0中設有效滯后距67 667.22m，步長大小為6 766.72m，步長數(shù)為10，模擬獲得不同閾值下的鹽分指示變異函數(shù)最佳模型及其參數(shù)。

從表中可以看出，土壤鹽分的空間相關性總體上隨著鹽分閾值的增大而增強，在第50百分位數(shù)之前均呈中等強度的空間相關性，鹽分從第60百分位數(shù)開始大多呈現(xiàn)強烈的空間相關性。RSS都很小，尤其是在第50分位數(shù)之前，在第30分位數(shù)處達最小。模型的決定系數(shù)R2在第50百分位數(shù)之前在0.758～0.995之間，模型擬合精度高，鹽分第60百分位數(shù)之后R2在0.220～0.386之間，擬合精度較低。

表1 1m以內土體各層及各層平均土壤鹽分含量的經典統(tǒng)計

表2 不同鹽分閾值下土壤鹽分的指示變異函數(shù)模型

表3 土壤鹽分指示閾值與指示變異函數(shù)的關系

3.3 指示閾值與指示變異函數(shù)的關系

由表3可見，當鹽分閾值由小逐漸增大時，指示變異函數(shù)值有一種由小到大、再到小的趨勢，轉折點在第50分位數(shù)（鹽分閾值0.13%）處，可見中位數(shù)閾值附近的指示變異函數(shù)結構性最差，其反應的指示概率的空間變化性最大，由此推測其指示概率的預測誤差也將最大。

3.4 指示閾值與指示概率函數(shù)的關系

指示概率均值隨著鹽分閾值的增大逐漸減少，由83.4%減少到27.6%，從第60分位數(shù)（鹽分閾值0.16）開始減小速度變慢；均方根預測誤差（RMSE）有先增加后減小的趨勢，由38.4%迅速增加至45.9%，后隨著閾值增大緩慢減小至38.7%，在第50分位數(shù)（鹽分閾值0.133）處概率預測誤差最大，與上文提及的指示變異函數(shù)在此值處最大一致，采用此閾值會降低風險評價與IK法估計的可靠性，在第20分位數(shù)處（鹽分閾值0.09）概率預測誤差最小?？梢?，單從IK法插值精度考慮，中位數(shù)鹽分閾值并不是最佳閾值。這點與楊奇勇等對山東禹城市[9]、劉全明等對內蒙古河套灌區(qū)[8]利用IK法進行的土壤鹽分空間變異性研究結果一致。但在楊與劉等的研究中，由于其研究區(qū)中位數(shù)以下的鹽分閾值小于0.1%，屬非鹽漬化，其概率分布圖對鹽漬化改良參考意義不大，因此建議研究區(qū)鹽分閾值應盡量在精度容許情況下選取超過中位數(shù)的偏大值。而本文的研究中，中位數(shù)以下的30和40分位數(shù)處(鹽分閾值分別是0.10%和0.11%)均屬輕度鹽漬化范疇，尤其是30分位數(shù)處，半變異函數(shù)模型擬合精度最高，土壤鹽分呈中等強度的空間相關性；超過中位數(shù)閾值的RMSE相對較小的第70、80百分位數(shù)處，土壤鹽分呈強烈的空間相關性，而其半變異函數(shù)的擬合精度較低。

3.5 指示概率分布圖

在鹽堿土壤改良過程中，更感興趣的往往不是某一點處土壤鹽分的具體值，而是鹽分大于某一閾值在空間上的分布概率。用IK法做出不同鹽分閾值下的指示概率分布圖如圖2。

在圖2中，隨著鹽分閾值的增大，高概率面積逐漸減小，低概率面積逐漸增大，導致平均概率逐漸由鹽分閾值0.09%的83.1%減小到鹽分閾值0.35%的27.6%，概率值分布的偏度由鹽分閾值0.09%的負偏（偏度為-0.9 775）過渡到鹽分閾值0.35%的正偏（偏度為0.6354），整體呈逐漸增大趨勢；概率值分布的峰度在鹽分閾值0.09%、0.10%、0.11%和0.13%處普遍偏大，在中位數(shù)閾值0.13%處達到最大，而在鹽分閾值0.16%和0.22%處急劇降低，在鹽分閾值0.35%處又開始回升，峰度隨著鹽分閾值的變化規(guī)律與研究區(qū)土壤鹽分的分布狀況密切相關。從概率值分布的偏度和峰度看，在中位數(shù)閾值0.13%處最接近正態(tài)分布（偏度-0.929，峰度3.252）。隨著鹽分閾值的增大，主要概率分布區(qū)間也由鹽分閾值0.09%的0.9～1.0過渡到鹽分閾值0.35%的0.0～0.1。當鹽分閾值0.10%時，高概率區(qū)主要分布在廣大的西部以及北部地區(qū)，濰坊北部是主要鹽場所在地，而西部則主要與廣泛分布的大棚蔬菜種植區(qū)導致的次生鹽漬化密切相關。在鹽分閾值0.22%處，鹽分高概率區(qū)主要分布在北部鹽田區(qū)以及膠萊河和白浪河河谷地帶。某一鹽分的高概率分布區(qū)域一般包含在小于它的上一個鹽分閾值的高概率分布范圍內，相應的低概率分布區(qū)域一般包含在大于它的上一個鹽分閾值的低概率分布范圍內。鹽分閾值可根據(jù)不同閾值目標下特定區(qū)域的概率與風險進行選擇。

萊州灣南岸地區(qū)鹽漬土的形成既有自然因素也有人為因素，較高的氣溫及較大的蒸發(fā)量使得受海水（或鹵水）浸染的礦化度較高的地下水通過蒸發(fā)作用及毛細作用將水中鹽分聚積在地表及地表下的土層中，形成了鹽漬土。除自然因素外，人類活動是造成土壤次生鹽漬化的重要原因。萊州灣是山東省糧食高產區(qū)和全國蔬菜市場和水產養(yǎng)殖的重要基地，自然條件優(yōu)越。但過大的人口密度（是全國人口平均密度的4倍多）及非農業(yè)用地比例的增加使得人口環(huán)境壓力十分突出。近年來隨著沿海改革開放戰(zhàn)略的實施，在各縣市經濟快速發(fā)展的同時，水資源的需求量猛增。由于過度開采地下水，使得地下水位以1～3m/yr的速度迅猛下降，形成大于2 500km2的地下水位負值區(qū)。此外萊州灣地區(qū)是我國重要的海鹽生產基地和海水養(yǎng)殖基地，海水養(yǎng)殖和抽鹵曬鹽等經濟活動把大量海水和鹵水引入陸地，人為地促進海水或咸水向內陸運移5～15km2，擴大了咸水入侵范圍。咸水入侵一方面直接浸漬土顆粒使其鹽漬化；另一方面，為了農田保產，人們不得不使用一部分含鹽量較高的地下半咸水進行灌溉，在沒有良好的排灌系統(tǒng)的情況下長期使用高礦化度的水進行灌溉，鹽分不斷在土壤表層聚集，使耕層土壤鹽分含量增加，必然導致了次生鹽漬化的發(fā)生。據(jù)各縣市區(qū)的統(tǒng)計資料，目前由海水入侵造成的次生鹽漬化面積已達2.5×104hm2，受害土地面積達1.01×105hm2。土壤鹽漬化的結果是土壤結構變差，理化性能降低，土地資源退化，嚴重的喪失生產能力，部分農田因地下水變咸而被棄耕成為荒地，栽培作物逐漸為鹽生植被所替代，從而進一步加劇了耕地資源的危機，也增加了區(qū)域生態(tài)恢復和生態(tài)開發(fā)的難度。

4 結論

萊州灣南岸濰坊北部的廣大平原地區(qū)土壤鹽分呈偏態(tài)分布，利用IK法，選擇不同分位數(shù)鹽分閾值，對土壤鹽分空間結構、鹽分閾值與指示概率和指示概率均方根預測誤差的關系以及IK指示概率插值圖進行了研究。

隨著土壤鹽分閾值的增大，鹽分空間相關性總體呈增強趨勢，擬合精度降低，指示變異函數(shù)值有一種由小到大、再到小的趨勢，轉折點在第50分位數(shù)處，可見中位數(shù)閾值附近的指示變異函數(shù)結構性最差，其反應的指示概率的空間變化性最大，導致其指示概率均方根預測誤差也最大，采用此閾值會降低風險評價與IK法估計的可靠性。中位數(shù)以下的30和40分位數(shù)處(鹽分閾值分別是0.10%和0.11%）均屬輕度鹽漬化范疇，尤其是30分位數(shù)處，半變異函數(shù)模型擬合精度最高，土壤鹽分呈中等強度的空間相關性；超過中位數(shù)閾值的RMSE相對較小的第70、80百分位數(shù)處，土壤鹽分呈強烈的空間相關性，而其半變異函數(shù)的擬合精度較低。

從不同鹽分閾值下的指示概率分布圖可見，隨著鹽分閾值的增大，高概率面積逐漸減小，低概率面積逐漸增大，鹽分預測平均條件概率隨著閾值的增大逐漸減少。主要概率分布區(qū)間由鹽分閾值0.09%的0.9～1.0過渡到鹽分閾值0.35%的0.0～0.1。概率值分布的偏度由鹽分閾值0.09%的負偏過渡到鹽分閾值0.35%的正偏，峰度先增大在中位數(shù)閾值0.13%處達到最大，后急劇降低，在鹽分閾值0.35%處又開始回升，在中位數(shù)閾值0.13%處最接近正態(tài)分布。某一鹽分的高概率分布區(qū)域一般包含在小于它的上一個鹽分閾值的高概率分布范圍內，相應的低概率分布區(qū)域一般包含在大于它的上一個鹽分閾值的低概率分布范圍內。鹽分閾值可根據(jù)不同閾值目標下特定區(qū)域的概率與風險進行選擇。