陳桂良，劉忠妹，許木果，黎小清，丁華平，楊春霞

（云南省熱帶作物科學(xué)研究所，云南景洪 666100）

0 引言

云南是中國種植面積最大、產(chǎn)膠最多、單產(chǎn)最高的優(yōu)質(zhì)天然橡膠生產(chǎn)基地[1]。云南橡膠園為典型的山地膠園，地形地貌復(fù)雜，小氣候復(fù)雜多變，成土母質(zhì)豐富，土壤養(yǎng)分空間差異明顯[2-3]，植膠初期土壤肥力較高，但在長期植膠后，土壤退化嚴(yán)重[4-5]。土壤有機(jī)質(zhì)是反映橡膠園土壤養(yǎng)分狀況的重要指標(biāo)，其含量變化受到廣泛關(guān)注[6-8]。傳統(tǒng)化學(xué)分析方法存在周期長、成本高、污染環(huán)境等不足，而高光譜技術(shù)憑借快速、簡便、無污染等優(yōu)勢(shì)，為土壤有機(jī)質(zhì)含量的快速測(cè)定開辟了新的途徑[9-10]。土壤有機(jī)質(zhì)含量快速、高精度估算模型的建立，可以更好地指導(dǎo)云南山地膠園精細(xì)化生產(chǎn)管理，對(duì)云南天然橡膠產(chǎn)業(yè)的高質(zhì)量發(fā)展意義重大。

國內(nèi)外學(xué)者普遍認(rèn)為土壤光譜信息中存在有機(jī)質(zhì)敏感波段，利用土壤有機(jī)質(zhì)敏感波長建立土壤有機(jī)質(zhì)含量估算模型，可以簡化模型，保持甚至提升模型精度[11-13]。為突出土壤有機(jī)質(zhì)光譜信號(hào)，通常會(huì)對(duì)原始光譜反射率進(jìn)行適當(dāng)變換，然后進(jìn)行特征波長優(yōu)選，再建立相應(yīng)的預(yù)測(cè)模型，以提高模型的穩(wěn)定性和精準(zhǔn)度[14-17]。有關(guān)光譜反演土壤有機(jī)質(zhì)含量的研究中，模型的建立大多數(shù)以經(jīng)驗(yàn)統(tǒng)計(jì)方法為主，由于不同區(qū)域土壤類型多樣且成因復(fù)雜，使得高光譜特征和估算模型的差異很大，難以建立統(tǒng)一的有機(jī)質(zhì)估算模型[18]。為了提高模型的有效性，針對(duì)特定區(qū)域建立區(qū)域化的有機(jī)質(zhì)光譜預(yù)測(cè)模型是常用手段[19]。很多學(xué)者對(duì)國內(nèi)典型區(qū)域進(jìn)行了土壤有機(jī)質(zhì)高光譜估算研究，如東北黑土[20]、干旱區(qū)典型綠洲[21]、紅壤區(qū)[22]、華南地區(qū)[23]、黃河三角洲[24]等。目前還未見利用高光譜技術(shù)估算云南山地膠園土壤有機(jī)質(zhì)的研究報(bào)道。

筆者以云南典型山地膠園東風(fēng)農(nóng)場(chǎng)為研究區(qū)，基于獲取的土壤有機(jī)質(zhì)和光譜反射率數(shù)據(jù)，在優(yōu)選光譜變換模式和特征波長的基礎(chǔ)上，采用多元線性回歸、偏最小二乘回歸和支持向量回歸，構(gòu)建云南山地膠園土壤有機(jī)質(zhì)含量的高光譜最優(yōu)估算模型，以期實(shí)現(xiàn)云南山地膠園土壤有機(jī)質(zhì)含量快速檢測(cè)。

1 材料與方法

1.1 研究區(qū)與土壤樣品采集

云南省景洪市東風(fēng)農(nóng)場(chǎng)位于景洪市西南部，在21°30′00″—21°46′12″N、100°34′48″—100°49′48″E 之間，海拔600～900 m 之間，成土母質(zhì)主要是花崗巖、千枚巖、砂頁巖和老沖積物，土壤主要為磚紅壤和磚紅壤性紅壤，下轄6 個(gè)分場(chǎng)，現(xiàn)有天然橡膠1.02 萬hm2。2019 年9 月，在東風(fēng)農(nóng)場(chǎng)選取了103 個(gè)GPS 定位采樣點(diǎn)，根據(jù)地形的差異，在每個(gè)GPS定位采樣點(diǎn)采集1～3個(gè)橡膠樹保護(hù)帶土壤樣品(0～20 cm)，共采集到土壤樣品225個(gè)。圖1為研究區(qū)位置及采樣點(diǎn)分布圖。

圖1 研究區(qū)位置及采樣點(diǎn)分布圖

1.2 土壤有機(jī)質(zhì)含量與光譜反射率測(cè)定

土壤樣品帶回實(shí)驗(yàn)室后，將每個(gè)土壤樣品混合均勻，剔除樹根、石塊等雜質(zhì)，經(jīng)室內(nèi)自然風(fēng)干后研磨并過0.15 mm 孔徑篩，用于土壤有機(jī)質(zhì)和光譜反射率測(cè)定。采用重鉻酸鉀容量法測(cè)定土壤有機(jī)質(zhì)[25]。采用FieldSpec4光譜儀（美國ASD公司產(chǎn)）進(jìn)行土壤光譜反射率測(cè)定，整個(gè)測(cè)定工作在一個(gè)同等于暗室的實(shí)驗(yàn)室進(jìn)行。用培養(yǎng)皿（半徑10 cm、深1.5 cm）裝滿土壤樣品，并將表面刮平，放在反射率近似為0的黑色橡膠墊上。光源為50 W鹵素?zé)?，光源入射角度?5°，距土壤樣品表面70 cm。傳感器探頭采用25°視場(chǎng)角，位于土壤樣品表面垂直上方15 cm處。光譜測(cè)定之前先進(jìn)行白板校正，測(cè)定期間每隔15 min 進(jìn)行白板校正，每個(gè)土壤樣品采集10條光譜曲線，算術(shù)平均后得到該土壤樣品的光譜反射率數(shù)據(jù)。

1.3 研究方法

1.3.1 光譜預(yù)處理及建模集劃分 剔除噪聲較大的350～399 nm波段，剩余的400～2500 nm光譜反射率數(shù)據(jù)用于本研究的建模研究。將光譜反射率重采樣為5 nm間隔(400、405、410、…、2495、2500 nm)。 采用Kennard-Stone(KS)算法對(duì)225 個(gè)土壤樣本進(jìn)行優(yōu)選，選出150個(gè)作為校正集，剩余75個(gè)作為驗(yàn)證集，校正集和驗(yàn)證集的土壤樣品有機(jī)質(zhì)含量統(tǒng)計(jì)見表1。

表1 用于模型校正和驗(yàn)證的土壤樣品有機(jī)質(zhì)含量統(tǒng)計(jì)

1.3.2 光譜變換為了消除光譜數(shù)據(jù)無關(guān)信息和噪聲，提高建模精度，首先采用3 種方法對(duì)重采樣后的光譜反射率R進(jìn)行光譜變換處理，即倒數(shù)對(duì)數(shù)變換[log(1/R)]、多元散射校正(MSC)、標(biāo)準(zhǔn)正態(tài)變換(SNV)，然后對(duì)光譜反射率R以及3 種變換形式光譜數(shù)據(jù)進(jìn)行Savitzky-Golay(SG)平滑或?qū)?shù)變換，以得到最佳的光譜變換模式。

1.3.3 特征波長篩選方法

（1）競(jìng)爭(zhēng)適應(yīng)重加權(quán)采樣。競(jìng)爭(zhēng)適應(yīng)重加權(quán)采樣(competitive adaptive reweighted sampling，CARS)，通過優(yōu)選偏最小二乘回歸模型回歸系數(shù)絕對(duì)值相對(duì)較大的波長變量，剔除模型回歸系數(shù)絕對(duì)值相對(duì)較小的波長變量，獲得模型交叉驗(yàn)證均方根誤差最小的優(yōu)選變量子集[26-27]。

（2）連續(xù)投影算法。連續(xù)投影算法(successive projections algorithm，SPA)，是在光譜矩陣中應(yīng)用變量投影操作尋找含有冗余信息最低、共線性最少的光譜特征變量組，最大程度地避免光譜信息重疊、簡化模型結(jié)構(gòu)、提高建模的速度和效率[28]。

（3）CARS-SPA。將CARS與SPA 2種方法結(jié)合用于特征波長篩選，首先用CARS 方法獲得優(yōu)選的變量子集，然后再應(yīng)用SPA 方法對(duì)CARS 方法獲得的變量子集進(jìn)一步優(yōu)選。

1.3.4 建模方法

（1）多元線性回歸。多元線性回歸(multiple linear regression，MLR)是由多個(gè)自變量的最優(yōu)組合通過回歸預(yù)測(cè)因變量的一種回歸方法，在已知多組自變量和所對(duì)應(yīng)的單組因變量時(shí)，可以對(duì)它們之間的關(guān)系進(jìn)行很好的擬合和預(yù)測(cè)[29]。

（2）偏最小二乘回歸。偏最小二乘回歸(partial least squares regression，PLSR)廣泛用于高光譜數(shù)據(jù)回歸建模，是建立穩(wěn)健線性光譜定量校正模型的通用方法。PLSR 集主成分分析、線性回歸分析和典型相關(guān)分析的優(yōu)點(diǎn)于一身，適于處理自變量存在多重共線性的回歸問題[30]。

（3）支持向量回歸。支持向量回歸(support vector regression，SVR) 是支持向量機(jī)(support vector machine，SVM)的重要應(yīng)用分支，用于回歸模型的構(gòu)建。SVM 通過核函數(shù)將低維非線性問題轉(zhuǎn)換成高維的線性問題，通過結(jié)構(gòu)風(fēng)險(xiǎn)最小化原則提高泛化能力，在保證最小化樣本的同時(shí)縮小了模型泛化誤差的上界，廣泛應(yīng)用于分類與回歸[31-33]。

1.3.5 模型驗(yàn)證模型驗(yàn)證采用決定系數(shù)(R2)、均方根誤差(RMSE)和相對(duì)分析誤差(RPD)，R2越接近1，RMSE越小，RPD越大，說明模型預(yù)測(cè)效果越好。1.5＜RPD＜2 表明模型只能對(duì)樣品進(jìn)行粗略估算，2.0＜RPD＜2.5表明模型具有較好的估算能力，2.5＜RPD＜3.0 表明模型具有很好的估算能力，RPD＞3.0 表明模型具有極好的估算能力[18]。

2 結(jié)果與分析

2.1 土壤有機(jī)質(zhì)含量與光譜反射率

參照GB/T 29570—2013 附錄C 膠園土壤養(yǎng)分含量正常指標(biāo)[34]，將所有樣本按有機(jī)質(zhì)含量高低劃分為高(＞25 g/kg)、中(20～25 g/kg)、低(＜20 g/kg)3 種有機(jī)質(zhì)含量水平，并計(jì)算每種有機(jī)質(zhì)含量水平土壤樣品的平均反射率光譜，發(fā)現(xiàn)400～2500 nm 波段范圍內(nèi)反射率與有機(jī)質(zhì)含量負(fù)相關(guān)，且不同有機(jī)質(zhì)含量水平土壤的反射率光譜曲線有相似的反射特征，總體呈現(xiàn)先增加后降低的拋物線型（圖2）。在可見光波段，受土壤發(fā)色團(tuán)和有機(jī)質(zhì)本身黑色的影響，存在較寬的吸收波段，反射率總體較低，隨著波長的增加而急劇上升；在900 nm 附近出現(xiàn)典型的Fe3+的吸收谷；在1400、1900、2200 nm 附近有3 個(gè)明顯的吸收峰，深度略有差別，這可能是黏土礦物中含有的水分子和羥基的吸收帶[35]。

圖2 3種有機(jī)質(zhì)含量水平土壤的平均原始反射率光譜

2.2 光譜變換模式篩選

對(duì)重采樣后的光譜反射率R進(jìn)行光譜變換處理后，對(duì)光譜反射率R以及3 種變換形式光譜數(shù)據(jù)進(jìn)行SG 平滑或?qū)?shù)變換模式優(yōu)選。優(yōu)選方法如下：基于SG 平滑或?qū)?shù)變換后的全波段光譜數(shù)據(jù)與有機(jī)質(zhì)含量數(shù)據(jù)，采用PLSR 及留一交叉驗(yàn)證方法建模，其中，SG 平滑或?qū)?shù)變換的濾波窗口在3～101 的奇數(shù)中篩選，多項(xiàng)式次數(shù)1～9 中篩選，導(dǎo)數(shù)階數(shù)在0、1、2 中優(yōu)選，PLSR主成分?jǐn)?shù)在1～20中篩選，按照交叉驗(yàn)證均方根誤差(root mean square error of cross-validation，RMSECV)選擇最佳的光譜變換模式。

log(1/R)、MSC 和SNV 變換在Unscrambler 9.7 中實(shí)現(xiàn)，SG 平滑或?qū)?shù)變換和PLSR 模型的建立通過MATLAB 2015 編程實(shí)現(xiàn)。研究發(fā)現(xiàn)MSC 和SNV 變換對(duì)模型影響不大，log(1/R)變換能夠顯著提升模型性能，log(1/R)變換后再進(jìn)行SG平滑或?qū)?shù)變換，可以進(jìn)一步提升模型預(yù)測(cè)能力，但不同的導(dǎo)數(shù)階數(shù)、濾波窗口大小對(duì)模型的預(yù)測(cè)能力影響較大（圖3）。log(1/R)結(jié)合SG平滑變換是最佳的光譜變換模式，其中，SG平滑變換模式為導(dǎo)數(shù)階數(shù)0，SG 濾波窗口5，多項(xiàng)式次數(shù)2 或3（表2）。

表2 不同光譜變換和不同導(dǎo)數(shù)階數(shù)下最優(yōu)模型預(yù)測(cè)效果

圖3 不同光譜變換下各階導(dǎo)數(shù)不同濾波窗口大小對(duì)應(yīng)的最小交叉驗(yàn)證均方根誤差

采用PLSR 模型和留一交叉驗(yàn)證法，分別基于全波段的光譜反射率R和最佳光譜變換光譜數(shù)據(jù)，建立土壤有機(jī)質(zhì)高光譜估算模型，利用驗(yàn)證集對(duì)模型進(jìn)行驗(yàn)證，相比光譜反射率R，基于最佳光譜變換光譜數(shù)據(jù)建立的模型精度顯著提高，RPD從2.433 增加到2.795（圖4）。表明應(yīng)用本研究的最佳光譜變換模式對(duì)光譜反射率R進(jìn)行光譜變換處理，能夠有效突出土壤有機(jī)質(zhì)光譜信號(hào)，有利于提高模型的穩(wěn)定性和精準(zhǔn)度。

圖4 基于全波段的光譜反射率R(a)和最佳光譜變換光譜數(shù)據(jù)(b)建立的PLSR模型預(yù)測(cè)效果

2.3 特征波長篩選

為進(jìn)一步提高模型的穩(wěn)定性和精準(zhǔn)度，基于篩選得到的最佳光譜變換光譜數(shù)據(jù)與土壤有機(jī)質(zhì)含量數(shù)據(jù)，選擇CARS[36]、SPA[37]、CARS-SPA 3 種特征波長篩選方法，提取得到3個(gè)特征波長組合（表3）。

表3 不同特征波長篩選方法提取的特征波長

2.4 模型構(gòu)建與驗(yàn)證

基于校正集的土壤有機(jī)質(zhì)含量與最佳光譜變換模式下的光譜數(shù)據(jù)，分別選用表3所示的3種特征波長組合對(duì)應(yīng)的光譜特征變量作為自變量，土壤有機(jī)質(zhì)含量作為因變量，并采用MLR、PLSR和SVR 3種方法構(gòu)建土壤有機(jī)質(zhì)高光譜估算模型。其中，MLR 和PLSR 模型在Unscrambler 9.7 中實(shí)現(xiàn)，結(jié)合留一交叉驗(yàn)證方法進(jìn)行建模。SVR 模型利用LIBSVM 軟件[38]實(shí)現(xiàn)，選取SVR 類型為epsilon-SVR，核函數(shù)類型為radial basis function，通過網(wǎng)格搜索和留一交叉驗(yàn)證進(jìn)行參數(shù)優(yōu)選。利用驗(yàn)證集對(duì)所建模型進(jìn)行精度驗(yàn)證，模型的預(yù)測(cè)效果如圖5 所示。從模型預(yù)測(cè)效果來看，CARSSVR最優(yōu)。

圖5 不同模型的預(yù)測(cè)效果

3 結(jié)論

本研究嘗試對(duì)土壤光譜反射率R進(jìn)行多種光譜變換后，再結(jié)合SG平滑或?qū)?shù)變換模式篩選，得到最佳的光譜變換模式為log(1/R)結(jié)合SG 平滑變換，其中，SG 平滑變換模式為導(dǎo)數(shù)階數(shù)0、SG 濾波窗口5、多項(xiàng)式次數(shù)2 或3?；谧罴压庾V變換光譜數(shù)據(jù)與土壤有機(jī)質(zhì)含量數(shù)據(jù)，選擇CARS、SPA、CARS-SPA 3 種方法進(jìn)行特征波長篩選，并采用MLR、PLSR和SVR 3種方法構(gòu)建土壤有機(jī)質(zhì)高光譜估算模型。結(jié)果表明，MLR、PLSR、SVR 3 種模型類型中，最佳模型分別為CARS-MLR、SPA-PLSR、CARS-SVR，RPD分別為2.745、2.617、2.947。CARS-SVR 模型估算精度最高，R2、RMSE、RPD分別為0.897、3.990 g/kg、2.947，該模型RPD位于2.5～3.0之間，具有很好的估算能力。本研究應(yīng)用高光譜技術(shù)，對(duì)云南省景洪市東風(fēng)農(nóng)場(chǎng)典型山地膠園土壤有機(jī)質(zhì)含量進(jìn)行估算，為云南山地膠園土壤有機(jī)質(zhì)含量的快速估算提供了參考。

4 討論

本研究采用3種方法[log(1/R)、MSC、SNV]對(duì)5 nm間隔光譜反射率R(400、405、…、2495、2500 nm)進(jìn)行光譜變換處理，并對(duì)光譜反射率R以及3 種變換形式光譜數(shù)據(jù)進(jìn)行SG 平滑或?qū)?shù)變換模式優(yōu)選，得到最佳的光譜變換模式為log(1/R)結(jié)合SG 平滑變換，其中，SG 平滑變換模式為導(dǎo)數(shù)階數(shù)0、SG 濾波窗口5、多項(xiàng)式次數(shù)2 或3。研究發(fā)現(xiàn)，log(1/R)變換能夠明顯提升土壤有機(jī)質(zhì)估算模型性能，這與方少文等[22]的研究一致。log(1/R)變換后再進(jìn)行SG 平滑或?qū)?shù)變換，可以進(jìn)一步提升模型預(yù)測(cè)能力，但不同的導(dǎo)數(shù)階數(shù)、濾波窗口大小對(duì)模型的預(yù)測(cè)能力影響較大。

基于全波段的最佳光譜變換光譜數(shù)據(jù)建立土壤有機(jī)質(zhì)PLSR模型，并利用驗(yàn)證集對(duì)模型進(jìn)行驗(yàn)證，RPD達(dá)2.795?；诤Y選得到的最佳光譜變換光譜數(shù)據(jù)與土壤有機(jī)質(zhì)含量數(shù)據(jù)，選擇CARS、SPA、CARS-SPA 3種方法進(jìn)行特征波長篩選，并采用MLR、PLSR 和SVR 3 種方法構(gòu)建土壤有機(jī)質(zhì)高光譜估算模型，利用驗(yàn)證集對(duì)各模型進(jìn)行驗(yàn)證，結(jié)果表明，MLR、PLSR、SVR 3 種模型類型中，最佳模型分別為CARS-MLR、SPA-PLSR、CARS-SVR，RPD分別為2.745、2.617、2.947。研究發(fā)現(xiàn)，PLSR 作為一種建立穩(wěn)健線性光譜定量校正模型的通用方法，在全波段模型中具有一定優(yōu)勢(shì)，但在特征波段模型中，不如SVR，甚至不如MLR，與王濤等[39]的研究結(jié)果一致，這可能也與特征波長提取方法有關(guān)。SVR作為一種非線性建模方法，特別是基于特征波長建模，在土壤有機(jī)質(zhì)高光譜估算研究中具有較大潛力。

為了盡可能地提高光譜數(shù)據(jù)與土壤有機(jī)質(zhì)含量之間的相關(guān)性，有效突出土壤有機(jī)質(zhì)光譜信號(hào)，本研究嘗試對(duì)土壤光譜反射率R進(jìn)行多種光譜變換后，再結(jié)合SG平滑或?qū)?shù)變換模式篩選，這為土壤有機(jī)質(zhì)含量高光譜估算模型的建立提供了新思路。本研究采集的土壤樣品包括花崗巖、千枚巖、砂頁巖和老沖積物等成土母質(zhì)，也是云南山地膠園主要母質(zhì)類型[5]。由于建立的有機(jī)質(zhì)估算模型主要針對(duì)云南山地膠園，為了模型的區(qū)域普適性，并未按母質(zhì)分類建模，而是將全部土壤樣品混合建模，但仍然取得了較好的效果。本研究僅用同一農(nóng)場(chǎng)的土壤樣本數(shù)據(jù)對(duì)模型進(jìn)行了驗(yàn)證，在后續(xù)研究中，還需在云南山地膠園開展更為廣泛的模型驗(yàn)證及優(yōu)化。