史鑫蕊 梁 浩 周 豐 胡克林
(1.中國農(nóng)業(yè)大學土地科學與技術學院,北京100193;2.農(nóng)業(yè)農(nóng)村部華北耕地保育重點實驗室,北京100193;3.河海大學農(nóng)業(yè)工程學院,南京210098;4.北京大學城市與環(huán)境學院,北京100871;5.教育部地表過程分析與模擬重點實驗室,北京100871)
土壤-作物系統(tǒng)模型通過計算機程序和數(shù)學方程描述作物生長過程及碳氮元素在土壤-作物-大氣連續(xù)體中的運移和轉化規(guī)律[1]。與田間試驗相比,不僅可以節(jié)省大量的人力物力,而且可以通過情景分析模擬不同的田間管理、氣候變化、土壤屬性以及作物品種等因子對作物生長、水分消耗以及碳氮運移和轉化過程的影響,在水肥優(yōu)化管理、品種改良以及應對氣候變化等方面發(fā)揮了重要作用[2]。目前比較成熟的土壤-作物系統(tǒng)模型主要有EPIC、WOFOST、DNDC、RZWQM、APSIM 和DSSAT 等,在小麥、玉米、水稻等作物上已進行了廣泛的驗證和應用,但是這些模型應用于不同地區(qū)時均需要對模型參數(shù)進行校準和驗證[3]。由于輸入?yún)?shù)眾多,模型參數(shù)率定過程困難。另外,參數(shù)校驗過程中會產(chǎn)生“異參同效”現(xiàn)象,大大降低了模型校驗的可靠性[4-5]。敏感度分析能夠從眾多輸入?yún)?shù)中識別出敏感度高的參數(shù),通過簡化或固定低敏感參數(shù),可以減小模型校驗的工作量、降低不確定性[6-8]。
模型敏感度分析主要有局部和全局敏感度分析兩種方法。局部敏感度法每次只改變一個參數(shù),評價單個參數(shù)變化對模擬結果的響應,忽略了參數(shù)間的交互作用對輸出結果的影響,僅適用于過程簡單的線性模型[9]。全局敏感度法綜合考慮各個參數(shù)及參數(shù)間的交互作用對輸出結果的影響,更適合于分析具有非線性的土壤-作物模型參數(shù)的敏感度[10]。常用的全局敏感度法主要有Morris、EFAST和Sobol'法等。DEJONGE 等[11]用Morris 和Sobol'兩種方法研究了不同灌溉條件下CERES-Maize 模型的參數(shù)敏感性,兩種方法得出的結果具有很強的相關性。宋明丹等[1]比較了Morris 和EFAST 兩種方法得出的CERES-Wheat 模型作物參數(shù)的敏感度,兩者結果基本一致,但Morris 方法計算工作量明顯較小。EFAST 和Sobol'均是基于方差分析的敏感度法,精度高且穩(wěn)定性強,但運算量大[12]。Morris法可以快速地從大量的模型輸入?yún)?shù)中識別出部分敏感度較高的參數(shù),運算量較小,已被廣泛應用于水文、生態(tài)和環(huán)境類模型的敏感度分析中[13-15]。俞雙恩等[16]采用Morris 法分析了DRAINMOD-S 模型6 個參數(shù)對土壤剖面含鹽量的的靈敏度。高穎會等[17]采用Morris 方法篩選出對SWMM 模型區(qū)域洪峰流量和徑流系數(shù)比較敏感的水文水力參數(shù)。JABLOUN 等[18]通過Morris 方法得到了Daisy 模型在模擬小麥-玉米輪作系統(tǒng)中作物產(chǎn)量和氮素淋洗的敏感參數(shù)。前人對各種模型的參數(shù)敏感度進行了大量研究,但大部分研究主要集中在旱地作物系統(tǒng)上,對于水田系統(tǒng)的研究還較少。
LIANG 等[19]在參考和借鑒國外主流模型的基礎上開發(fā)了適用于水稻的WHCNS_Rice(Soil water heat carbon nitrogen simulator for rice)模型,該模型考慮了中國的氣候、土壤、環(huán)境和農(nóng)業(yè)管理特點,相比國外模型更適用于中國高強度集約化的管理模式,已成功用于模擬長江中游地區(qū)常規(guī)淹水和覆膜旱作條件下的水稻生長及稻田水氮遷移過程[20]。水田與旱地的環(huán)境條件及作物生長過程機理不同,水田考慮的過程更為復雜,輸入?yún)?shù)更多,而且具有很強的非線性,這就導致模型參數(shù)校準過程更為困難,大大限制了水田模型的推廣應用。此外,水田與旱地系統(tǒng)模型參數(shù)敏感性的異同點尚不明確。因此,本文以WHCNS_Rice 模型為例,應用長江中游地區(qū)2 年的水稻田間試驗數(shù)據(jù),基于Morris 和Sobol'兩種方法對該模型進行參數(shù)敏感度分析,旨在為提高水田系統(tǒng)模型的校準效率及其應用提供技術指導。
試驗地位于長江中游地區(qū)的湖北省荊州農(nóng)業(yè)氣象試驗站(30°21'N,112°9'E),種植制度為水稻-油菜輪作。研究區(qū)屬濕潤季風氣候,年平均氣溫16℃,年均降雨量1 095 mm。于2017 年和2018 年水稻生長季進行了田間試驗,對水稻生長、水分消耗及氮素去向等指標進行了動態(tài)觀測。水稻品種為“豐兩優(yōu)香一號”,2 年移栽日期分別為6 月6 日和6 月9 日,收獲日期分別為9 月14 日和9 月18 日,種植密度為9 萬株/hm2。小區(qū)面積150 m2(25 m×6 m),3次重復。根據(jù)當?shù)剞r(nóng)民管理習慣,施用純氮230.5 kg/hm2,純磷26.6 kg/hm2,純鉀40.4 kg/hm2。氮肥分3 次施用,分施比例為5∶3∶1,其中基肥和第2 次追肥施用復合肥,第1 次追肥施用尿素。水稻生長季內(nèi)水分管理采取淹水-曬田-再次淹水-乳熟期后自然落干至收獲[21]。
試驗前開挖土壤剖面,取樣測定0 ~90 cm 土壤剖面各層的基本理化性質。試驗地土壤類型為水耕人為土,pH 值8.1,土壤有機碳質量比7.9 g/kg,全氮質量比1.2 g/kg。具體的土壤物理與水力學性質見表1。水稻生育期內(nèi)白天每2 h 測定一次田面水高度。通過流量計記錄每次發(fā)生灌水和徑流時的水量,取水樣采用流動分析儀(AA3 型,Bran +Luebbe,德國)測定銨態(tài)氮和硝態(tài)氮含量。2017 年采用動態(tài)箱法測定氨揮發(fā)量,水稻關鍵生長期取植株樣干燥測定秸稈和籽粒干物質質量,收獲時測產(chǎn),同時用碳氮元素分析儀(Flash2000 型,Thermo,美國)測定植株全氮含量。其他的測定指標及方法見文獻[21]。
表1 土壤剖面基本物理與水力學性質Tab.1 Soil basic physical and hydraulic properties in soil profile
WHCNS_Rice 模型以天為步長,由氣象數(shù)據(jù)和作物生物學參數(shù)驅動,可用于模擬水稻生長、稻田水分運動、氮素運移和轉化以及熱傳導等過程。模型采用Penman-Monteith 公式估算參考作物蒸散量。采用Green-Ampt 模型和Richards 方程分別模擬非飽和條件下土壤水分入滲和再分布過程,飽和條件下采用達西定律。每天的田面水高度通過水分平衡方程計算,公式為
式中 Pn、Pn-1——第n 天、第n-1 天田面水高度,mm I
n——第n 天灌溉量,mm
P′n——第n 天降雨量,mm
ETn——第n 天蒸散量,mm
Rn——第n 天徑流量,mm
Fn——第n 天水分入滲通量,mm
稻田徑流的計算與旱地不同,稻田通常有一定高度的田埂,將此高度作為最大田面水允許高度。在強降雨條件下,需先將水稻田池蓄滿,超過田埂之后才會產(chǎn)生徑流。
土壤有機質周轉以及碳氮循環(huán)直接參考Daisy模型。土壤氮素運移過程使用對流-彌散方程描述,源匯項考慮了有機質礦化、生物固持、尿素水解、氨揮發(fā)、硝化、反硝化和作物吸收等過程。采用改進的PS123 作物模型對作物生長發(fā)育進程、干物質生產(chǎn)和分配及作物產(chǎn)量形成過程進行模擬。詳細的模型原理見文獻[19 -20]。
1.3.1 Morris 方法
Morris 是一種基于篩選分析的全局敏感度分析方法,利用OAT(One factor at a time)搜索途徑,通過計算輸入?yún)?shù)對輸出結果的基礎效應進行靈敏度分析。該方法有效平衡了計算的效率和準確性,對于輸入?yún)?shù)多且運算負荷大的模型具有很好的適用性。每個輸入?yún)?shù)對輸出結果的影響程度Ei[22]計算公式為
其中
式中 y——參數(shù)i 變化前的輸出值
y*——參數(shù)i 變化后的輸出值
Δi——參數(shù)i 的變化量Δ 介于1/(p -1)和1 -1/(p -1)之間,p 表示參數(shù)的水平。
由于Morris 方法參數(shù)取值的隨機性,易造成取值誤差,所以需進行r 次重復,模型總運行次數(shù)為r(L+1)次,L 表示參數(shù)個數(shù)。根據(jù)Morris 方法建議,采樣因子r 設置為10,輸入?yún)?shù)個數(shù)L 為28,共運行模型290 次。最后計算每個參數(shù)對模擬值的影響程度Ei的平均值μ 和標準差σ。μ 越大,說明參數(shù)對輸出結果越敏感,而σ 表示的是參數(shù)間的交互作用,σ 越大表示該參數(shù)與其他參數(shù)交互作用越強。
1.3.2 Sobol'方法
Sobol'是一種基于方差的全局敏感性分析法,通過分解輸出變量的方差來定量化評價各輸入?yún)?shù)和參數(shù)間交互作用對輸出變量的影響程度[23]。假設y 為模型輸出變量;Xi為模型輸入?yún)?shù)(共m個)。y 的方差D(y)可分解為
式中 Di——參數(shù)i 產(chǎn)生的方差
Dij——參數(shù)i 和j 相互作用產(chǎn)生的方差
Dijk——參數(shù)i、j 和k 相互作用產(chǎn)生的方差
D1,2,…,m——m 個參數(shù)共同作用產(chǎn)生的方差
對于參數(shù)i,可用一階敏感度指數(shù)Si反映單個參數(shù)的敏感性,全階敏感度指數(shù)STi表示單個參數(shù)與其他所有參數(shù)的共同影響,具體公式為
式中 D~i——除了參數(shù)i 之外的參數(shù)的方差
SALTELLI 等[23]建議,Sobol'方法取樣個數(shù)設置為n(m+2),其中m 為參數(shù)個數(shù)(28 個),n 取值應大于100。本研究計算了n 為1 024 時的Sobol'敏感度結果,需運行模型29 696 次。
1.3.3 參數(shù)選擇、分布和范圍
WHCNS_Rice 模型輸入?yún)?shù)主要分為土壤水力學參數(shù)、氮素轉化參數(shù)和作物參數(shù)3 類。首先根據(jù)實測的土壤理化性質數(shù)據(jù),基于土壤傳遞函數(shù)獲得水力學參數(shù)初始值[24],作物參數(shù)和氮素轉化參數(shù)初始值參考相關文獻數(shù)據(jù)或者直接采用默認值。結合文獻中其他模型敏感度分析結果[10],選擇28 個參數(shù)進行敏感性分析(表2),其中水力學參數(shù)每層5個共15 個、氮素轉化參數(shù)5 個、作物參數(shù)8 個。由于各參數(shù)的分布概率未知,假定所有的參數(shù)服從均勻分布。模型參數(shù)取值范圍采用初始值的±10%作為其區(qū)間的上限和下限,采用蒙特卡洛方法對每個參數(shù)隨機采樣。輸入?yún)?shù)組合的產(chǎn)生以及最后的結果分析均采用SimLab 2010 軟件完成[25]。
表2 模型輸入?yún)?shù)Tab.2 Model input parameters
根據(jù)模型敏感度分析結果,本文基于2017 年的大田實測數(shù)據(jù),采用試錯法對WHCNS_Rice 模型中的敏感參數(shù)進行率定。首先對土壤水力學參數(shù)進行調整,使模型模擬的田面水高度、蒸散量、徑流量等與實測值基本吻合。然后率定氮素轉化參數(shù),使氨揮發(fā)量、反硝化量、氮淋洗量、徑流損失量等模擬值與實測值吻合,若缺少實測值應使模擬值接近文獻參考值。最后,調整作物參數(shù),使產(chǎn)量、干物質質量、吸氮量等指標的模擬值與實測值盡量吻合。由于參數(shù)間存在交互作用,不同參數(shù)調整過程會相互影響,因此該過程需要周而復始,最終大部分模擬值與實測值誤差最小時模型率定過程結束。模型率定結束后,固定所有輸入?yún)?shù),根據(jù)2018 年實際田間管理和氣象數(shù)據(jù)運行模型,用實測的田面水高度、蒸散量、作物吸氮量和干物質質量等來驗證模型。若模型模擬值與實測值誤差在可接受范圍內(nèi),說明模型輸入?yún)?shù)合理,可用來模擬該地區(qū)水稻生長及稻田水分運動和氮素轉化過程。
2.1.1 水分輸出項對不同參數(shù)的敏感度
根據(jù)模型輸入?yún)?shù)初始值,各參數(shù)設定±10%的變化范圍,基于Morris 敏感度分析方法,通過SimLab 軟件共產(chǎn)生了290 組模型組合參數(shù)。不同參數(shù)組合下分別運行WHCNS_Rice 模型,得到了290 組不同的輸出結果,最后采用SimLab 軟件進行參數(shù)敏感度分析。
圖1 為模型模擬的蒸散量(ET)、水分滲漏量、徑流量和田面水高度對各輸入?yún)?shù)的敏感度結果。ET 對中期作物系數(shù)Kmid最敏感,其次為Kend和Kini,但是σ 較低,說明與其他參數(shù)交互作用較弱(圖1a)。Tsum是第4 敏感參數(shù),通過影響作物生育期天數(shù)間接影響ET,與其他參數(shù)有較強的交互作用。其余參數(shù)μ 均小于10 mm,但是F、Smax、θs等參數(shù)的σ 均較高,說明部分土壤水力學參數(shù)和作物參數(shù)有較強的交互作用。
圖1 ET、滲漏量、田面水高度及徑流量對WHCNS_Rice 輸入?yún)?shù)的Morris 敏感度均值和標準差Fig.1 Average Morris mean effect and spread of WHCNS_Rice parameter with respect to ET,drainage,ponding water depth and runoff
對于田面水高度、滲漏量和徑流量,敏感度最高的參數(shù)均為θs1,其μ 和σ 均較高(圖1)。其次為F1,具有較高的μ,但σ 相對θs1偏低。第3 敏感參數(shù)是Ks2,其μ 較小,但σ 較高。表層土壤飽和含水率和田間持水率很大程度上決定了土壤的蓄水能力,Ks2為犁底層土壤飽和導水率,決定了土壤水分的入滲速率。土壤蓄水率和導水率共同決定了土壤中水分的運動,從而對田面水高度、滲漏量和徑流量產(chǎn)生顯著影響?;谒制胶庠恚珽T、土壤儲水量、田面水高度、滲漏量和徑流量,這些水分平衡項之間存在相互影響、此消彼長的關系。土體儲水量的增加必然導致滲漏量降低,ET 會影響田面水高度的變化,進而對徑流量產(chǎn)生影響。
2.1.2 氮素損失項對不同參數(shù)的敏感度
WHCNS_Rice 模型稻田氮素損失項主要包括氮淋洗、反硝化、氨揮發(fā)和徑流損失等4 個氮素損失項。從模型敏感度分析結果來看,不同土層飽和含水率和田間持水率均對氮淋洗有較大的影響,氮淋洗對θs1和F1最敏感(圖2a)。其余敏感度較高的幾個參數(shù)按μ 值由大到小依次為θs2、θs3、F3、F2、W2,可以看出氮淋洗對各個土層的飽和含水率及田間持水率均比較敏感,這些參數(shù)值越大表示土壤蓄水能力越強,可減少水分滲漏,從而減小氮素隨水分下移造成的淋洗損失。其余參數(shù)敏感度均低于5 kg/hm2,對氮淋洗的影響較小。
圖2 氮淋洗、反硝化、徑流和氨揮發(fā)對WHCNS_Rice 輸入?yún)?shù)的Morris 敏感度均值和標準差Fig.2 Average Morris mean effect and spread of WHCNS_Rice parameter with respect to nitrogen leaching,denitrification,runoff and ammonia volatilization
對反硝化作用影響最大的參數(shù)是F1,但其σ 較低(圖2b)。其次是θs1和θs2,且具有較高的σ。反硝化作用是在厭氧條件下發(fā)生的,水力學參數(shù)通過影響土壤中的水分狀況進而對反硝化產(chǎn)生影響。反硝化經(jīng)驗系數(shù)Ad是第4 敏感參數(shù),與其他參數(shù)交互作用較弱,Ad越大,通過反硝化損失的氮素越多,對反硝化輸出項有最直接的影響。
氮徑流損失對輸入?yún)?shù)的敏感性分析結果與氮淋洗相似,θs1和F1敏感度最大,其次是θs3和θs2(圖2c)。高蓄水能力下可以減少降雨和灌溉時產(chǎn)生的徑流損失。下層土壤的飽和含水率、田間持水率、飽和導水率等參數(shù)敏感度均值μ 相對較低,但非線性效應較強,主要通過影響水分的下滲速率對徑流產(chǎn)生間接作用。
氨揮發(fā)對θs1和θs2最敏感,且具有較強的非線性效應,上層土壤的飽和含水率通過影響土壤中的銨態(tài)氮濃度間接影響氨揮發(fā)(圖2d)。對氨揮發(fā)影響較大的第3 個敏感參數(shù)是Kv,Kv越大,氨揮發(fā)速率越快。Smax是第4 敏感參數(shù),Smax決定了覆蓋度,覆蓋度影響了NH3從水面向大氣中的運動過程。不僅如此,覆蓋度是影響土面或水面溫度的重要原因,而氨揮發(fā)與溫度密切相關,高溫會加快氨逸出。
2.1.3 作物生長指標對不同參數(shù)的敏感度
從干物質質量、產(chǎn)量、吸氮量和葉面積指數(shù)(LAI)對模型輸入?yún)?shù)的敏感度結果來看,這些作物生長指標均對作物參數(shù)Tsum和Smax具有較高的敏感度(圖3)。Tsum控制了作物生長的物候期,Tsum變大作物生育期延長,反之生育期變短,Smax決定了作物冠層覆蓋度,因此這兩個作物參數(shù)對作物生長及產(chǎn)量形成均有較大影響(圖3a、3b)。土壤水力學參數(shù)θs1和θs2對干物質質量、產(chǎn)量和吸氮量影響較大,并與其他參數(shù)有較強的交互作用。由于水稻土一般具有犁底層,決定了水分入滲和養(yǎng)分淋失速率,因此Ks2對水稻生長也有一定的影響。土壤水力學參數(shù)通過影響稻田土壤中的水分狀況及氮素供應水平間接影響作物生長過程。
圖3 干物質質量、產(chǎn)量、吸氮量和LAI 對WHCNS_Rice 輸入?yún)?shù)的Morris 敏感度均值和標準差Fig.3 Average Morris mean effect and spread of WHCNS_Rice parameter with respect to dry matter,yield,crop uptake and LAI
對LAI 影響較大的3 個參數(shù)依次是Smax、Tsum、Smin,且σ 較高(圖3d)。Smax和Smin是與LAI 直接相關的參數(shù),直接影響LAI。Tsum通過控制生育期的長短對LAI 產(chǎn)生作用,三者共同決定了LAI。LAI 對Kmid的敏感度較低,但具有較高的σ,與其他參數(shù)交互作用較強。對干物質質量影響明顯的4 個參數(shù)依次為Tsum、Smax、θs1、θs2,其敏感度均值μ 均大于1 000 kg/hm2,θs2、Tsum和θs1的σ 較大。與前4 個參數(shù)相比,F(xiàn)1、Ks2、F3、Smin對干物質質量的影響相對較小,其μ 介于300 ~500 kg/hm2,干物質質量對其余參數(shù)不敏感。
對作物產(chǎn)量影響較大的參數(shù)依次為θs1、θs2、Tsum、F1、Smax、Ks2、F3(圖3b)。產(chǎn)量對θs1和θs2敏感度最高,且交互作用強。從前面的結果可以看出,θs1和θs2對各項水分消耗及氮素損失均有較大的影響,進而會對作物生長期間的水分和氮素供應狀況產(chǎn)生作用,最終影響產(chǎn)量形成。作物吸氮量的敏感參數(shù)與產(chǎn)量基本相同,依次是θs1、Tsum、θs2、Smax、F1、Ks2、F3、W2。
Sobol'方法共產(chǎn)生了29 696 組模型參數(shù)組合,分別運行WHCNS_Rice 模型輸出相應的結果,然后用SimLab 2010 軟件進行全局敏感度分析。圖4 為28 個輸入?yún)?shù)對不同作物生長、水分及氮素平衡輸出項的一階及全階敏感度指數(shù)。一階敏感度指數(shù)劃定超過閾值0.1 的為敏感參數(shù)??梢钥闯觯州敵鲰椫姓羯⒘康拿舾袇?shù)分別為Kmid和θr1,滲漏量、徑流量及田面水高度的敏感參數(shù)均為θs1和F1。氮素損失項中,氮淋洗對F1敏感,氮徑流和氮反硝化對θs1和F1敏感,氨揮發(fā)對θr1和Ks1敏感。作物生長指標中,干物質質量、吸氮量和LAI 均對Tsum和Smax敏感,而產(chǎn)量只對Tsum敏感,Smax對其影響不大。整體上Sobol'一階敏感度指數(shù)得出的敏感參數(shù)與Morris 方法基本一致,但是氮素平衡項中氮淋洗和氨揮發(fā)無法獲得對其影響較大的參數(shù),因此改用手動調參,發(fā)現(xiàn)參數(shù)Kv和Ad分別對氨揮發(fā)和氮反硝化輸出結果有一定的影響,由此可見,Morris 得到的結果準確可靠。
從全階敏感度指數(shù)看,大部分輸出項很難篩選出明顯的敏感參數(shù)。不同參數(shù)對蒸散量、滲漏量、徑流量和田面水高度的全階敏感度指數(shù)分別介于0.81 ~1.00、0.55 ~0.89、0.27 ~0.71、0.44 ~0.85之間,差異均不明顯。氮素損失項中,θs1和F1對氮徑流和氮反硝化的全階敏感度指數(shù)明顯高于其余參數(shù),而各參數(shù)對氨揮發(fā)和氮淋洗的全階敏感度指數(shù)差異不明顯。但作物生長輸出項的敏感參數(shù)比較明顯,除Tsum和Smax,θs1和θs2也有較高的敏感度,其次為F1和F2,這也與Morris 的結果一致(圖3)。
某一參數(shù)全階敏感度指數(shù)與一階敏感度指數(shù)的差值即為該參數(shù)與其余參數(shù)交互作用對模型輸出項的影響程度。差值越大,表明該參數(shù)交互作用越強,這與Morris 方法中的σ 相對應,解釋了作物生長指標中部分土壤水力學參數(shù)敏感度較高的原因(圖3)。
圖4 WHCNS_Rice 模型參數(shù)對不同輸出項的一階敏感度指數(shù)及全階敏感度指數(shù)Fig.4 First-order and total sensitivity of WHCNS_Rice model parameters for different objective functions
與Morris 方法相比,Sobol'方法較好地考慮了參數(shù)間的交互作用,這是導致兩種方法結果略有不同的主要原因,但是Sobol'分析需要大批量運算才能得出有效的結果。另外,本文同時也計算了參數(shù)組合分別為928 和3 712 次時的敏感度結果,其敏感度指數(shù)出現(xiàn)大量的負值,無法獲得準確的敏感參數(shù)。此外,Sobol'產(chǎn)生的大量參數(shù)組合中個別參數(shù)在模型計算中無法輸出結果,導致結果為空。這種情況下,采用與之相近的上一組模型結果作為輸出值,這也造成了Sobol'分析結果的誤差。因此,綜合考慮計算量、可靠性及精度,Morris 不僅可以有效篩選出敏感參數(shù),而且計算量較小,是適合于稻田土壤-作物系統(tǒng)WHCNS_Rice 模型參數(shù)敏感度分析的有效方法。
綜上模型敏感度分析結果,輸入?yún)?shù)的變化會對模擬結果造成很大的不確定性,合理的參數(shù)值對模擬結果的準確度十分重要。在模型調參過程中,應重點關注土壤水力學參數(shù)中的θs、F 和Ks,這幾個參數(shù)對水分運動、氮素運移和轉化以及作物生長均有顯著影響。其次是作物參數(shù),應重點關注Tsum、Smax和作物系數(shù)K。最后是氮素轉化參數(shù)中的Ad和Kv。以上參數(shù)是對模型敏感度較高的參數(shù),由于參數(shù)存在很強的非線性效應,且參數(shù)間存在交互作用,需要反復調整才能達到比較理想的效果。調參過程中對于敏感度較低的參數(shù),比如Rmax等,應根據(jù)實際情況進行設置,會使模型模擬效果更好。
基于敏感度分析結果,使用2017 年實測數(shù)據(jù)校準模型。固定非敏感參數(shù),通過“試錯法”重點校準高敏感參數(shù),使模型模擬的干物質質量、吸氮量、田面水高度和蒸散量與實測值盡量吻合,最后采用2018 年實測數(shù)據(jù)驗證模型,校驗后的參數(shù)值見表2。
從圖5 可以看出,模型校準年份(2017 年)的干物質質量、作物吸氮量、田面水高度和ET 的模擬值與實測值的相關系數(shù)r 分別為0.998、0.979、0.502和0.852,與2018 年模型驗證相對應的相關系數(shù)分別為0.998、0.985、0.879 和0.856。其中2017 年田面水高度模擬值與實測值的相關系數(shù)較低,主要是因為移栽后15 ~27 d 田面水高度模擬值和實測值相差較大所致。由于該時段并未發(fā)生大的降雨,只有第22、23 天有少量降雨,故有可能是人為觀測誤差造成的,但此誤差在可接受范圍內(nèi)。校準(驗證)的線性回歸方程的斜率分別為1.024 3、1.072 1、0.938 3 和0.975 0(1.000 1、1.120 5、1.099 9 和1.086 5),均接近1,表明模擬值與實測值顯著相關(P <0.01)(圖5)??傮w來說,模型模擬的各項指標與實測值之間具有很好的一致性,校驗后的WHCNS_Rice 模型可用于模擬該地區(qū)水稻生長過程、水分消耗以及稻田氮素的遷移和轉化過程。
圖5 校準和驗證年份模擬和實測的干物質質量、作物吸氮量、田面水高度和ET 的相關性Fig.5 Measured/simulated values with linear regression line for dry matter,crop nitrogen uptake,ponding water depth and ET in 2017 and 2018
土壤-作物系統(tǒng)模型是實現(xiàn)數(shù)字化農(nóng)田作物生長監(jiān)測、預測和決策支持的重要手段[12]。模型的參數(shù)校準是模型應用的重要前提和關鍵步驟[26],也是模型功能拓展的基礎,同時輸入?yún)?shù)的合理性和區(qū)域代表性決定了模擬結果的準確性。由于模型參數(shù)較多,在缺乏先驗知識的情況下,模型調參過程非常困難,限制了模型的推廣應用。敏感性分析可以從大量的參數(shù)中識別出敏感參數(shù),減少模型校驗的工作量和不確定性,進而提高模型校準的效率和精度。
土壤-作物模型輸入?yún)?shù)大致可以分為3 大類:土壤水力學參數(shù)、作物參數(shù)和碳氮轉化參數(shù),不同類型的參數(shù)對模型輸出結果的影響有明顯差異。梁浩等[10]對WHCNS 模型小麥和玉米作物的參數(shù)敏感度進行了分析,結果表明土壤含水率、土壤硝態(tài)氮含量、作物產(chǎn)量和LAI 的綜合響應對作物參數(shù)最敏感,其次是土壤水力學參數(shù),而氮素轉化參數(shù)敏感度較低,這與本研究參數(shù)敏感度結果略有不同。本研究中除LAI 和干物質質量外,表現(xiàn)為土壤水力學參數(shù)敏感度高于作物參數(shù),尤其是水氮各輸出項。這可能是由水田和旱地土壤性質及環(huán)境條件的顯著差異造成的。水田與旱地土壤中的水氮運移有很大不同,土壤水力學參數(shù)是影響土壤中水分運動和氮素遷移的主要參數(shù)[27]。淹水稻田生育期內(nèi)很少出現(xiàn)水分脅迫情況,但是由于生育期內(nèi)頻繁降雨和大量灌溉,很容易發(fā)生大量的徑流和水分滲漏損失,其損失量就取決于土壤的水力學性質,尤其是犁底層飽和導水率,決定了稻田水分下滲速率。不僅如此,水分運動直接影響了土壤中氮素的遷移與轉化過程,從而會影響土壤的供氮水平。DEJONGE 等[11]的研究結果表明,在充分灌溉條件下,CERES-Maize 模型模擬的產(chǎn)量、蒸散量和LAI 對作物參數(shù)最敏感,而虧缺灌溉下田間持水率和萎蔫含水率是主要的敏感參數(shù),由此可見,土壤中的水分狀況是影響參數(shù)敏感性的重要因素。
關于作物參數(shù)對不同模型輸出項的敏感度研究已經(jīng)有大量報道。ZHAO 等[28]研究了APSIM 模型輸出的產(chǎn)量、生物量、開花期和成熟期對作物參數(shù)的敏感度,發(fā)現(xiàn)產(chǎn)量和生物量均對積溫比較敏感。VANUYTRECHT 等[6]采用AquaCrop 模型模擬分析了4 個國家不同作物產(chǎn)量對參數(shù)的敏感度,結果表明水稻產(chǎn)量對積溫、最大覆蓋度敏感度較高。本研究中蒸散量、干物質質量、產(chǎn)量、吸氮量及LAI 對作物參數(shù)敏感度較高,特別是Tsum和Smax兩個參數(shù),這與上述研究結果一致。除Tsum和Smax外,本研究作物參數(shù)中作物系數(shù)是第3 敏感參數(shù),尤其是對ET有較大影響。由于WHCNS_Rice 模型中潛在作物蒸散量(ETp)是由參考作物蒸散量(ET0)與作物系數(shù)乘積得到的[29],因此不同時期作物系數(shù)的大小會對ET 產(chǎn)生直接影響。這與LIANG 等[30]的研究結果不太一致,其結果中對春玉米產(chǎn)量影響較大的第3 個參數(shù)是最大同化率Amax,而本研究中水稻生長對Amax不敏感,這主要是由于C3 和C4 作物光合作用機理不同造成的,在模型中采用了不同的光合作用公式。
相對于土壤水力學參數(shù)和作物參數(shù),模型輸出項對氮素轉化參數(shù)的敏感度相對較低,只有氨揮發(fā)量和反硝化量分別對Kv和Ad敏感。而梁浩等[10]的研究結果表現(xiàn)為最大硝化速率Vn和硝化半飽和系數(shù)Kn敏感度較高,這與本研究的結果不一致。主要是由于水田和旱地的差異造成的,大量研究表明旱地硝化作用比較強烈,土壤硝態(tài)氮含量高,而銨態(tài)氮含量較低,因此與硝化作用相關的參數(shù)敏感度高。而水田土壤中硝態(tài)氮含量明顯低于銨態(tài)氮[31],硝化能力較弱,氨揮發(fā)和氮反硝化作用比較強烈,導致兩者的敏感參數(shù)正好相反。水田與旱地由于土壤水分條件不同,造成了氮素遷移和轉化過程差異較大,相應的敏感參數(shù)也有較大差異,尤其是氮素轉化參數(shù),因此模型重點校準參數(shù)也應根據(jù)土壤水分狀況調整。
不同的作物品種、環(huán)境條件及田間管理均會對模型參數(shù)的敏感度產(chǎn)生較大影響。LIANG 等[30]研究表明不同水肥管理條件下模型參數(shù)敏感度有較大的變化。ASSENG 等[32]發(fā)現(xiàn)不同的氣候條件對產(chǎn)量模擬值造成了很大的不確定性。ZHAO 等[28]研究指出不同的施肥量比氣候條件和土壤類型的改變對模型參數(shù)敏感度的影響更大。本文僅研究了某一特定條件下的參數(shù)敏感度,沒有分析不同的土壤類型、氣候條件及水肥管理等差異造成的不確定性。為了擴大WHCNS_Rice 模型應用范圍,應進一步全面系統(tǒng)地考慮環(huán)境條件、作物品種、管理措施等不同情況下的參數(shù)敏感度及不確定性。
(1)采用Morris 和Sobol'兩種全局敏感度分析方法,分別針對作物生長、水分運動及氮素損失3 部分輸出項,具體分析了稻田土壤-作物系統(tǒng)WHCNS_Rice 模型輸出項對土壤水力學參數(shù)、氮素轉化參數(shù)和作物參數(shù)的敏感性。結果表明,整體上Sobol'一階敏感度指數(shù)得出的敏感參數(shù)與Morris 方法基本一致,但是綜合考慮計算量及精度,Morris 方法可以快速并有效地篩選出模型敏感性參數(shù)。模型各個輸出項均對土壤水力學參數(shù)有較高的敏感度,尤其是對飽和含水率、田間持水率和飽和導水率的敏感度較高。其次是作物參數(shù),其中作物生育期總有效積溫、最大比葉面積和作物系數(shù)對產(chǎn)量、干物質質量、葉面積指數(shù)和蒸散量有較大影響,對水分運動和氮素遷移過程影響較小。相反,作物生長過程各輸出項對土壤水力學參數(shù)和作物參數(shù)均有較高的敏感性。各輸出項對氮素轉化參數(shù)敏感度相對最低,只有氨揮發(fā)一階動力學系數(shù)和反硝化經(jīng)驗系數(shù)分別對氨揮發(fā)和反硝化過程有一定的影響。
(2)水田與旱地水分條件的不同,導致水田系統(tǒng)土壤水力學參數(shù)敏感度高于作物參數(shù),而旱地正好相反。另外,旱地系統(tǒng)由于硝化作用強烈,因此與硝化反應相關的參數(shù)(最大硝化速率Vn和硝化半飽和系數(shù)Kn)對氮素輸出項有較大影響,而在水田中與氨揮發(fā)和氮反硝化過程相關的參數(shù)對輸出結果影響較大。
(3)基于模型參數(shù)敏感度篩選結果,根據(jù)稻田實測數(shù)據(jù)對WHCNS_Rice 模型參數(shù)進行了校驗。結果表明,模型模擬值與實測值具有較好的一致性,該模型可用于模擬長江中游地區(qū)水稻生長及稻田水氮遷移過程。該方法大大提高了模型校準效率,為土壤-作物系統(tǒng)模型的推廣應用提供了技術支持。