許寧，宋建新，王淼，李佳，梁中欽，顧曉鶴，鄧光遠(yuǎn)，單偉東，范鳳翠，毛婭楠

（1.河北省農(nóng)業(yè)技術(shù)推廣總站，河北 石家莊 050011；2.河北省畜牧站，河北 石家莊 050000；3.平山縣農(nóng)牧局，河北 平山050400；4.國(guó)家農(nóng)業(yè)信息化工程技術(shù)研究中心，北京 100097；5.高邑縣大營(yíng)農(nóng)業(yè)技術(shù)推廣區(qū)域站，河北 高邑 051330；6.河北省農(nóng)林科學(xué)院農(nóng)業(yè)信息與經(jīng)濟(jì)研究所，河北 石家莊 050051；7.河北省農(nóng)作物引育種中心，河北 石家莊 050000）

耕地是重要的農(nóng)業(yè)資源之一，也是糧食生產(chǎn)最基本的資源條件，對(duì)保障國(guó)家糧食安全以及維護(hù)社會(huì)經(jīng)濟(jì)穩(wěn)定與可持續(xù)發(fā)展具有重要意義。2018年“中央1 號(hào)”文件明確指出，要大力實(shí)施耕地質(zhì)量保護(hù)與提升行動(dòng)。河北省是農(nóng)業(yè)大省，耕地面積約633.33 萬hm2，其中糧食作物占用耕地400 萬hm2。推進(jìn)耕地質(zhì)量保護(hù)與提升，是農(nóng)業(yè)管理部門的工作重點(diǎn)之一。開展耕地質(zhì)量遙感評(píng)估[1～4]，掌握中低產(chǎn)田的空間分布現(xiàn)狀，對(duì)于指導(dǎo)河北省耕地質(zhì)量保護(hù)提升、糧食生產(chǎn)可持續(xù)發(fā)展、水肥高效精準(zhǔn)管理具有重要意義[5～7]。

傳統(tǒng)的耕地質(zhì)量評(píng)價(jià)主要是通過人為選取采樣點(diǎn)，然后對(duì)樣點(diǎn)的指標(biāo)數(shù)據(jù)進(jìn)行整理和分析，以點(diǎn)代面得到某個(gè)區(qū)域的耕地質(zhì)量狀況。該方法存在工作周期長(zhǎng)、空間分布廣等特點(diǎn)，導(dǎo)致數(shù)據(jù)較為粗糙，影響了成果在生產(chǎn)實(shí)踐中的應(yīng)用[8～12]。隨著我國(guó)衛(wèi)星觀測(cè)體系的日益完善和遙感技術(shù)的深入應(yīng)用，高空間、高時(shí)間、高光譜的衛(wèi)星遙感數(shù)據(jù)將有助于解決耕地質(zhì)量綜合評(píng)價(jià)的多個(gè)瓶頸問題[13，14]。遙感技術(shù)能夠快速監(jiān)測(cè)耕地的土壤養(yǎng)分，動(dòng)態(tài)監(jiān)測(cè)多年時(shí)間序列的作物收獲產(chǎn)量和全生育期的長(zhǎng)勢(shì)信息[15～17]，可為快速、動(dòng)態(tài)、準(zhǔn)確、高效監(jiān)測(cè)耕地綜合質(zhì)量提供強(qiáng)有力的技術(shù)支持。

1 研究思路

研究區(qū)域?yàn)榭h域尺度的大面積冬小麥-夏玉米輪作耕地。基于耕地質(zhì)量?jī)?nèi)涵剖析，將耕地質(zhì)量分為耕地養(yǎng)分質(zhì)量和耕地生產(chǎn)力質(zhì)量。對(duì)于耕地養(yǎng)分質(zhì)量，采用高光譜影像篩選敏感波段和變換算法，基于多元線性回歸法構(gòu)建耕地有機(jī)質(zhì)、全氮、速效鉀、有效磷含量的遙感監(jiān)測(cè)模型。對(duì)于耕地生產(chǎn)力質(zhì)量，利用歷史多年時(shí)間序列的低分辨率遙感數(shù)據(jù)，開展典型作物種植類型遙感分類；利用當(dāng)年實(shí)測(cè)單產(chǎn)樣本，構(gòu)建典型作物遙感估產(chǎn)模型，推算長(zhǎng)時(shí)間序列的典型作物單產(chǎn)，基于多年平均產(chǎn)量水平評(píng)價(jià)耕地生產(chǎn)力（圖1）。在此基礎(chǔ)上，構(gòu)建耕地質(zhì)量綜合評(píng)價(jià)模型，實(shí)現(xiàn)縣域尺度的耕地質(zhì)量空間制圖。

2 數(shù)據(jù)獲取

2.1 衛(wèi)星遙感數(shù)據(jù)獲取

2014年6月10 日，采用環(huán)境減災(zāi)小衛(wèi)星超光譜影像1 張，評(píng)價(jià)安平縣的耕地養(yǎng)分質(zhì)量。通過http：//edcdaac.usgs.gov/main.asp 獲 取2005 ～2014年 時(shí) 間序列MODIS 影像360 張，評(píng)價(jià)該縣的耕地生產(chǎn)力質(zhì)量。

2.2 試驗(yàn)數(shù)據(jù)獲取

獲取研究區(qū)域內(nèi)2014年的測(cè)土配方樣本100 個(gè)，樣本數(shù)據(jù)指標(biāo)包括有機(jī)質(zhì)含量、全氮含量、速效磷含量和速效鉀含量。獲取研究區(qū)域內(nèi)當(dāng)季冬小麥單產(chǎn)實(shí)際測(cè)量樣本20 個(gè)。涉及的樣本點(diǎn)均附有經(jīng)緯度位置信息。

3 研究方法

3.1 數(shù)據(jù)預(yù)處理

利用ENVI 軟件對(duì)獲取的影像進(jìn)行輻射定標(biāo)、大氣校正和幾何糾正等標(biāo)準(zhǔn)化數(shù)據(jù)處理，獲得統(tǒng)一空間坐標(biāo)系統(tǒng)的遙感數(shù)據(jù)。

3.1.1 輻射定標(biāo) 對(duì)遙感影像數(shù)據(jù)的DN 值進(jìn)行輻射定標(biāo)，即把無量綱的DN 值轉(zhuǎn)換為有量綱的分譜輻射亮度值的過程。根據(jù)公式（1），將圖象的輻射亮度轉(zhuǎn)化為反射率[18]：

式中，γ：反射率（%）；d：太陽(yáng)距地面的距離（km）；Esun：太陽(yáng)平均輻射強(qiáng)度〔J/（cm2·min）〕；θ：太陽(yáng)高度角（°）。

3.1.2 大氣校正 利用ENVI 中的FLAASH 模型對(duì)影像進(jìn)行大氣校正。校正內(nèi)容包括設(shè)置尺度轉(zhuǎn)換因子、圖像中心坐標(biāo)、傳感器高度、圖像分辨率、海拔高度、數(shù)據(jù)獲取日期和衛(wèi)星過境時(shí)間、大氣模型、水氣反演和氣溶膠模型等參數(shù)。

3.1.3 幾何糾正 采用二次多項(xiàng)式中的雙線性內(nèi)插法作為幾何精校正數(shù)學(xué)模型對(duì)遙感影像進(jìn)行重采樣。通過目視選取圖像上易分辨且較精細(xì)的特征點(diǎn)，如道路交叉點(diǎn)、河流岔口、建筑邊界、農(nóng)田界限等作為地面控制點(diǎn)，建立圖像坐標(biāo)和參考坐標(biāo)之間的多項(xiàng)式校正模型。

3.2 耕地質(zhì)量遙感評(píng)價(jià)方法

從耕地質(zhì)量?jī)?nèi)涵看，耕地質(zhì)量包括耕地養(yǎng)分質(zhì)量和耕地生產(chǎn)力質(zhì)量。因此，從耕地養(yǎng)分質(zhì)量和土壤生產(chǎn)力質(zhì)量2 個(gè)方面對(duì)耕地質(zhì)量進(jìn)行綜合評(píng)價(jià)。對(duì)于耕地養(yǎng)分質(zhì)量評(píng)價(jià)，采用高光譜影像篩選敏感波段和換算方法，基于多元線性回歸構(gòu)建耕地4 個(gè)養(yǎng)分指標(biāo)的遙感監(jiān)測(cè)模型；對(duì)于耕地生產(chǎn)力質(zhì)量評(píng)價(jià)，利用歷年時(shí)間序列的低分辨率遙感數(shù)據(jù)對(duì)作物種植類型進(jìn)行分類，利用當(dāng)年實(shí)測(cè)作物單產(chǎn)樣本數(shù)據(jù)構(gòu)建典型作物遙感估產(chǎn)模型，估算時(shí)間序列的作物單產(chǎn)，根據(jù)評(píng)價(jià)區(qū)域耕地的產(chǎn)量水平，評(píng)價(jià)耕地生產(chǎn)力質(zhì)量。

3.2.1 耕地養(yǎng)分含量遙感評(píng)價(jià) 以高光譜影像的光譜反射率及多種數(shù)學(xué)變化形式作為特征向量，與土壤樣點(diǎn)的養(yǎng)分（有機(jī)質(zhì)、全氮、有效磷、速效鉀）含量進(jìn)行相關(guān)性分析，分別篩選4 種養(yǎng)分指標(biāo)的敏感光譜波段和變換形式，進(jìn)而采用多元回歸分析方法構(gòu)建耕地養(yǎng)分遙感監(jiān)測(cè)模型，結(jié)合高光譜影像進(jìn)行區(qū)域耕地的養(yǎng)分（有機(jī)質(zhì)、全氮、有效磷、速效鉀）含量空間制圖。

（1）土壤樣點(diǎn)高光譜影像反射率的提取。根據(jù)土壤樣點(diǎn)的經(jīng)緯度位置，將土壤樣點(diǎn)與高光譜影像空間疊加分析，提取每個(gè)土壤樣點(diǎn)的光譜反射率信息。

（2）耕地養(yǎng)分敏感波段的篩選。不同的數(shù)學(xué)變換有助于篩選耕地各種養(yǎng)分的敏感波段，變換形式分別為反射率、反射率的倒數(shù)、反射率的對(duì)數(shù)、反射率對(duì)數(shù)的倒數(shù)、反射率一階微分、反射率倒數(shù)的一階微分、反射率對(duì)數(shù)的一階微分、反射率對(duì)數(shù)的倒數(shù)的一階微分。將土壤樣本的有機(jī)質(zhì)含量、全氮含量、有效磷含量、速效鉀含量分別與光譜反射率及其變化形式進(jìn)行逐波段相關(guān)性分析，得到每個(gè)波段與4 種養(yǎng)分含量的相關(guān)系數(shù)（R），根據(jù)相關(guān)系數(shù)的高低篩選出有機(jī)質(zhì)含量、全氮含量、有效磷含量、速效鉀含量的敏感波段及其最佳數(shù)學(xué)變換形式。

（3）耕地養(yǎng)分遙感監(jiān)測(cè)模型的構(gòu)建。根據(jù)耕地養(yǎng)分指標(biāo)所對(duì)應(yīng)的敏感波段及數(shù)學(xué)變換形式，通過多元線性回歸方法建立預(yù)測(cè)模型（2）。從所有土壤樣本中隨機(jī)抽取2/3 樣本量建立多元線性回歸模型，其余1/3樣本用于監(jiān)測(cè)模型的精度驗(yàn)證，通過決定系數(shù)和相對(duì)均方根誤差（RMSE）進(jìn)行評(píng)價(jià)。

式中，SNi：耕地4 個(gè)養(yǎng)分指標(biāo)的實(shí)測(cè)值；PSNi：耕地4 個(gè)養(yǎng)分指標(biāo)的預(yù)測(cè)值；n：測(cè)土配方的土壤樣本數(shù)量（個(gè)）。

（4）耕地養(yǎng)分指標(biāo)的空間填圖。將高光譜影像敏感波段的反射率及數(shù)學(xué)變換圖層代入上述監(jiān)測(cè)模型，即可計(jì)算得到評(píng)價(jià)區(qū)域耕地養(yǎng)分指標(biāo)的空間分布圖。

3.2.2 耕地生產(chǎn)力質(zhì)量遙感評(píng)價(jià) 影響耕地生產(chǎn)力水平的因素較多，因此，某一年的作物產(chǎn)量并不能代表評(píng)價(jià)區(qū)域耕地生產(chǎn)力質(zhì)量水平。選擇＞10 a 的時(shí)間序列遙感影像開展耕地作物類型遙感分類，構(gòu)建基于群體長(zhǎng)勢(shì)的作物遙感估產(chǎn)模型，實(shí)現(xiàn)長(zhǎng)時(shí)間序列的耕地生產(chǎn)力評(píng)價(jià)。以研究區(qū)域內(nèi)小麥-玉米輪作區(qū)的小麥單產(chǎn)水平代表耕地生產(chǎn)質(zhì)量，利用冬小麥生長(zhǎng)周期的NDVI 信息與野外實(shí)測(cè)單產(chǎn)樣本，構(gòu)建冬小麥單產(chǎn)遙感評(píng)價(jià)模型，分析長(zhǎng)時(shí)間序列的耕地生產(chǎn)力水平變化情況，實(shí)現(xiàn)耕地生產(chǎn)力質(zhì)量空間制圖。

（1）冬小麥的提取。選取3月25日～4月30日（冬小麥起身期至開花期）研究區(qū)域的MODIS影像，計(jì)算低分辨率時(shí)間序列遙感影像的歸一化植被指數(shù)[19]，形成時(shí)間序列NDVI 影像，根據(jù)公式（3）進(jìn)行S-G 濾波處理。以NDVI ＝0 為臨界閾值，當(dāng)NDVI＞0 時(shí)提取為冬小麥，當(dāng)NDVI≤0 時(shí)識(shí)別為非小麥。

（2）冬小麥遙感估產(chǎn)模型的構(gòu)建。在當(dāng)年野外實(shí)測(cè)單產(chǎn)樣本的支持下，獲取冬小麥整個(gè)生長(zhǎng)季內(nèi)的MODIS 16 d 最大化合成NDVI 數(shù)據(jù)，剔除越冬期的無效數(shù)據(jù)層，采用多元線性回歸法構(gòu)建冬小麥單產(chǎn)遙感評(píng)價(jià)模型，得到長(zhǎng)時(shí)間序列的逐年冬小麥單產(chǎn)空間分布圖，計(jì)算冬小麥多年平均單產(chǎn)。具體流程為：以冬小麥分布圖對(duì)長(zhǎng)時(shí)間序列MODIS NDVI 圖進(jìn)行空間裁剪，提取冬小麥分布范圍內(nèi)的NDVI 數(shù)據(jù)集；將當(dāng)年冬小麥實(shí)測(cè)單產(chǎn)樣本與當(dāng)年NDVI 圖像進(jìn)行空間疊加分析，提取每個(gè)樣本地塊對(duì)應(yīng)的各層NDVI 值；采用多元線性回歸方法對(duì)當(dāng)年樣本地塊的實(shí)測(cè)單產(chǎn)和多時(shí)相NDVI 值構(gòu)建冬小麥單產(chǎn)遙感評(píng)價(jià)模型（4），求取各期NDVI 自變量的回歸系數(shù)，并計(jì)算模型決定系數(shù)和RMSE；將當(dāng)年以及歷史各年的NDVI 序列圖像代入上述估產(chǎn)模型，計(jì)算目標(biāo)區(qū)的逐年冬小麥單產(chǎn)，得到每年冬小麥單產(chǎn)空間分布圖；將時(shí)間序列的冬小麥估產(chǎn)圖累加求平均，計(jì)算每個(gè)冬小麥像元的平均單產(chǎn)，表征該耕地像元的生產(chǎn)力水平。

式中，y：冬小麥單位面積產(chǎn)量（kg/hm2）；xi：冬小麥生長(zhǎng)全程中的MODIS 16 d 最大化合成的NDVI數(shù)據(jù)（越冬期冬小麥生長(zhǎng)緩慢，可認(rèn)為該階段NDVI對(duì)冬小麥最終產(chǎn)量貢獻(xiàn)忽略不計(jì)）；an：估產(chǎn)模型各個(gè)自變量的回歸系數(shù)。

3.3 耕地質(zhì)量綜合評(píng)價(jià)

根據(jù)耕地養(yǎng)分質(zhì)量和生產(chǎn)力質(zhì)量對(duì)耕地綜合質(zhì)量的表征能力，構(gòu)建耕地綜合遙感監(jiān)測(cè)評(píng)價(jià)模型，實(shí)現(xiàn)耕地質(zhì)量的綜合評(píng)價(jià)。

3.3.1 評(píng)價(jià)指標(biāo)權(quán)重確定與歸一化處理 各種評(píng)價(jià)因子對(duì)耕地綜合質(zhì)量的貢獻(xiàn)不同，采用特爾菲法[20]確定各種耕地質(zhì)量評(píng)價(jià)因子的權(quán)重。根據(jù)公式（5）和（6），計(jì)算Ej和δj。通過公式（7），對(duì)數(shù)據(jù)進(jìn)行歸一化處理，使各個(gè)指標(biāo)的變化范圍都為［0，1］。

式中，m：專家人數(shù)（人）；aij：第i 位專家對(duì)因子j 的評(píng)分值（分）；Ej：因子j 的平均值；δj：因子j 的方差。

式中，X：各個(gè)指標(biāo)含量；Xmin：該指標(biāo)的最小值；Xmax：該指標(biāo)的最大值。

3.3.2 耕地質(zhì)量綜合評(píng)價(jià)與分級(jí) 以耕地地塊為基本單元，根據(jù)耕地養(yǎng)分質(zhì)量和生產(chǎn)力質(zhì)量的遙感監(jiān)測(cè)結(jié)果，構(gòu)建基于多源遙感信息的耕地質(zhì)量綜合評(píng)價(jià)模型（8）。統(tǒng)計(jì)養(yǎng)分含量平均值（mean）和標(biāo)準(zhǔn)差（SD），采用正態(tài)（偏正態(tài)）統(tǒng)計(jì)理論的雙閾值劃分策略，以mean-2×SD、mean 和mean+2×SD 為耕地綜合質(zhì)量的4個(gè)等級(jí)的劃分閾值，根據(jù)公式（8），在遙感軟件中建立評(píng)價(jià)模型計(jì)算得到目標(biāo)區(qū)域的耕地質(zhì)量綜合評(píng)價(jià)值（Qi）。參照相關(guān)的耕地分等定級(jí)和土地利用土地信息資料對(duì)閾值做出適當(dāng)調(diào)整，實(shí)現(xiàn)評(píng)價(jià)區(qū)域的耕地質(zhì)量分級(jí)。具體劃分標(biāo)準(zhǔn)為：當(dāng)某塊耕地的Qi＞mean+2×SD 時(shí)，判定為一等耕地；當(dāng)mean＜Qi≤mean+2×SD時(shí)，判定為二等耕地；當(dāng)mean-2×SD＜Qi≤mean 時(shí)，判定為三等耕地；當(dāng)Qi≤mean-2×SD 時(shí)，判定為四等耕地。

式中，Q：耕地綜合質(zhì)量；N：耕地養(yǎng)分質(zhì)量；P：耕地生產(chǎn)力質(zhì)量，即長(zhǎng)時(shí)序的冬小麥單產(chǎn)均值（kg/hm2）；ai：各個(gè)評(píng)價(jià)指標(biāo)的權(quán)重；N1：耕地有機(jī)質(zhì)含量（%）；N2：耕地全氮含量（%）；N3：耕地速效鉀含量（%）；N4：耕地有效磷含量（%）。

4 耕地質(zhì)量遙感綜合評(píng)價(jià)方法的應(yīng)用

安平縣位于河北省中南部，地處太行山前沖積扇前緣，境內(nèi)多為滹沱河沖積平原。地勢(shì)平緩，略顯西高東低。屬于半干旱半濕潤(rùn)大陸性季風(fēng)氣候區(qū)，四季分明，干濕交替。耕地面積3.31 萬hm2，土壤母質(zhì)類型以河流沖積物為主，土層深厚。表層土壤質(zhì)地以砂壤、輕壤質(zhì)為主，土壤顏色發(fā)灰。土壤養(yǎng)分包括有機(jī)質(zhì)、全氮、速效磷、速效鉀和微量元素等。

運(yùn)用上述方法，對(duì)安平縣耕地質(zhì)量進(jìn)行遙感綜合評(píng)價(jià)，結(jié)果（圖2）顯示，安平縣一等和二等耕地?cái)?shù)量較少，且分布零散，主要集中在研究區(qū)域的西南部、東北部和北部地區(qū)；三等和四等耕地?cái)?shù)量較多，分布相對(duì)廣泛。

圖2 河北省衡水市安平縣耕地質(zhì)量綜合評(píng)價(jià)圖Fig.2 Comprehensive evaluation map of cultivated land quality in Anping County，Hengshui City，Hebei Province

進(jìn)一步對(duì)三、四等耕地的數(shù)據(jù)進(jìn)行深入分析，發(fā)現(xiàn)在這些耕地中存在3 個(gè)問題： （1）未開發(fā)利用的土地中存在一定面積的鹽堿地。安平縣鹽堿地屬大陸鹽堿化濃縮型，在鹽漬土中以鹽化為主。經(jīng)過多年的鹽堿地綜合治理和農(nóng)業(yè)利用后，該縣耕地鹽堿化程度已經(jīng)不明顯，一般不會(huì)對(duì)農(nóng)業(yè)生產(chǎn)造成大的影響，但在未利用的土地中仍存在一定數(shù)量的鹽化土壤。（2）中低產(chǎn)田仍然存在。在經(jīng)過一系列的鹽化土改良行動(dòng)后，鹽堿地綜合治理取得了一定成效。但是，該研究區(qū)域的土壤為“漏砂型”，保水保肥能力低于其他土壤類型，導(dǎo)致該區(qū)域仍然存在一定數(shù)量的中低產(chǎn)田。（3）土壤環(huán)境質(zhì)量較差。研究區(qū)域的有關(guān)部門在積極促進(jìn)企業(yè)快速成長(zhǎng)過程中，忽略了環(huán)境保護(hù)問題，致使部分耕地受到污染。雖然導(dǎo)致這些問題出現(xiàn)的原因是多方面的，既有客觀原因，也有主觀原因，但總體來看，不良自然環(huán)境條件是形成中低產(chǎn)田的重要客觀原因，水肥不足是形成中低產(chǎn)田的直接因素，不合理的人類活動(dòng)加劇了土壤質(zhì)量問題。

5 對(duì)策建議

5.1 完善灌溉系統(tǒng)建設(shè)，保證小麥生產(chǎn)順利進(jìn)行

安平縣位于缺水干旱區(qū)，對(duì)于小麥生產(chǎn)而言，干旱是主要的災(zāi)害因素。完備的灌溉設(shè)施可為小麥高產(chǎn)穩(wěn)產(chǎn)提供有力保障。研究區(qū)的有關(guān)部門應(yīng)積極進(jìn)行水利設(shè)施基礎(chǔ)建設(shè)，根據(jù)實(shí)際情況，科學(xué)設(shè)計(jì)井間距和灌溉渠系。

5.2 加強(qiáng)土壤培肥措施的實(shí)施，提高土壤地力

根據(jù)安平縣土壤的肥力特征以及水、肥、鹽的相互關(guān)系，采取針對(duì)性措施培肥土壤，提高地力水平。（1）盡可能多地把土壤生產(chǎn)的有機(jī)質(zhì)還回土壤；（2）合理調(diào)整作物布局和耕作措施，減少土壤有機(jī)質(zhì)的消耗。