林陳捷，劉振華，張小媛，胡月明，4，5，劉 洛

(1.華南農(nóng)業(yè)大學(xué)資源環(huán)境學(xué)院，廣東廣州 510642；2.廣東省土地信息工程技術(shù)研究中心，廣東廣州 510642；3.廣東省土地利用與整治重點實驗室，廣東廣州 510642；4.廣州市華南自然資源科學(xué)技術(shù)研究院，廣東廣州 510642；5.海南大學(xué)熱帶作物學(xué)院，海南海口 570228)

耕地是經(jīng)過開墾用以種植農(nóng)作物并經(jīng)常耕耘的土地，是農(nóng)業(yè)生產(chǎn)的重要資源，對保障糧食安全、生態(tài)安全及可持續(xù)發(fā)展具有重要作用[1]。耕地質(zhì)量評價是保障耕地可持續(xù)利用的重要手段[2]，傳統(tǒng)的耕地質(zhì)量制圖主要是通過野外采樣、實驗室化驗，再借助數(shù)理統(tǒng)計方法，以點帶面來評價區(qū)域耕地質(zhì)量狀況[3-4]。該方法依賴大量監(jiān)測樣點，樣點布設(shè)主觀性強(qiáng)、費時費力、效率較低，難以滿足現(xiàn)代農(nóng)業(yè)高效管理的需要。隨著遙感技術(shù)的快速發(fā)展，給精確高效評價耕地質(zhì)量帶來了新思路。

目前，耕地質(zhì)量遙感評價研究主要是基于作物光譜的耕地質(zhì)量制圖，已有學(xué)者嘗試?yán)每梢姽?、近紅外植被指數(shù)進(jìn)行耕地質(zhì)量預(yù)測研究。如Liu等使用增強(qiáng)植被指數(shù)(EVI)和溫度植被干旱指數(shù)(TVDI)，結(jié)合坡度、道路通達(dá)度等構(gòu)建了基于PSR框架的耕地質(zhì)量評價模型，并與傳統(tǒng)指標(biāo)下的壓力—狀態(tài)—響應(yīng)(PSR)框架耕地質(zhì)量評價結(jié)果精度進(jìn)行對比，證明了遙感耕地質(zhì)量制圖的可行性[5]。李因帥等從影像中提取地表含水量指數(shù)(SWCI)和歸一化植被指數(shù)(NDVI)作為特征變量，構(gòu)建了小麥種植區(qū)的耕地質(zhì)量評價模型[6]。官炎俊等從Landsat影像中研究區(qū)的NDVI并結(jié)合研究區(qū)的糧食產(chǎn)量數(shù)據(jù)進(jìn)行修正，構(gòu)建了基于糧食產(chǎn)量修正的NDVI耕地質(zhì)量評價模型[7]。此外，一些學(xué)者使用植被光譜遙感產(chǎn)品進(jìn)行耕地質(zhì)量制圖研究。如馬佳妮等利用時間序列的MODIS數(shù)據(jù)結(jié)合氣象數(shù)據(jù)，通過光能利用率(VPM)模型獲取研究區(qū)作物的植被凈初級生產(chǎn)力(NPP)，進(jìn)而得到多年NPP的均值，并以此進(jìn)行耕地質(zhì)量評價制圖[8]。Zhu等使用經(jīng)驗貝葉斯克里金插值(EBK)算法對水稻不同生育期的MODIS總初級生產(chǎn)力(GPP)產(chǎn)品進(jìn)行降尺度，構(gòu)建了基于水稻生長期GPP的耕地質(zhì)量評價模型[9]。

盡管這些方法提供了寶貴的科學(xué)經(jīng)驗，但仍沒有統(tǒng)一耕地質(zhì)量的光譜響應(yīng)指標(biāo)，需要引入更多的植被光譜，獲取對耕地質(zhì)量更為敏感的光譜指標(biāo)。因此本研究引入植被紅邊指數(shù)，篩選對耕地質(zhì)量敏感的最佳紅邊指數(shù)，構(gòu)建基于紅邊指數(shù)的耕地質(zhì)量關(guān)聯(lián)模型，并以廣東省廣州市增城區(qū)為例進(jìn)行耕地質(zhì)量遙感制圖。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)位于廣東省廣州市增城區(qū)(23°05′～23°37′N，113°32′～114°00′E)，屬于南亞熱帶海洋性季風(fēng)氣候，炎熱多雨，長夏無冬，全年可種植作物。增城區(qū)是全國著名的絲苗米之鄉(xiāng)，2019年增城區(qū)統(tǒng)計年鑒顯示，研究區(qū)耕地的主要作物為水稻，因此，本研究選擇種植水稻的耕地作為研究對象。在研究區(qū)內(nèi)，根據(jù)分層抽樣的方法采集400個耕地樣點，并利用隨機(jī)抽樣按3 ∶1的比例對樣點劃分建模集(圖1-c中的紅點)和驗證集(圖1-c中的黑點)。樣點的耕地質(zhì)量數(shù)據(jù)從2018年增城區(qū)耕地質(zhì)量數(shù)據(jù)庫獲得。樣點布設(shè)過程中，充分考慮了耕地質(zhì)量的分布情況，并通過耕地利用現(xiàn)狀數(shù)據(jù)庫和影像目視解譯確定樣點耕地的利用類型。此外，為了驗證基于紅邊指數(shù)的耕地質(zhì)量制圖精度，另選取了100個樣點作為制圖精度驗證點(圖1-c中的黃點)。

1.2 影像數(shù)據(jù)及預(yù)處理

Sentinel-2是攜帶多光譜成像儀(MSI)的高分辨率多光譜成像衛(wèi)星，MSI可覆蓋13個光譜波段(從可見光和近紅外到短波紅外，具有不同的空間分辨率)，幅寬可達(dá)290 km(表1)。Sentinel-2影像是唯一在紅邊范圍含有3個波段的遙感數(shù)據(jù)，可高效監(jiān)測植被健康信息。根據(jù)前人的研究，水稻抽穗期的影像是評價耕地質(zhì)量的最佳時期[10]，因此，本研究選擇覆蓋研究區(qū)的抽穗期Sentinel-2影像(影像獲取于2018年10月7日)作為遙感數(shù)據(jù)來源。Sentinel-2影像來源于哥白尼開放訪問中心(https：//scihub.copernicus.eu/dhus/#/home)。

為了盡量減小大氣的影響，還需對影像進(jìn)行大氣校正。大氣校正在sen2cor插件上完成，得到消除大氣影響、較真實反映地物的Level-2A數(shù)據(jù)。Sentinel-2的紅邊波段分辨率為20 m，為了提高精度，本研究在snap上對Level-2A數(shù)據(jù)影像進(jìn)行重采樣，得到紅邊波段10 m分辨率的影像。然后將其各波段的像元亮度(DN)值除以10 000，還原為地表反射率。計算公式如下：

ρλ=Qλ/10 000。

(1)

式中：ρλ為λ波段地表反射率；Qλ為λ波段DN值。

1.3 研究方法

1.3.1 紅邊指數(shù) 紅邊是綠色植物在670～760 nm 之間反射率增高最快的點，也是一階導(dǎo)數(shù)光譜在該區(qū)間內(nèi)的拐點[11]。它與植被的各種理化參數(shù)是緊密相關(guān)的，間接反映著耕地狀態(tài)的信息，因此本研究從遙感紅邊波段出發(fā)，獲取耕地質(zhì)量光譜響應(yīng)的遙感指標(biāo)。紅邊指數(shù)的計算公式如表2所示。根據(jù)上述紅邊指數(shù)計算公式，本研究通過ENVI 5.3的band math模塊從Sentinel-2影像中獲取14個紅邊指數(shù)。

表1 Sentinel-2影像多光譜波段主要參數(shù)信息

1.3.2 耕地質(zhì)量最佳紅邊光譜指標(biāo)篩選 本研究利用梯度提升樹(GBDT)算法對14個紅邊指數(shù)進(jìn)行篩選。GBDT是一種基于Booting改進(jìn)的決策樹算法，其進(jìn)行特征選擇的主要原理是通過計算每個紅邊指數(shù)在單棵樹中重要度的平均值，將其作為紅邊指數(shù)的特征重要度(FI)，并根據(jù)重要度進(jìn)行特征選擇[21-22]。特征重要度的公式如下：

(2)

式中：Ij(Tm)是紅邊指數(shù)在單棵樹的重要度；M為樹的數(shù)量。單棵樹的特征重要度通過計算按這個特征分裂之后損失的減少值得出。公式如下：

(3)

式中：L為葉子節(jié)點的數(shù)量；即L-1為非葉子節(jié)點的數(shù)量；vt是和節(jié)點t關(guān)聯(lián)的特征；it是節(jié)點分裂之后損失值的平方。

表2 紅邊指數(shù)及計算公式

此外，本研究利用方差膨脹因子(VIF)進(jìn)一步篩選紅邊指數(shù)，從而消減變量間的共線性。當(dāng)0

1.3.3 耕地質(zhì)量預(yù)測模型構(gòu)建 本研究基于所篩選的最佳紅邊指數(shù)，采用偏最小二乘回歸(PLSR)法、嶺回歸(RR)法和BP神經(jīng)網(wǎng)絡(luò)(BPNN)構(gòu)建紅邊指數(shù)和耕地質(zhì)量的關(guān)聯(lián)模型。偏最小二乘回歸法是由伍德和阿巴諾等在1983年提出的一種新型的多元分析方法，該方法結(jié)合了多元線性回歸分析、典型相關(guān)分析以及主成分分析的特點[22]。偏最小二乘回歸算法通過對系統(tǒng)中的數(shù)據(jù)信息進(jìn)行分解和篩選，從而達(dá)到消除變量多重共線性的效果，因此偏最小二乘回歸可以用于解決很多普通多元回歸無法解決的問題。該算法的最終目標(biāo)是輸出一個線性模型，如公式(4)所示。

Y=Xβ+ε。

(4)

式中：Y為經(jīng)過歸一化的因變量(耕地利用等指數(shù))；X為經(jīng)過歸一化處理的自變量(紅邊指數(shù))；β為系數(shù)矩陣；ε為殘差矩陣。

嶺回歸是一種基于最小二乘法的改良算法，其通過放棄最小二乘的無偏性，以損失部分精度為代價，獲得普適性更高的回歸模型[23-24]，其模型可表示為

Y=Xγ(k)+e。

(5)

式中：Y為經(jīng)過歸一化的因變量(耕地利用等指數(shù))；X為經(jīng)過歸一化處理的自變量(紅邊指數(shù))；e為耕地利用等指數(shù)觀測誤差；γ(k)為嶺回歸系數(shù)；γ(k) 的解法如下：

γ(k)=(XTX+kI)-1XTY。

(6)

式中：Y為經(jīng)過歸一化的因變量(耕地利用等指數(shù))；X為經(jīng)過歸一化處理的自變量(紅邊指數(shù))；I是單位矩陣；k是參數(shù)，通過嶺跡圖確認(rèn)。

BP神經(jīng)網(wǎng)絡(luò)由Rumelhart和McClelland為首的科學(xué)家于1986年提出，該算法利用梯度下降法搜索權(quán)值，并使用反向傳播來傳播偏差反復(fù)調(diào)整網(wǎng)絡(luò)，直到模型的仿真性盡可能高時結(jié)束學(xué)習(xí)[25-26]。網(wǎng)絡(luò)主要包含輸入層、隱含層、輸出層3個部分，隱含層和輸出層之間的傳輸函數(shù)如公式(7)所示。

(7)

式中：Oi為輸入層信息，即紅邊指數(shù)；Oj為隱含層信息；ωjt表示輸入層對隱含層的權(quán)值；fi為輸入層對隱含層的傳遞函數(shù)，在本研究中選用Trainlm 函數(shù)；θj為隱含層閾值；m為輸入層節(jié)點數(shù)量。

(8)

式中：Ok為耕地利用等指數(shù)預(yù)測值；fi為隱含層對輸出層傳遞函數(shù)，在本研究中選用Purelin 函數(shù)；ωkj表示隱含層對輸出層的權(quán)值；θk為輸出層的閾值；w為隱含層節(jié)點數(shù)量。

如果隱含層神經(jīng)元數(shù)量過少，網(wǎng)絡(luò)不能具有必要的學(xué)習(xí)能力和信息處理能力。反之，會大大增加網(wǎng)絡(luò)結(jié)構(gòu)的復(fù)雜性，并且更易陷入局部極小點，使網(wǎng)絡(luò)的學(xué)習(xí)速度變得很慢。本研究所選用的隱含層神經(jīng)元數(shù)量由經(jīng)驗公式確定。

ni≤nh≤2ni+1。

(9)

式中：ni為紅邊指數(shù)數(shù)量；本研究為7；nh為隱含層神經(jīng)元數(shù)量。

若預(yù)測值與實測值相差較大，則進(jìn)行誤差前饋，由輸出層向輸入層逐層修正連接權(quán)值，從而減小均方誤差。如公式(10)所示：

(10)

式中：O為實測耕地利用等指數(shù)；Ok為耕地利用指數(shù)預(yù)測值；N為訓(xùn)練次數(shù)；n為樣點數(shù)量。

1.3.4 精度驗證 本研究利用決定系數(shù)(R2)和歸一化均方根誤差(NRMSE)評估模型預(yù)測效果。其中，R2[公式(11)]用于表明預(yù)測值與實測值之間的線性擬合程度，而NRMSE[公式(12)]用于衡量預(yù)測值與實測值之間的偏差。

(11)

(12)

2 結(jié)果與分析

2.1 耕地質(zhì)量光譜響應(yīng)的最佳紅邊指數(shù)

本研究利用Python 3.7軟件中scikit-learn模塊的網(wǎng)格調(diào)參(Grid Search)確定梯度提升樹算法的最佳參數(shù)后，對14個紅邊指數(shù)進(jìn)行特征選取，根據(jù)多次試驗結(jié)果對比分析，設(shè)定特征重要度的閾值為0.080 。當(dāng)紅邊指數(shù)的重要度<0.080時，該指數(shù)將被去除，然后使用VIF檢驗紅邊指數(shù)的共線性，檢驗標(biāo)準(zhǔn)為去除VIF>10的變量[29]，得到了7個與耕地質(zhì)量有高相關(guān)性的紅邊指數(shù)(最佳紅邊指數(shù))，分別為X1：MCARI2、X2：MTCI、X3：NDVIRE2、X4：NDVIRE3、X5：PSRI、X6：REP、X7：SRRE1，如圖2所示。

2.2 耕地質(zhì)量最佳預(yù)測模型

本研究將最佳紅邊指數(shù)作為自變量，耕地利用等指數(shù)作為因變量，使用PLSR、RR和BPNN構(gòu)建兩者的關(guān)系模型。其中，基于PLSR的耕地質(zhì)量關(guān)系模型為

Y=-0.053 2×X1-0.006 0×X2+0.128 7×X3-0.123 3×X4+0.078 2×X5-0.009 5×X6-0.050 3×X7。

(13)

基于RR的耕地質(zhì)量關(guān)系模型為

Y=-0.319 0×X1-0.103 8×X2-1.121 0×X3-1.012 4×X4+0.773 4×X5-0.089 0×X6-0.208 9×X7。

(14)

本研究中，BPNN包含了3個隱含層和9個神經(jīng)元節(jié)點，迭代次數(shù)為1 000次，模型的建立和運行在MATLAB R2017a上完成。3個模型預(yù)測結(jié)果如圖3所示。其中，基于BPNN模型所建立的預(yù)測模型表現(xiàn)出了最好的預(yù)測能力，其R2、NRMSE分別為0.70，10.00%；PLSR模型和RR模型的擬合精度較為接近，R2均小于0.45。

本研究使用100個樣點作為驗證集，驗證結(jié)果見圖4。BPNN模型(R2=0.63)的驗證效果明顯優(yōu)于其他2種模型，這與訓(xùn)練集的精度結(jié)果相似，進(jìn)一步說明了建模效果的可靠性。

由圖3和圖4可知，BPNN模型散點分布更為接近1 ∶1線，這說明在本研究中，BPNN模型具有較高的預(yù)測能力。PLSR模型和RR模型在耕地利用指數(shù)為1 500～1 700區(qū)間具有較好的擬合能力，而對于兩端的數(shù)據(jù)(利用指數(shù)小于1 500或大于 1 700)的偏差值較大。已有研究表明，BPNN在處理復(fù)雜的非線性建模時有著很好的擬合能力[23，30]，而PLSR和RR通常用來確定自變量和因變量之間的線性關(guān)系。本研究中，非線性模型(BPNN)的驗證精度明顯高于其他2種線性模型(PLSR和 RR)，其原因可能是耕地質(zhì)量和紅邊指數(shù)之間存在明顯的非線性的關(guān)聯(lián)關(guān)系。

2.3 研究區(qū)耕地質(zhì)量遙感空間制圖

本研究利用REI-BPNN模型對影像進(jìn)行耕地質(zhì)量遙感空間制圖，并根據(jù)農(nóng)用地分等定級規(guī)程(GB/T 28407—2012《農(nóng)業(yè)用地質(zhì)量分級規(guī)定》)將結(jié)果劃分為4個等級(圖5-a)。結(jié)果表明，研究區(qū)的西南部耕地等級較高，10等耕地主要集中在西南部，9等耕地分布較為均勻，研究區(qū)北部耕地等級相較于南部偏低，遙感耕地制圖的耕地等級分布趨勢與耕地質(zhì)量數(shù)據(jù)庫(圖5-b)較為一致。為了驗證使用BPNN模型進(jìn)行耕地質(zhì)量制圖的精度，本研究選取了100個樣點(圖1-c中的黃點)計算其R2和NRMSE。結(jié)果顯示，REI-BPNN模型的耕地質(zhì)量制圖精度R2為0.60，NRMSE為14.80%，表明利用REI-BPNN模型進(jìn)行耕地質(zhì)量遙感空間制圖具有可行性(圖6)。

3 討論

耕地質(zhì)量評價是保護(hù)耕地的重要措施，高效、快速監(jiān)測耕地質(zhì)量具有重要意義。針對當(dāng)前耕地質(zhì)量評價費時費力的問題，已有學(xué)者嘗試使用遙感技術(shù)進(jìn)行耕地質(zhì)量評價。當(dāng)前已有研究利用作物光譜進(jìn)行耕地質(zhì)量評價，本研究首次引入了紅邊指數(shù)進(jìn)行耕地質(zhì)量預(yù)測，其構(gòu)建的耕地質(zhì)量光譜響應(yīng)預(yù)測模型，NRMSE為10.00%，與前人的研究[4]相比，精度有所提高。

通過對比3種模型(PLSR、RR、BPNN)，BPNN模型表現(xiàn)出最佳的預(yù)測效果)，而RR和PLSR等2種線性模型，在建模精度上較為接近。這表明耕地質(zhì)量與光譜變量之間可能存在非線性關(guān)系，這與前人的研究[4，9]基本一致。

為了驗證模型在區(qū)域尺度預(yù)測耕地質(zhì)量的精度，將REI-BPNN模型運用與Sentinel-2影像相結(jié)合，對研究區(qū)進(jìn)行基于紅邊指數(shù)的耕地質(zhì)量制圖。對比前人的研究[4]，制圖精度有所提高(R2=0.60，NRMSE=14.80%)，紅邊指數(shù)的加入在一定程度上提高了耕地質(zhì)量預(yù)測的精度，說明紅邊指數(shù)在耕地質(zhì)量預(yù)測上具有一定優(yōu)勢。

由于BPNN模型的初始參數(shù)(如神經(jīng)元節(jié)點數(shù))具有很大不確定性，可能出現(xiàn)過擬合或欠擬合的情況，對模型的精度和穩(wěn)定性產(chǎn)生一定影響[31]。因此，下一步研究中，將考慮對BPNN模型進(jìn)行優(yōu)化改進(jìn)(如粒子群優(yōu)化)，以提高其精度和穩(wěn)定性。

盡管本研究證明了紅邊波段驅(qū)動下遙感耕地質(zhì)量評價具有極大的潛力，但僅以水稻種植區(qū)為研究對象，為了進(jìn)一步驗證紅邊波段和耕地質(zhì)量的關(guān)系，在下一步研究中，我們將把研究范圍擴(kuò)大到其他利用類型的耕地(如花生地、玉米地)中，建立普適性更強(qiáng)的耕地質(zhì)量評價模型。

4 結(jié)論

利用遙感技術(shù)實現(xiàn)快速、高效耕地質(zhì)量評價對于保護(hù)耕地具有重要意義，本研究以廣州市增城區(qū)為研究對象，利用3種算法(PLSR、RR、BPNN)構(gòu)建紅邊指數(shù)與耕地質(zhì)量的關(guān)聯(lián)模型，在此基礎(chǔ)上，利用Sentinel-2影像進(jìn)行耕地質(zhì)量制圖。結(jié)果表明：(1)本研究利用GBDT算法和VIF篩選的7個最佳紅邊指數(shù)是可靠的。(2)基于最佳紅邊指數(shù)，利用PLSR、RR、BPNN等3種算法構(gòu)建的模型中，BP神經(jīng)網(wǎng)絡(luò)構(gòu)建的模型表現(xiàn)出了最佳的精度(R2、NRMSE分別為0.70，10.00%)，反映了耕地質(zhì)量和紅邊指數(shù)之間存在非線性關(guān)聯(lián)關(guān)系。(3)基于耕地質(zhì)量遙感最佳光譜響應(yīng)模型，利用Sentinel-2影像進(jìn)行空間尺度的耕地質(zhì)量制圖，其制圖精度可靠(R2為0.60，NRMSE為14.80%)。本研究為耕地質(zhì)量提供了一種新的評價方式，為耕地質(zhì)量快速精準(zhǔn)評價提供新的思路。