崔 穎 藺宏宏 謝 云 劉素紅

北京師范大學(xué)地表過(guò)程與資源生態(tài)國(guó)家重點(diǎn)實(shí)驗(yàn)室, 北京 100875

水是作物生長(zhǎng)的重要因子, 提高作物水分利用效率, 既可降低農(nóng)作區(qū)水資源消耗, 彌補(bǔ)水資源短缺, 并可有效提高作物產(chǎn)量。作物生長(zhǎng)模型是研究水分利用效率和作物產(chǎn)量的重要工具, 能在單點(diǎn)尺度上定量描述給定作物在一定環(huán)境下整個(gè)生長(zhǎng)周期內(nèi)的生長(zhǎng)發(fā)育和產(chǎn)量形成過(guò)程[1-2], 進(jìn)而建立起作物生長(zhǎng)對(duì)土壤及氣象環(huán)境變化的響應(yīng)系統(tǒng)[3], 為精準(zhǔn)農(nóng)業(yè)生產(chǎn)提供切實(shí)可行的研究方法。目前應(yīng)用比較廣的作物模型主要有AquaCrop、WOFOST、CERES等[4], AquaCrop 雖然比后2 個(gè)模型開(kāi)發(fā)得晚, 但其核心部分的作物需水量計(jì)算方法早在20 世紀(jì)70 年代就已被廣泛應(yīng)用[5]。這些模型都是以日為步長(zhǎng)模擬作物的生長(zhǎng)過(guò)程與產(chǎn)量形成, 但模擬過(guò)程略有差異。AquaCrop 模型重點(diǎn)考慮了水分脅迫, 以水分利用效率為基礎(chǔ)模擬產(chǎn)量[6]。WOFOST 和CERES 模型都是先模擬環(huán)境適宜條件下的最大產(chǎn)量, 然后考慮熱量、水分和氮磷等肥力脅迫。不同的是WOFOST對(duì)呼吸作用考慮更加全面, 且與AquaCrop 模型類似,首先模擬生物量, 然后利用收獲指數(shù)與生物量相乘得到作物產(chǎn)量[7-8]。CERES 模型則是模擬作物的穗數(shù)、穗粒數(shù)和穗粒重后計(jì)算產(chǎn)量[9]。旱地農(nóng)業(yè)受水分脅迫為主, AquaCrop 模型強(qiáng)調(diào)水分脅迫機(jī)制、模擬過(guò)程相對(duì)簡(jiǎn)單、且所需輸入?yún)?shù)較少, 使其得到廣泛應(yīng)用[10]。

目前已有多位國(guó)內(nèi)外學(xué)者應(yīng)用AquaCrop 模型在不同地區(qū)進(jìn)行了作物生長(zhǎng)動(dòng)態(tài)模擬并開(kāi)展了模型適應(yīng)性評(píng)估研究: 劉興冉等[11]將AquaCrop 模型應(yīng)用于石家莊市欒城區(qū)夏玉米水分研究中, 證明AquaCrop 模型能夠較好地模擬夏玉米的產(chǎn)量、生物量和冠層發(fā)育過(guò)程以及表層土壤水儲(chǔ)量的動(dòng)態(tài)變化;劉琦等[12]在山西省壽陽(yáng)縣進(jìn)行了覆膜和露地春玉米種植對(duì)比試驗(yàn), 以驗(yàn)證AquaCrop 模型對(duì)土壤含水率、農(nóng)田蒸散和冠層覆蓋度的模擬精度, 認(rèn)為該模型能較好地模擬晉中地區(qū)旱作覆膜春玉米的耗水、生長(zhǎng)和產(chǎn)量形成過(guò)程; Iqbal 等[13]應(yīng)用AquaCrop 模型對(duì)石家莊市欒城區(qū)冬小麥的生長(zhǎng)過(guò)程進(jìn)行了模擬,分析了不同程度水分脅迫條件下的模擬精度, 發(fā)現(xiàn)在水分嚴(yán)重虧缺時(shí)模擬結(jié)果與實(shí)際產(chǎn)量出入較大,而當(dāng)參數(shù)率定后則能有效提高產(chǎn)量模擬精度?；诖? AquaCrop 模型在應(yīng)用前需對(duì)強(qiáng)敏感性參數(shù)進(jìn)行率定。模型參數(shù)敏感性分析是篩選強(qiáng)敏感性參數(shù)的主要方法, 其中OAT (one factor at a time)[14]方法是最快速有效的敏感性方法之一, 即每次只改變1 個(gè)參數(shù)并根據(jù)參數(shù)改變量與模擬結(jié)果變化量的相關(guān)程度進(jìn)行篩選。已有部分研究對(duì)其他作物生長(zhǎng)模型參數(shù)敏感性進(jìn)行探究, 但關(guān)于AquaCrop 模型參數(shù)敏感性的研究仍有待進(jìn)一步補(bǔ)充和完善: 劉剛等[15]使用OAT 方法在山東省禹城市對(duì)ALMANAC 模型參數(shù)進(jìn)行了敏感性分析, 為后續(xù)應(yīng)用模型模擬黃淮海平原地區(qū)冬小麥和夏玉米的生長(zhǎng)提供了科學(xué)依據(jù); 宋明丹等[16]采用Morris 和EFAST 的敏感性分析方法,分析了陜西省楊凌區(qū)CERES-Wheat 模型中對(duì)開(kāi)花期、成熟期、產(chǎn)量、地上生物量4 個(gè)輸出變量的敏感性參數(shù), 發(fā)現(xiàn)2 種方法篩選出的強(qiáng)敏感性的作物參數(shù)基本一致, 其中Morris 法所需的計(jì)算量更小,操作性更強(qiáng), 而EFAST 法精確度更高, 可定量分析各個(gè)參數(shù)的敏感性, 在進(jìn)行參數(shù)的不確定性分析方面更有優(yōu)勢(shì)。可見(jiàn), 簡(jiǎn)單有效的敏感性分析方法即可實(shí)現(xiàn)強(qiáng)敏感性參數(shù)篩選的目的, 從而有針對(duì)性的進(jìn)行參數(shù)率定, 這對(duì)于提高AquaCrop 模型的模擬精度和區(qū)域適用能力而言十分必要。

東北黑土區(qū)是世界重要黑土帶之一, 土壤肥沃,氣候適宜, 主要出產(chǎn)玉米和大豆等作物, 是國(guó)家糧食安全戰(zhàn)略性保障基地[17-18]。特定的土壤和自然地理?xiàng)l件使得其雨養(yǎng)農(nóng)業(yè)生產(chǎn)模式有別于其他糧食主產(chǎn)區(qū), 大氣降水成為該地作物產(chǎn)量的主要影響因子之一[19]。為能應(yīng)用AquaCrop 模型模擬東北黑土區(qū)作物產(chǎn)量的水響應(yīng)機(jī)制, 彌補(bǔ)該模型在東北地區(qū)應(yīng)用研究的缺失, 本研究采用OAT 敏感性分析方法在分析模型輸入?yún)?shù)敏感性的基礎(chǔ)上, 基于作物試驗(yàn)觀測(cè)數(shù)據(jù)對(duì)敏感參數(shù)進(jìn)行率定, 然后利用大田觀測(cè)資料驗(yàn)證參數(shù)率定后的AquaCrop 模型在東北黑土區(qū)的適用性, 以便通過(guò)控制和減小環(huán)境變量效應(yīng),調(diào)整和優(yōu)化作物種植措施, 為該地區(qū)農(nóng)業(yè)高產(chǎn)可持續(xù)發(fā)展提供科學(xué)依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)位于黑龍江省嫩江縣鶴山農(nóng)場(chǎng), 地處小興安嶺西南麓、由小興安嶺向松嫩平原過(guò)渡的鶴北流域,地理位置(48°43′—49°03′N, 124°56′—126°21′E)。該區(qū)海拔310～390 m, 地形起伏較小, 有坡長(zhǎng)坡緩的特征[20], 坡度大多在3°～6°之間, 坡長(zhǎng)多為2000～4000 m。屬寒溫帶大陸性季風(fēng)氣候, 冬季寒冷干燥, 夏季溫?zé)岫嘤? 年均溫約0.4oC; 多年平均年降水量約534 mm, 6 月至8 月降水量約占全年的66.6%[21]。區(qū)內(nèi)土壤以黑土和草甸土為主, 屬于東北典型黑土區(qū), 因有機(jī)質(zhì)含量豐富而呈黑色。黑土土壤剖面表層為20～40 cm 的黑土層, 土壤質(zhì)地比較黏重, 土壤容重為1.39 g cm–3, 田間持水量為0.28 g g–1, 凋萎濕度為0.18 g g–1。向下為過(guò)渡層逐漸到達(dá)母質(zhì)層, 包括黃黏土、黃沙土或礫石等。研究區(qū)開(kāi)墾于20 世紀(jì)40 年代, 流域內(nèi)農(nóng)業(yè)用地占總土地面積的62.4%, 采用大型機(jī)械化作業(yè)和集中農(nóng)業(yè)管理方式, 主要農(nóng)作物為大豆和玉米每年輪作, 一般于每年5 月初播種至當(dāng)年9 月底收獲。玉米生長(zhǎng)期為150 d, 大豆生長(zhǎng)期為130 d。

1.2 AquaCrop 模型簡(jiǎn)介

AquaCrop 模型是由聯(lián)合國(guó)糧食及農(nóng)業(yè)組織(Food and Agriculture Organization of the United Nations, FAO)第33 號(hào)灌排文件Yield response to water[5]中的作物-水分響應(yīng)方程演變而來(lái)。該方程假設(shè), 作物產(chǎn)量(yield,Y) 的高低是對(duì)作物蒸散量(evapotranspiration,ET)大小的響應(yīng)結(jié)果:

式中,xY、aY分別為作物最大產(chǎn)量(kg m–2)和實(shí)際產(chǎn)量(kg m–2)。ExT 、EaT 分別為作物潛在蒸散量(mm)和實(shí)際蒸散量(mm)。潛在蒸散量指充分供水條件下的蒸散量, 首先根據(jù) FAO 第 56 號(hào)文件[22]中的Penman-Monteith 公式計(jì)算出參考蒸散量ET0[23], 然后再乘以作物系數(shù)得到。實(shí)際蒸散量是指實(shí)際供水情況下的蒸散量, 會(huì)受水分供給不足導(dǎo)致的土壤水分脅迫影響。ky為作物產(chǎn)量對(duì)土壤水分虧缺響應(yīng)的敏感系數(shù), 隨作物生育期變化。建立的AquaCrop 作物生長(zhǎng)模型以日為模擬步長(zhǎng), 通過(guò)輸入氣候、作物、土壤, 田間管理如灌溉、施肥及除草等數(shù)據(jù), 可以模擬作物的光合、呼吸、蒸騰等過(guò)程, 最終輸出冠層覆蓋度、生物量、產(chǎn)量等結(jié)果。

作物產(chǎn)量是作物成熟時(shí)的生物量(biomass, Bx)與收獲指數(shù)(harvest index, HI)相乘得到。作物生物量則是由歸一化的水分生產(chǎn)力即單位耗水量累積的生物量(normalized water productivity, WP*)乘以作物蒸騰量與參考蒸散量的比值得到:

式中,Bi為日地上生物量(t hm–2); Tri為日蒸騰量(mm), ETo,i為日參考蒸散量(mm)。WP*為歸一化水分生產(chǎn)力(gram m–2), 是將作物水分生產(chǎn)力除以標(biāo)準(zhǔn)蒸散條件下的作物蒸散量。實(shí)際計(jì)算時(shí), 可根據(jù)大氣實(shí)際CO2濃度、C3或C4作物類型、以及作物生長(zhǎng)期和產(chǎn)量成熟期進(jìn)行參數(shù)調(diào)整。

利用Penman-Monteith 公式計(jì)算出參考蒸散量后, 根據(jù)冠層覆蓋度最大時(shí)的作物蒸騰系數(shù)和實(shí)際生長(zhǎng)過(guò)程的冠層覆蓋度(canopy coverage, CC), 可以將作物生長(zhǎng)過(guò)程中的作物蒸騰量(transpiration,Tr)從蒸散量中分離出來(lái), 剔除了土壤蒸發(fā)的非生產(chǎn)性消耗性用水(evaporation, E)干擾。此外隨著冠層擴(kuò)張,因遮陰和冠層對(duì)空氣平流的影響都會(huì)影響作物蒸騰,需要對(duì)CC 進(jìn)行修訂, 用CC*表示:

式中, CC*為考慮了冠層遮陰和空氣流動(dòng)影響后的冠層覆蓋度(%); Kcbx為充分供水條件下參考作物蒸散的蒸騰系數(shù)。實(shí)際作物蒸騰受土壤水分脅迫影響,分別表現(xiàn)為缺水情況下的氣孔導(dǎo)度脅迫, 用脅迫系數(shù)Kssto表示, 以及多水情況下的土壤通氣性脅迫,用脅迫系數(shù)Ksaer表示。二者都是0～1 的無(wú)量綱參數(shù),表示脅迫程度, 值為1 時(shí)無(wú)脅迫。

冠層覆蓋度隨作物生長(zhǎng)發(fā)生變化, 直接影響作物的光合作用。模型采用指數(shù)形式表達(dá)CC 的增衰變化[24], 包括指數(shù)增長(zhǎng)期和平穩(wěn)增長(zhǎng)期(式4), 衰落期(式5)。

式中, CC 為冠層覆蓋度(%);t為從出苗開(kāi)始累積的時(shí)間; CC0為初始冠層覆蓋度(%), 一般取90%出苗時(shí)的平均幼苗覆蓋度; CCx為冠層覆蓋度達(dá)到最大時(shí)的值(%); CGC 為冠層覆蓋度增加速率(canopy growth coefficient), 表示單位生長(zhǎng)度日冠層覆蓋度的增加量(%, GDD–1); CDC 為冠層覆蓋度衰減速率(canopy decline coefficient), 表示單位生長(zhǎng)度日冠層覆蓋度的減少量(%, GDD–1)。如果受到土壤水分脅迫, 需要對(duì)CGC 和CDC 進(jìn)行修正, 修正公式如式(6)～(7)。

式中, CGCadj和CDCadj分別表示受水分脅迫影響的CGC 和CDC, Ksexp為水分對(duì)冠層覆蓋增加的脅迫系數(shù); Kssen為水分對(duì)冠層覆蓋減少的脅迫系數(shù)。此外, 收獲指數(shù)也會(huì)受到水分脅迫的影響, 缺水時(shí)利用Kssto調(diào)整; 滯水時(shí)利用調(diào)整。

植物水分依靠根系吸水, 因此模型考慮了根區(qū)中的根系生長(zhǎng), 及其對(duì)水分脅迫的響應(yīng)。給定土壤水分上限和下限, 當(dāng)根區(qū)水分消耗達(dá)到或超過(guò)上限時(shí),會(huì)影響冠層覆蓋增加, 但不會(huì)影響根系加深, 因?yàn)樗置{迫對(duì)冠層的影響強(qiáng)于對(duì)根系的影響。如果水分脅迫低于葉片氣孔關(guān)閉的閾值, 根系受到脅迫, 不會(huì)加深, 否則可以加深。根系增長(zhǎng)速率計(jì)算公式如下:

式中,Z為t時(shí)刻的有效根深(m);Zini為根系初始深度(m);Zx為最大有效根深(m);t0為90%出苗的時(shí)間;tx為根系達(dá)到最大有效根深的時(shí)間;t為根系生長(zhǎng)時(shí)間。

AquaCrop 模型涉及氣候、土壤、作物、管理等因素, 同時(shí)考慮了充分供水和水分脅迫對(duì)作物生長(zhǎng)發(fā)育和產(chǎn)量形成的影響。但其中的控制方程和參數(shù)都是基于試驗(yàn)觀測(cè)獲得的經(jīng)驗(yàn)方程和不同地區(qū)的參數(shù)?？紤]到不同地區(qū)之間氣候、土壤、作物品種、管理等的差異, 將模型應(yīng)用于某地區(qū)時(shí), 應(yīng)進(jìn)行參數(shù)的率定和模型的驗(yàn)證, 使其適用于當(dāng)?shù)貤l件。

1.3 數(shù)據(jù)來(lái)源

研究收集了所在地鶴山農(nóng)場(chǎng)氣象站觀測(cè)數(shù)據(jù),包括日最高/低氣溫(°C)、日降水量(mm)、太陽(yáng)輻射量(MJ m–2d–1)、水汽壓(kPa), 并利用FAO 發(fā)布的ET0Calculator 工具計(jì)算ET0(mm d–1), 建立氣象數(shù)據(jù)庫(kù)。通過(guò)測(cè)定各試驗(yàn)點(diǎn)土壤剖面的容重、田間持水量和凋萎濕度等土壤屬性數(shù)據(jù), 建立了東北黑土區(qū)土壤特性數(shù)據(jù)庫(kù), 部分參數(shù)如表1 所示。

表1 試驗(yàn)點(diǎn)的土壤物理參數(shù)表Table 1 Soil physical parameters of the experimental locations

參數(shù)率定所用數(shù)據(jù)來(lái)自2016—2017 年位于黑龍江省嫩江縣鶴山農(nóng)場(chǎng)的北京師范大學(xué)九三水土保持試驗(yàn)站作物試驗(yàn)。種植的玉米品種為德美亞1 號(hào), 大豆品種為1467。分別以鶴山農(nóng)場(chǎng)玉米實(shí)際平均種植密度9.25 萬(wàn)株 hm–2和當(dāng)?shù)卮筇锎蠖狗N植密度45 萬(wàn)株 hm–2為標(biāo)準(zhǔn)密度, 按標(biāo)準(zhǔn)密度的60%、80%、100%、120%、140%共5 個(gè)密度, 3 個(gè)重復(fù), 2 種作物各15 個(gè)小區(qū)種植, 單個(gè)小區(qū)面積為5 m × 8 m, 小區(qū)間隔0.5 m。試驗(yàn)播種時(shí)按當(dāng)?shù)胤N植方式每公頃施肥量為0.4 kg, 3 種底肥尿素、磷酸二銨和硫酸鉀的配比為1.0∶1.3∶0.3, 并依據(jù)管理標(biāo)準(zhǔn)按時(shí)中耕、除草以保證作物的正常生長(zhǎng)。每2d 觀測(cè)1 次作物生育及土壤水分動(dòng)態(tài), 記錄各個(gè)小區(qū)作物出苗期、開(kāi)花期、成熟期等關(guān)鍵生育期的出現(xiàn)日期。玉米成熟后(約150 d), 采生長(zhǎng)狀況接近平均水平的植株12 株風(fēng)干后將所有植株籽?；旌暇鶆蚝蟛蓸臃湃?5℃烘箱烘干至恒重后稱重。大豆成熟后(約130 d), 各小區(qū)采2 m × 2 m 的大豆樣本地上部分用恒溫鼓風(fēng)烘干箱105℃殺青2 h 后轉(zhuǎn)80℃烘干至恒重后稱重。

模型驗(yàn)證使用2011—2018 年鶴山農(nóng)場(chǎng)選取的16 塊農(nóng)田監(jiān)測(cè)點(diǎn)觀測(cè)的數(shù)據(jù)。在生長(zhǎng)季內(nèi)(5 月至9月)每10 d 觀測(cè)1 次監(jiān)測(cè)點(diǎn)的作物生長(zhǎng)動(dòng)態(tài), 觀測(cè)指標(biāo)包括株高、地上生物量和收獲指數(shù)等。對(duì)試驗(yàn)小區(qū)和農(nóng)田實(shí)測(cè)作物數(shù)據(jù)進(jìn)行整理與分析, 建立鶴北小流域作物特征參數(shù)數(shù)據(jù)庫(kù)。

1.4 模型參數(shù)敏感性分析與率定

AquaCrop 模型部分參數(shù)不隨地理位置、種植時(shí)間和管理措施變化而被稱為保守參數(shù), 無(wú)需調(diào)整[25]。為了認(rèn)識(shí)非保守參數(shù)中對(duì)模擬結(jié)果的影響程度, 采用敏感性分析的方法評(píng)估模型參數(shù)對(duì)模擬結(jié)果的影響程度[26]。本文采用OAT 法, 參考FAO 作物參數(shù)手冊(cè)[27]及實(shí)際作物特性確定擬進(jìn)行率定的模型參數(shù)(表2), 使用2011—2018 年大田數(shù)據(jù), 在AquaCrop模型中逐一改變某特征參數(shù)為標(biāo)注范圍極值或默認(rèn)值±10%及±20%的調(diào)整值, 保持其余變量為原數(shù)據(jù)集值不變, 對(duì)比作物產(chǎn)量的模型模擬結(jié)果, 分析產(chǎn)量對(duì)不同變量的敏感性。敏感性分析結(jié)果由相對(duì)敏感度(relative sensitivity, RS)[28]表示(式9), 其值經(jīng)過(guò)了標(biāo)準(zhǔn)化處理,可以進(jìn)行參數(shù)之間敏感性大小的比較[29]。

式中, RS 為相對(duì)敏感度,Y為模擬的作物產(chǎn)量,X為變量, ΔX為變量的調(diào)整值。按RS 值由大到小可將變量對(duì)產(chǎn)量影響的敏感程度劃分為: 極敏感(RS>0.6)、一般敏感(0.6>RS>0.1)和低敏感(RS<0.1)。

本研究采用2016—2017 年不同種植密度下的作物田間試驗(yàn)觀測(cè)數(shù)據(jù)率定模型參數(shù), 對(duì) RS 較低的參數(shù)使用默認(rèn)值, 僅對(duì) RS 值較高的參數(shù)進(jìn)行率定。

表2 AquaCrop 模型待敏感性分析參數(shù)介紹Table 2 Introduction of parameters to be sensitive analysis in AquaCrop model

1.5 模型精度評(píng)價(jià)

基于建立的氣候、作物、土壤數(shù)據(jù)庫(kù)和率定后的參數(shù), 采用2011—2018 年大田監(jiān)測(cè)點(diǎn)實(shí)測(cè)的18組玉米和16 組大豆產(chǎn)量對(duì)模型進(jìn)行驗(yàn)證。模型模擬精度評(píng)價(jià)采用決定系數(shù)(coefficient of determination,R2)、均方根誤差(root mean square error, RMSE)、標(biāo)準(zhǔn)均方根誤差(normalized root mean square error,NRMSE)和模擬效率(model efficiency, ME)作為精度評(píng)定指標(biāo), 評(píng)估預(yù)測(cè)結(jié)果與實(shí)測(cè)結(jié)果的差距, 計(jì)算公式如式(10)～(13)所示,R2和ME 越接近1、RMSE和NRMSE 越接近0 說(shuō)明預(yù)測(cè)精度越高。

式中,n為總觀測(cè)數(shù);Oi和Si分別為第i次觀測(cè)的實(shí)測(cè)值和預(yù)測(cè)值, t hm–2;和分別為平均實(shí)測(cè)值和平均預(yù)測(cè)值, t hm–2。

表3 AquaCrop 模型部分參數(shù)相對(duì)敏感度Table 3 Relative sensitivity of parameters in AquaCrop model

2 結(jié)果與分析

2.1 模型參數(shù)的相對(duì)敏感度

應(yīng)用AquaCrop 模型對(duì)東北黑土區(qū)主要作物生長(zhǎng)過(guò)程進(jìn)行模擬時(shí), 兩種作物的參數(shù)相對(duì)敏感度計(jì)算結(jié)果如表3 所示(相對(duì)敏感度為0 的參數(shù)表中已略)。對(duì)不同作物而言, 同一參數(shù)對(duì)產(chǎn)量的敏感性有較大差異。玉米產(chǎn)量對(duì)以下參數(shù)表現(xiàn)為極強(qiáng)的敏感性: 冠層形成和枯萎前的作物系數(shù)(crop coefficient before canopy formation and senescence, KcTr,x)、WP*、參考收獲指數(shù)(reference harvest index, HI0)、CDC、最大有效根深(maximum effective rooting depth,Zx)。大豆產(chǎn)量的極敏感性參數(shù)較少, 包括限制冠層伸展的水分脅迫系數(shù)曲線的形狀因子(shape factor for water stress coefficient for canopy expansion,Pexshp)、KcTr,x、根區(qū)根系伸展曲線的形狀因子(shape factor describing root zone expansion, Rexshp)。部分參數(shù)的相對(duì)敏感度對(duì)兩種作物而言基本一致, 如限制冠層伸展的土壤水分消耗上限閾值(soil water depletion threshold for canopy expansion-upper threshold, Pexp-up)、因衰老、氮元素虧缺導(dǎo)致的作物系數(shù)下降速率(decline of crop coefficient as a result of ageing, nitrogen deficiency, DeKcTr,x)、產(chǎn)量形成期的歸一化水分生產(chǎn)力(water productivity normalized for ET0and CO2during yield formation, WP*yf)、CCx、CGC 和KcTr,x。RS 結(jié)果描述了參數(shù)變化對(duì)最終產(chǎn)量值的影響, 在一定程度上忽略了各參數(shù)間的相互影響及不同生長(zhǎng)發(fā)育階段的影響程度變化, 如CGC 的減小會(huì)導(dǎo)致最大冠層覆蓋日的延遲到來(lái); CGC 降低或最大冠層覆蓋度過(guò)大均會(huì)使HI0降低等。盡管該方法得到的敏感性分析結(jié)果存在一定的不確定性, 但已達(dá)到了篩選敏感性參數(shù)的目的, 可服務(wù)于參數(shù)率定。

2.2 典型作物參數(shù)率定

基于敏感性分析結(jié)果, 玉米及大豆的敏感性參數(shù)差異較大, 需對(duì)兩種作物分別進(jìn)行參數(shù)率定。玉米作物主要率定了HI0、Zx、WP*等參數(shù), 率定值與默認(rèn)值最大差距 230%, 率定后保持不變的有Pexp-up 和WP*yf。大豆作物主要率定了CGC、CDC、Zx、WP*yf 等參數(shù), 率定值與默認(rèn)值最大差距166%,率定后保持不變的有KcTr,x和Rexshp。表4 展示了AquaCrop 模型模擬東北黑土區(qū)玉米及大豆生長(zhǎng)的率定后參數(shù)。對(duì)水分脅迫參數(shù)調(diào)整較大, 尤其是Pexp-up、限制冠層伸展的土壤水分消耗下限閾值(soil water depletion threshold for canopy expansion-lower threshold, Pexp-lw)、Pexshp 等限制冠層生長(zhǎng)的水分脅迫參數(shù)在率定后使模擬結(jié)果精確度有了很大提高, 更適用于研究區(qū)土壤及氣候大環(huán)境。率定后的AquaCrop 模型參數(shù)為精確模擬東北黑土區(qū)玉米和大豆產(chǎn)量奠定了基礎(chǔ), 也為后續(xù)該地區(qū)作物研究提供了參考。

表4 AquaCrop 模型模擬東北黑土區(qū)典型作物生長(zhǎng)的參數(shù)率定結(jié)果Table 4 Calibrated parameters for simulating black soil area in Northeast China

2.3 模型在東北黑土區(qū)的適應(yīng)性

利用大田監(jiān)測(cè)點(diǎn)的18 組玉米和16 組大豆的作物實(shí)測(cè)產(chǎn)量, 對(duì)率定參數(shù)前后的實(shí)測(cè)產(chǎn)量與模型預(yù)測(cè)產(chǎn)量進(jìn)行對(duì)比, 結(jié)果見(jiàn)圖1。圖中可見(jiàn): (1)參數(shù)率定前后, 模型預(yù)測(cè)精度提升顯著。玉米的預(yù)測(cè)產(chǎn)量回歸系數(shù)由0.34 提升至0.89, 回歸方程截距由7.50 t hm–2降至1.42 t hm–2,R2由0.4797 提升至0.7745,RMSE 由3.271 t hm–2降至1.076 t hm–2, NRMSE 由0.293 降至0.097, ME 由1.333 降至0.747; 大豆的預(yù)測(cè)產(chǎn)量回歸系數(shù)由0.80 提升至0.88, 回歸方程截距由0.34 t hm–2降至0.29 t hm–2,R2由0.5811 提升至0.7789, RMSE 由0.941 t hm–2降至0.299 t hm–2,NRMSE 由0.561 降至0.178, ME 由1.674 降至0.747。參數(shù)率定修正了預(yù)測(cè)值在趨勢(shì)上的偏離和系統(tǒng)性誤差, 提高了模型的適用能力。(2)玉米和大豆的實(shí)測(cè)產(chǎn)量與率定后預(yù)測(cè)產(chǎn)量的R2分別為0.775 和0.779,RMSE 分別為1.076 t hm–2和0.299 t hm–2, NRMSE分別為0.097 和0.178, ME 分別為0.747 和0.730, 均表明率定后的模型具有較好的預(yù)測(cè)效果。從玉米和大豆的預(yù)測(cè)結(jié)果可見(jiàn), 針對(duì)不同作物的不同產(chǎn)量狀態(tài), AquaCrop 模型的預(yù)測(cè)效果存在差異。例如玉米低產(chǎn)時(shí)的模擬精度較低, 回歸方程截距較大, 但相較于公頃產(chǎn)量而言誤差在可接受范圍內(nèi); 而大豆無(wú)論是高產(chǎn)還是低產(chǎn)時(shí), 其模擬精度基本保持無(wú)偏。

3 討論

相對(duì)敏感性分析結(jié)果發(fā)現(xiàn), 對(duì)玉米和大豆產(chǎn)量比較敏感的參數(shù)可分為3 類: 一是與冠層相關(guān)的參數(shù), 如Pexshp 和CDC; 二是與根系相關(guān)的參數(shù), 如Rexshp 和Zx; 三是與作物遺傳特征相關(guān)的參數(shù), 如KcTr,x、WP*和HI0。葉片既是光合作用的場(chǎng)所又是植物蒸騰的器官, 植物所需水分來(lái)自根系吸水, 作物系數(shù)水分生產(chǎn)力決定了作物需水量和水分生產(chǎn)效率,收獲指數(shù)則決定了生物量中經(jīng)濟(jì)產(chǎn)量的高低, 因此上述參數(shù)的微小變化會(huì)使產(chǎn)量發(fā)生較大變化。然而玉米和大豆卻對(duì)3 類參數(shù)的敏感度有所差異。在與冠層相關(guān)的參數(shù)中, 玉米對(duì)冠層衰減系數(shù)更為敏感,而大豆則對(duì)影響冠層生長(zhǎng)的水分脅迫參數(shù)更為敏感。從根系參數(shù)看, 由于大豆根系一般比玉米根系分布淺, 大豆產(chǎn)量對(duì)根系生長(zhǎng)的形狀參數(shù)更為敏感,而玉米產(chǎn)量則對(duì)最大有效根深更為敏感。由于玉米需水量大且對(duì)水分變化敏感, 其對(duì)影響作物需水量大小的作物系數(shù)和水分生產(chǎn)力更為敏感, 大豆則表現(xiàn)為一般敏感。兩種作物都對(duì)影響經(jīng)濟(jì)產(chǎn)量的收獲指數(shù)非常敏感。綜上, 敏感性分析不僅反映出模型模擬作物產(chǎn)量的關(guān)鍵影響過(guò)程, 還揭示了對(duì)不同作物的差異, 具有很好的適用性。實(shí)際應(yīng)用中, 應(yīng)確保這些敏感參數(shù)輸入的準(zhǔn)確性。

本文采用的OAT 方法屬于局部敏感性分析方法的一種, 計(jì)算效率較高, 分析的是單個(gè)參數(shù)變化對(duì)模型輸出變量的影響, 一定程度上忽略了參數(shù)間的相關(guān)性。與之相比, 以擴(kuò)展傅里葉幅度檢驗(yàn)法(extended fourier amplitude sensitivity test, EFAST)為代表的全局敏感性分析方法[30-31]考慮了參數(shù)間交叉作用和作物在不同生長(zhǎng)階段的參數(shù)變化, 更為全面。已有研究采用EFAST 方法對(duì)AquaCrop 模型進(jìn)行了敏感性分析[32], 結(jié)果表明在北京地區(qū)對(duì)冬小麥產(chǎn)量的強(qiáng)敏感性參數(shù)包括: Pexp-up、HI0、KcTr,x、CDC、WP*。與本文OAT 方法得出的玉米產(chǎn)量敏感參數(shù)具有很高的一致性, 與大豆產(chǎn)量敏感參數(shù)也存在部分重疊, 也說(shuō)明采用 OAT 方法分析的AquaCrop 模型敏感性參數(shù)是合理的。AquaCrop 模型在松嫩平原預(yù)測(cè)春小麥產(chǎn)量時(shí)[33], 模擬產(chǎn)量與實(shí)測(cè)產(chǎn)量的R2在0.887～0.969 之間, RMSE 為0.11 t hm–2, 模擬效果好于本研究, 說(shuō)明該模型在東北黑土區(qū)具有良好的適用性和發(fā)展前景, 應(yīng)選擇更多的地區(qū)驗(yàn)證玉米和大豆的產(chǎn)量模擬結(jié)果。

AquaCrop 模型中許多假設(shè)條件均基于理想情況, 諸如土壤層分布均勻、無(wú)病蟲(chóng)害、無(wú)雜草、無(wú)養(yǎng)分限制等, 未考慮自然或人為影響造成的突發(fā)情況, 導(dǎo)致模擬結(jié)果與實(shí)際產(chǎn)量存在偏差。對(duì)敏感性高的參數(shù)進(jìn)行率定可以提高模型的模擬精度。作物生長(zhǎng)機(jī)理模型在單點(diǎn)上已經(jīng)有了較高的模擬精度, 將其擴(kuò)展到更大的區(qū)域范圍進(jìn)行連續(xù)模擬時(shí), 會(huì)由于單點(diǎn)數(shù)據(jù)的隨機(jī)性、模擬過(guò)程的復(fù)雜性以及變量的自相關(guān)等, 導(dǎo)致模擬結(jié)果具有不確定性。如何融合單點(diǎn)觀測(cè)和多源觀測(cè)如遙感數(shù)據(jù), 是將作物模型擴(kuò)展到更大空間尺度, 并提高模擬精度的發(fā)展方向。

遙感衛(wèi)星數(shù)據(jù)[34]能夠快速獲取大尺度地表信息而在農(nóng)業(yè)監(jiān)測(cè)與作物估產(chǎn)中有巨大應(yīng)用空間。集合卡爾曼濾波(ensemble kalman filtering, EnKF)[35]等數(shù)據(jù)同化算法能以實(shí)時(shí)觀測(cè)數(shù)據(jù)對(duì)先驗(yàn)的模型數(shù)據(jù)結(jié)果進(jìn)行再分析, 改進(jìn)了系統(tǒng)的可靠性, 從而為作物產(chǎn)量預(yù)測(cè)提供更大參考價(jià)值。研究還需結(jié)合多源數(shù)據(jù)和同化方法進(jìn)一步優(yōu)化AquaCrop 模型在東北黑土區(qū)的作物產(chǎn)量模擬結(jié)果, 控制和減小模型模擬不確定性對(duì)產(chǎn)量預(yù)測(cè)造成的誤差。

4 結(jié)論

AquaCrop 模型的參數(shù)敏感性分析結(jié)果表明水分脅迫等環(huán)境變量對(duì)作物產(chǎn)量有較大影響, 主要集中在以下幾個(gè)參數(shù): Pexshp、Rexshp、KcTr,x、WP*、HI0、CDC 和Zx。不同作物的參數(shù)敏感性差異較大。精準(zhǔn)的田間觀測(cè)可實(shí)現(xiàn)模型參數(shù)的率定, 有效彌補(bǔ)了模型在客觀生物理化意義上的缺失, 降低了作物模型中敏感性參數(shù)引起的誤差, 可提高模型的預(yù)測(cè)精度。率定的參數(shù)越豐富, 模型模擬效果越好, 率定參數(shù)的個(gè)數(shù)與模擬精度呈正相關(guān)。玉米和大豆的實(shí)測(cè)產(chǎn)量與參數(shù)率定后AquaCrop 模型模擬的產(chǎn)量之間差異較小, 其R2分別為0.775 和0.779, RMSE 分別為1.076 t hm–2和0.299 t hm–2, NRMSE 分別為0.097 和0.178, ME 分別為0.747 和0.730。表明參數(shù)率定后的AquaCrop 模型可較準(zhǔn)確地估算東北黑土區(qū)玉米和大豆產(chǎn)量, 適用于作物品種和環(huán)境條件相似的作物產(chǎn)量模擬與預(yù)測(cè)。