常熟地區(qū)水稻土飽和導水率的間接方法研究

1、常熟地區(qū)水稻土飽和導水率的間接方法研究基金項目： 研究工作得到所知識創(chuàng)新工程領域前沿項目（ISSASIP0719），支撐計劃課題（2006BAD10A06-03），973計劃課題（2005CB121103），863專題課題（2006AA10Z208）的資助。作者簡介：張均華(1983-)，博士研究生，主要從事區(qū)域土壤水分溶質運移模擬研究。Email:zhangjh*通訊作者：劉建立，Email:jlliu張均華1,2 劉建立1* 朱安寧1（1. 中國科學院南京土壤研究所，南京210008；2. 中國科學院研究生院，北京100049）摘要：作為直接試驗的一種替代方法，利用土壤基本物理性質通過土壤

2、轉換函數預測飽和導水率簡便易行，成本低廉，并且預測精度能滿足實際研究的需要。本研究利用目前得到較多應用的9種基于多元回歸分析建立的轉換函數來構建、校正預測土壤飽和導水率的經驗公式，并與人工神經網絡方法相比較。結果表明，人工神經網絡總體預測效果要優(yōu)于基于多元回歸分析建立的轉換函數，并且Cosby(1984)在輸入參數較少的基礎上預測飽和導水率精度最高。本文以Cosby(1984)預測常熟水稻土壤飽和導水率，進一步利用GIS的空間描述能力與函數的定量分析能力，得到區(qū)域尺度飽和導水率的分布狀況，為該地區(qū)區(qū)域尺度數值模擬的運行提供基礎參數支持。關鍵詞：土壤飽和導水率；多元回歸函數；人工神經網絡；GIS

3、中圖分類號：S152.7 文獻標識碼：A飽和水力傳導率（Ks）是重要的土壤水力學參數，是定量分析土壤水分入滲、地表產流等過程的基礎。土壤特別是表層土壤飽和導水率空間變異性非常強烈，需要采集大量樣本或進行田間現(xiàn)場測定才能得到一個較大空間尺度上的飽和導水率空間分布1。但是通常認為室內填裝樣本的測定結果無法表征田間原始狀態(tài)。田間原位測定土壤飽和傳導率的方法主要有Guelph滲透儀、流速滲透儀、圓盤滲透儀和雙管法。由于時間，人力成本等問題，原位測定僅能用于有限的點，無法滿足較大空間尺度（如區(qū)域上）的研究。已有的研究表明，土壤粒徑組成，容重，有機質含量等基本土壤性質與飽和導水率之間存在一定關系，可以通過

4、常規(guī)土壤調查中積累的大量基本土壤性質資料間接估算或預測區(qū)域尺度飽和導水率的空間分布2-5。并且在許多實際問題中，間接方法的預測結果也能提供足夠精度的Ks預測值。但到目前為止，在長江中下游典型水稻土上的應用很少有報道。本文即由此出發(fā)，以江蘇常熟市為研究區(qū)域，利用目前得到較多應用的基于多元回歸分析的土壤轉換函數(Pedo-transfer functions，簡稱PTFs)來構建、校正預測土壤飽和導水率的經驗公式，并與人工神經網絡方法(ANN)相比較。進一步利用GIS的空間描述能力與函數的定量分析能力，對區(qū)域尺度的飽和導水率做出定量描述與評價，為區(qū)域尺度數值模擬的運行提供基礎參數支持。1材料與方法

5、1.1數據來源土壤轉換函數方法普遍存在明顯的地域局限性。為克服轉換函數的區(qū)域適用性問題，建立PTFs時除采用45個常熟水稻土樣本外，還利用地跨河南封丘、長垣、延津、原陽四縣的天然文巖渠流域潮土樣本169個，以及隨機選取UNSODA數據庫6中各種質地類型的土壤樣本289個，共計503個樣本數據。水稻土和潮土飽和導水率由ALFranklin Engineers生產的圓盤滲透儀(CSIRO disc permeameter)現(xiàn)場測定。土壤基本性質包括粒徑組成、容重、有機質含量等，常規(guī)方法室內測定。土壤特征水分含量數據(如飽和含水量和-33 kPa含水量)根據土壤轉換函數預測7。1.2分析方法1.2.

6、1土壤轉換函數選用上世紀八十年代以來得到較多應用的9種PTFs，從503個樣本數據中隨機選取336個構建PTFs。通過與實測值比較，重新優(yōu)化函數系數。剩余167個樣本檢驗函數的預測能力，以期得到一個飽和導水率預測值與實測值最為接近的PTFs。在對函數系數進行優(yōu)化處理時采用了自行編制的基于改進的Levenberg-Marquard算法的FORTRAN程序。ANN方法通過非線性迭代可以獲得輸入基本性質與輸出飽和導水率之間（理論上）的最優(yōu)關系。本研究應用無內部反饋的誤差反向傳播前向神經網絡（BP網絡）對隨機選取的336個樣本數據進行訓練，剩余167個樣本數據進行驗證。進而對ANN與PTFs進行比較。

7、研究表明，Ks多服從正態(tài)分布，因此本文在對Ks實測值進行對數轉換的基礎上進行Ks的預測、比較。表1 用于本研究的預測飽和導水率的土壤轉換函數Table 1 The Pedo-transfer functions in this study函數來源Function Source自變量個數Independent Variable輸入變量Input variableCosby(1984)82Sand, Clay,Puckett(1985)9-101ClaySaxton(1986)9-103Sand, Clay, qSDanePuckett(1992)111ClayRawls(1992)123Sand

8、, Clay, PorosityCampbellShiozawa(1994)9132Sand, ClayWosten(1997)144Silt, Clay, BD, OMWosten(1999)154BD, OM, Sand, ClayKC方程(1999)164Sand, Clay, qS, q33注：Sand土壤砂粒含量（%）；Silt土壤粉粒含量（%）； Clay土壤黏粒含量（%）；BD土壤容重(g/cm3)；OM土壤有機質（%）；Porosity土壤孔隙度（%）； q33 土壤水勢為-33kPa時土壤水分含量（cm3/cm3）；qS 土壤飽和水分含量（cm3/cm3）。1.2.2 轉換函

9、數預測精度評價本文根據均方根誤差、確定性系數以及Akaike信息準則評價土壤轉換函數的預測精度。(1)均方根誤差(Root mean square error, RMSE)，計算Ks預測值的平均精度： （1）(2)確定性系數(Coefficient of Determination,R2),反映預測值和實測值相關性的重要指標： （2）（3）Akaike信息準則（Akaikes information criterion, AIC）11,評價方法的實用性。一種方法中需要涉及的參數越多，其實用性就會下降，最好的方法應該有最小的AIC。 （3）公式(1)至（3）中，表示實測值，表示預測值，N表示樣本

10、數，np為涉及的參數數量。1.2.3飽和導水率的區(qū)域尺度擴展在對比ANN與PTFs基礎上，利用預測Ks效果最好的方法預測常熟水稻土耕層土壤204個采樣點的Ks。結合GIS的空間描述能力與函數的定量分析能力，進行點數據的區(qū)域尺度擴展，以期為區(qū)域尺度數值模擬的運行提供基礎參數支持。圖1為常熟市204個采樣點的空間分布1 張福鎖,朱兆良. 農田生態(tài)系統(tǒng)氮素行為與氮肥高效利用。圖1 常熟市204個采樣點空間分布Fig.1 The spatial distribution of Changshu in 204 samples2結果與分析2.1基于多元回歸分析建立的土壤轉換函數預測飽和導水率精度通過9種P

11、TFS預測Ks對比實測值分析發(fā)現(xiàn)，各函數預測精度總體不高。AIC值平均2596.141，RMSE平均2.62（cm.d-1），R2平均0.49. 可以看出幾乎所有函數模型均用到粒徑分布數據，這些土壤物理性質比較容易獲得，所以為預測Ks提供可能。但也有不少函數預測值與實測值差距較大，需要重新進行數據擬合，優(yōu)化函數系數。2.2基于多元回歸分析的土壤轉換函數優(yōu)化函數的可靠性通過對不同于構建函數時的獨立序列的預測能力來評價17。本研究在兼顧各種質地屬性的基礎上，以隨機選取的336個樣本構建PTFs，剩余167個樣本檢驗其可靠性。這樣通過數據擬合得到的優(yōu)化系數能更好的反映現(xiàn)有土壤飽和導水率。KC方程在公

12、式中沒有系數，不用進行優(yōu)化設置。表2 調整系數后多元回歸分析預測167個樣本點飽和導水率與神經網絡方法比較Table 2 Comparison of Predictive value between adjusted PTFs and ANN in 167 samples飽和導水率saturated hydraulic conductivity多元回歸分析優(yōu)化系數后比較Comparison of Ks by adjusted PTFs人工神經網絡方法比較Comparison of Ks by ANN AICRMSE(cm.d-1)R2AICRMSE(cm.d-1)R2Cosby(1984)71

13、7.7740.6560.651592.7360.5870.794Puckett(1985)727.3660.6790.577618.8120.6620.708Saxton(1986)717.9010.6820.590603.0850.6240.736DanePuckett(1992)727.3660.6990.498618.8120.6620.708Rawls(1992)697.6720.6140.600587.2240.5950.727CampbellShiozawa(1994)717.7740.6770.641592.7360.5870.794Wosten(1997)702.1800.61

14、50.495581.4370.5780.764Wosten(1999)1072.6591.8750.241610.2100.5760.746KC方程(1999)1195.9631.7290.598603.6430.6250.744表2為優(yōu)化系數后PTFs生成的Ks預測值與ANN方法生成的Ks對比實測值分析。可以看出，PTFs調整系數后AIC比原來有很大程度的降低，線性相關性也有很大提高，Ks預測值與實測值之間差異性顯著下降。Cosby(1984), Rawls(1992) 、CampbellShiozawa(1994)、Wosten(1997)這四個函數預測Ks效果較好。Wosten(1997

15、)輸入參數比較多，理論上參數越多，預測效果越理想。但隨著輸入參數的增多，采用土壤轉換函數的主要優(yōu)勢（即實用性）也會減弱。所以本研究選擇Cosby(1984)作為預測Ks的最佳函數。圖1-a為Cosby(1984)飽和導水率實測值與預測值之間的對比分析。圖2 Cosby(1984)多元回歸函數與人工神經網絡對飽和導水率的預測值與實測值比較Fig. 2 Comparison of Ks between ANN and Multiple regression analysis of Cosby(1984)通過與優(yōu)化后PTFs預測能力比較，可以發(fā)現(xiàn)ANN方法總體上精度更高。其中以Cosby(1984)

16、，Rawls(1992)，Wosten(1997)參數為輸入變量的ANN預測效果比較好。Cosby(1984)輸入參數少，在保證預測精度下能有較高的實用性，所以本研究選擇Cosby(1984)的輸入參數作為ANN方法的輸入變量，即根據砂粒含量與黏粒含量預測Ks。圖1-b為ANN方法根據砂粒含量與黏粒含量建立的Ks預測值對比實測值分析。2.3土壤基本性質對飽和導水率影響性分析通過相關性分析，可以識別出影響Ks的關鍵因素。表3為土壤基本性質與Ks相關性分析。可以看出Ks多與砂粒含量、孔隙度成正相關，而與粉粒、黏粒含量，-33Kpa時水分含量成負相關，這是因為大孔隙占飽和水運動的大多數。土壤質地和結

17、構與導水率也有直接關系，砂質土壤通常比細質土壤具有更高的飽和導水率，同樣，具有穩(wěn)定團粒結構的土壤，比不穩(wěn)定團粒結構的土壤，傳導水分要快得多，后者在潮濕時結構被破壞，黏粒和粉砂粒阻塞較大孔隙的連接通道。表3轉換函數所用參數與飽和導水率相關性分析Table 3 Correlation Analysis between parameters and Ks土壤基本性質soil physical properties土壤飽和導水率saturated hydraulic conductivitySand0.668*Silt-0.645*Clay-0.467*BD-0.183*Porosity0.428*O

18、M-0.047q33-0.423*qS-0.091注：*表示置信度達95%，*表示置信度達到99%由表3還可以看出Ks與土壤粒徑分布數據相關性較好，與容重，孔隙度等數據相關程度次之。所以本研究用Cosby(1984)預測Ks在實際應用中是可行的。圖3為以Cosby(1984)的輸入參數作為ANN輸入變量預測常熟水稻土15個采樣點Ks與實測值比較。RMSE為0.192，R2為0.937,可見以砂粒含量與黏粒含量預測Ks在較大尺度研究中能滿足精度需要。圖3 神經網絡方法-Cosby(1984)預測常熟15個樣本點飽和導水率Fig. 3 Predicted value of Ks in Changs

19、hu by ANN-Cosby(1984) in 15 samples2.4飽和導水率的區(qū)域尺度擴展圖4 常熟市飽和導水率空間分布Fig. 4 Distribution of Saturated Hydraulic Conductivity in Changshu通過上述分析，可以發(fā)現(xiàn)ANN方法預測Ks效果總體上要優(yōu)于PTFs。其中根據砂粒含量與黏粒含量預測的Ks與實測值最為接近。所以本研究以蘇州常熟市范圍內耕層土壤204個采樣點的土壤粒徑分布數據為基礎，以Cosby(1984)的輸入參數作為ANN輸入變量，預測采樣點Ks。通過對采樣點Ks進行地統(tǒng)計學插值，將點源數據擴展到區(qū)域尺度。圖4為常熟

20、地區(qū)耕層Ks的空間分布狀況。從圖中可以看出粒徑分布的空間狀況，西部砂粒含量高于東部，黏粒含量與之相反，而Ks的空間分布規(guī)律與砂粒相似，這與表3相關分析中所揭示的結果相符合。3 結論土壤轉換函數具有易獲取輸入數據，節(jié)省調查采樣成本，計算簡單等特點。本研究通過構建、校正PTFs，并與ANN方法比較，探求最大區(qū)域適用性的土壤轉換函數。進一步的利用GIS的空間描述能力與函數的定量分析能力，對區(qū)域尤其是水稻土的Ks做出定量描述與評價，這對于該地區(qū)區(qū)域尺度數值模擬運行中基礎參數的獲得有非常重要的實際意義。通過分析發(fā)現(xiàn)PTFs與ANN中預測Ks效果最好的均為Cosby（1984）。相關性分析中，Ks與土壤粒

21、徑分布數據相關程度最高，所以本研究選擇Cosby（1984）預測不同質地屬性的Ks能滿足精度需要。Cosby（1984）輸入參數較易獲得，這就為轉換函數的實際應用提供了可能。而ANN對Ks的預測效果總體優(yōu)于PTFS，這是由于ANN通過非線性迭代可以獲得輸入基本性質與輸出Ks之間的最優(yōu)關系。 參考文獻1 Klute A, Dirken C. Hydraulic conductivity and diffusivity: Laboratory methods. In: Klute A(ed.) Methods of soil analysis. Part 1. Physical and miner

22、alogical methods. 2nd edition. Soil Science Society of America, Madison, Wis. 1986, 687734.2黃元仿, 李韻珠. 土壤水力性質的估算土壤轉換函數J. 土壤學報, 2002, 39(4): 517-523.3 van Genuchten, M.T., 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soilsJ. Soil Sci. Soc. Am. J. 44, 892898.4朱安

23、寧, 張佳寶, 陳效民等. 封丘地區(qū)土壤傳遞函數研究J. 土壤學報, 2003, 40(1): 53-58.5劉建立, 徐紹輝, 劉 慧. 估計土壤水分特征曲線的間接方法研究進展J. 水利學報. 2004, 2: 68-76.6 Minasny, B., A.B. McBratney. 2001. The Australian soil texture boomerang: a comparison of the Australian and USDA/FAO soil particle-size classification systems. Aust. J. Soil Res., 2001

24、, 39, 14431451.7 Saxton, K.E. and W.J. Rawls. 2004. Soil water characteristic estimates by texture and organic matter for hydrologic solutions. Proc. , Soil Sci. Soc. Amer. International Conference, Seattle, WA Nov. 1-3, 2004.8 Cosby, B.J., G.M. Hornberger, R.B. Clapp, T.R. Ginn. 1984. A statistical

25、 exploration of the relationship of soil moisture characteristics to the physical properties of soils. Water Resources research, 20(6):682-690.9 Minasny, B., A.B. McBratney. 2000. Evaluation and development of hydraulic conductivity pedotransfer functions for Australian soil. Australian Journal of S

26、oil Research, 38: 905-926.10 Puckett, W.E., J.H. Dane, B.F. Hajek. 1985. Physical and mineralogical data to determine soil hydraulic properties. SSSAJ, 49:831-836.11 Minasny, B., A.B. McBratney. 2000. Evaluation and development of hydraulic conductivity pedotransfer functions for Australian soil. Au

27、stralian Journal of Soil Research, 38: 905-926.12 Rawls, W.J., L.R. Ahuja, D.L. Brakensiek. 1992. Estimating soil hydraulic properties from soils data. In: M. Th.van genuchten, F.J. Leij, and L.J. Lund (eds) Indirect Methods for Estimating the Hydrauli Properties of Unsaturated Soils. University of

28、California, Riverside,329-340.13 Campbell, G.S., S. Shiozawa. 1992. Prediction of hydraulic properties of soils using PSD and bulk density data. pp. 317-328. In M. Th.van Genuchten et al. (eds) Proceedings of Inter. Workshop on Indirect Methods for Estimating the Hydraulic properties of Unsaturated

29、Soils.14 Wosten, J.H.M. 1997. Pedotransfer functions to evaluate soil quality. In Gregorich E.G., Carter M.R. (eds), Soil Quality for Crop Production and Ecosystem Health. Development in Soils Science, Vol. 25, Elsevier, Amsterdam, pp. 221-245.15 Wosten, J.H.M., A. Lilly, A. Nemes, C. Le Bas. 1999.

30、Development and use of a database of hydraulic properties of European soils. Geoderma, 90:169-185.16 Pachepsky, Y.A., D.J. Timlin, L.R. Ahuja. 1999. Estimating saturated soil hydraulic conductivity using water retention data and neural networks. Soil Science, 164(8): 552-560.17 Wosten J H M, Pacheps

31、ky Y A, Rawls W J. 2001. Pedo-transfer functions: bridging the gap between available basic soil data and missing soil hydraulic characteristics. Journal of Hydrology . 251: 123150.Indirect methods to estimate saturated soil hydraulic conductivity in paddy soil of ChangshuZHANG Jun-hua1,2, LIU Jian-li1*, ZHU An-ning1(1. Institute of soil science, Chinese academy of sciences, Nanjing 210008,China; 2. Graduate University of Chinese Academy of Sciences, Beijing 100049,China)Abstract: As an alternative approach of direct experiment, using soil physical