曹佳萍,張黎明,邱龍霞,邢世和,馬 丹

(福建農(nóng)林大學(xué)資源與環(huán)境學(xué)院/土壤生態(tài)系統(tǒng)健康與調(diào)控福建省高校重點(diǎn)實(shí)驗(yàn)室 福州 350002)

土壤有效磷是表征土壤磷素營(yíng)養(yǎng)豐缺和環(huán)境質(zhì)量?jī)?yōu)劣的重要因子,其空間分布受到外界環(huán)境(如氣象、地形和植被等)和土壤內(nèi)部理化性質(zhì)(如土壤pH 等)的影響而存在變異性,快速而準(zhǔn)確地預(yù)測(cè)區(qū)域表層土壤中有效磷含量和制作空間分布圖對(duì)于肥力監(jiān)測(cè)、農(nóng)業(yè)可持續(xù)發(fā)展和生態(tài)環(huán)境保護(hù)具有重要意義。

傳統(tǒng)的土壤有效磷含量分析與制圖需要采集大量的土壤樣點(diǎn),往往費(fèi)時(shí)、費(fèi)力、成本高,還會(huì)增加生態(tài)環(huán)境風(fēng)險(xiǎn),而且面向較大區(qū)域的制圖精度較低。衛(wèi)星遙感圖像可反映土壤表層信息而被用于數(shù)字土壤制圖研究,很多衛(wèi)星遙感圖像可免費(fèi)獲取,且遙感數(shù)據(jù)空間覆蓋具有連續(xù)性,作為二次數(shù)據(jù)源可提高稀疏樣點(diǎn)土壤屬性的制圖效率,利用衛(wèi)星遙感圖像聯(lián)合其他環(huán)境變量成為稀疏土壤樣點(diǎn)的較大區(qū)域數(shù)字土壤制圖的重要手段。Cambule 等僅通過(guò)60組實(shí)測(cè)樣本和Landsat TM、氣象、地形、母質(zhì)等數(shù)據(jù)信息,實(shí)現(xiàn)Limpopo 國(guó)家公園(10 410 km)土壤有機(jī)碳儲(chǔ)量估算; 龐龍輝等利用MODIS 等數(shù)據(jù)和205 個(gè)稀疏樣點(diǎn)實(shí)現(xiàn)青海省(72.23×10km)全氮、有機(jī)碳等土壤屬性預(yù)測(cè); Henderson 等借助Landsat TM 影像、地形、氣象數(shù)據(jù),僅通過(guò)2124 個(gè)樣點(diǎn)實(shí)現(xiàn)澳大利亞全境(769.2×10km)表層土壤有效磷含量預(yù)測(cè)。與常用于土壤屬性預(yù)測(cè)制圖的Landsat 數(shù)據(jù)相比,Sentinel-2 數(shù)據(jù)具有更高的空間分辨率和3個(gè)紅邊波段。隨著遙感數(shù)據(jù)空間分辨率的提高,土壤屬性制圖精度越高,而且紅邊波段在土壤屬性如有機(jī)質(zhì)和全氮預(yù)測(cè)方面具有明顯的優(yōu)勢(shì)。因此,利用Sentinel-2 數(shù)據(jù)可能會(huì)提高稀疏樣點(diǎn)的土壤屬性制圖精度,利用衛(wèi)星遙感數(shù)據(jù)輔助有限土壤樣點(diǎn)的數(shù)字土壤制圖可能是未來(lái)土壤屬性制圖的重要研究方向之一。

國(guó)內(nèi)外利用遙感數(shù)據(jù)間接獲取與土壤屬性密切的地表光譜信息,結(jié)合地形、氣象等成土環(huán)境因子輔助預(yù)測(cè)表層土壤屬性的相關(guān)研究,主要集中在土壤有機(jī)碳或全氮的預(yù)測(cè),常用的預(yù)測(cè)方法包括線性回歸模型、地統(tǒng)計(jì)模型和機(jī)器學(xué)習(xí)模型。目前土壤磷元素制圖多采用地統(tǒng)計(jì)插值方法,且研究區(qū)多為基于充足采樣點(diǎn)的平坦地區(qū),如Mohamed等通過(guò)普通克里格插值實(shí)現(xiàn)埃及西北海岸土壤有效磷含量預(yù)測(cè)(采樣密度1.54 個(gè)點(diǎn)?km),Shen 等通過(guò)聯(lián)合植被和地形變量的協(xié)同克里格插值實(shí)現(xiàn)黑龍江4 個(gè)黑土區(qū)(采樣密度為0.99～64.52 個(gè)點(diǎn)?km)土壤全磷預(yù)測(cè)。地統(tǒng)計(jì)模型應(yīng)用于土壤屬性制圖較廣,但其預(yù)測(cè)結(jié)果依賴于樣本數(shù)量和空間自相關(guān)性,當(dāng)實(shí)際樣本不滿足二階平穩(wěn)假設(shè)或內(nèi)蘊(yùn)假設(shè)時(shí),會(huì)造成預(yù)測(cè)結(jié)果置信度降低,限制了其在樣本稀疏的復(fù)雜地貌區(qū)土壤屬性制圖的應(yīng)用。而隨機(jī)森林(Random Forest,RF)等機(jī)器學(xué)習(xí)模型通過(guò)充分挖掘土壤屬性與環(huán)境因子之間的關(guān)系,在采樣點(diǎn)稀疏、空間自相關(guān)弱的情況下也能做出穩(wěn)定預(yù)測(cè),且RF相對(duì)于大多數(shù)統(tǒng)計(jì)模型具有處理高維度數(shù)據(jù)、共線性不敏感和防止過(guò)度擬合等優(yōu)勢(shì)。王茵茵等研究證明了RF 在復(fù)雜地貌區(qū)土壤有機(jī)質(zhì)預(yù)測(cè)的有效性; 姜賽平等通過(guò)聯(lián)合地形因子、歸一化植被指數(shù)和土壤類型等輔助因子的RF 模型,實(shí)現(xiàn)基于160個(gè)采樣點(diǎn)的熱帶地形復(fù)雜區(qū)海南島(32 900 km)土壤有機(jī)質(zhì)含量空間預(yù)測(cè)。然而,目前運(yùn)用協(xié)同多源環(huán)境變量的RF 方法預(yù)測(cè)土壤磷素研究較少,尤其是在稀疏采樣的中國(guó)南方丘陵地區(qū)土壤有效磷數(shù)字制圖效果仍不清楚。

中國(guó)南方丘陵地區(qū)地形起伏較大,耕地田塊面積小,分布比較零散,較小區(qū)域內(nèi)或鄰近田塊的土壤有效磷含量的差異較大,相較于平坦地區(qū)采樣點(diǎn)布設(shè)難度大,獲取數(shù)量易受到限制,導(dǎo)致該區(qū)域耕地土壤有效磷制圖精度較低,繪制難度大。在樣點(diǎn)有限的情況下,如何快速準(zhǔn)確地實(shí)現(xiàn)南方丘陵地區(qū)耕地土壤有效磷制圖成為一個(gè)重要的研究課題。因此,本文以典型南方丘陵區(qū)福建省縣域耕地面積最大的建甌市為研究區(qū),利用96 個(gè)稀疏土壤樣點(diǎn)屬性、Sentinel-2 光譜數(shù)據(jù)、氣象和地形數(shù)據(jù)作為輸入變量,通過(guò)皮爾遜相關(guān)分析及后向剔除法篩選環(huán)境因子,構(gòu)建并對(duì)比分析不同環(huán)境變量組合和土壤有效磷之間的RF 模型,篩選最優(yōu)組合以預(yù)測(cè)建甌市耕地表層土壤有效磷含量與空間分布圖,以期為南方丘陵地區(qū)耕地土壤磷元素的精準(zhǔn)管理與環(huán)境評(píng)估提供科學(xué)依據(jù)。

1 材料與方法

1.1 研究區(qū)域概況

建甌市位于福建省北部(圖1),總面積為4233 km,2018年初耕地面積為40 470.7 hm。屬中亞熱帶海洋性季風(fēng)氣候,無(wú)霜期平均為276.6 d,年均溫為18.7 ℃,年均降水量為1662.9 mm,年均日照為1812.7 h。主要地形為東南沿海典型低山丘陵區(qū),地勢(shì)東南高、西南低,四周環(huán)山,中部為河谷平原、丘陵和山間盆谷,海拔為76～1821 m。土壤類型有水稻土、紅壤、黃壤、潮土和紫色土,其中全市耕地總面積的97.65%為水稻土。主要糧食作物為水稻(),種植制度以一年兩熟為主。

圖1 研究區(qū)范圍及土壤采樣點(diǎn)分布圖Fig.1 Map of the study area and soil sampling points distribution in the study area

1.2 數(shù)據(jù)來(lái)源與預(yù)處理

研究區(qū)97 個(gè)耕地表層(0～20 cm)土壤樣點(diǎn)的理化屬性數(shù)據(jù)來(lái)源于國(guó)家農(nóng)業(yè)農(nóng)村部2017年末耕地質(zhì)量監(jiān)測(cè)與評(píng)價(jià)樣點(diǎn)數(shù)據(jù),樣點(diǎn)采集按照密度控制、均勻性和代表性原則,在選定樣點(diǎn)的耕地地塊均勻隨機(jī)地采集10～15 個(gè)耕層土樣混合后用四分法留取1.5 kg 樣品,采樣時(shí)記錄樣點(diǎn)的地理坐標(biāo)、高程等信息。土壤有效磷采用碳酸氫鈉浸提-鉬銻抗比色法測(cè)定,土壤pH 采用酸度計(jì)法測(cè)定。剔除遙感影像上微量云層影響的1 個(gè)采樣點(diǎn),共96 個(gè)土壤樣點(diǎn)參與建模(圖1)。其他空間數(shù)據(jù)包括:1)研究區(qū)土壤pH 柵格數(shù)據(jù),從2017年福建省3634 個(gè)土壤pH 實(shí)測(cè)數(shù)據(jù)插值獲取的10 m×10 m 柵格數(shù)據(jù)中掩膜提取,數(shù)據(jù)來(lái)源及測(cè)定方法同上; 2)研究區(qū)耕地矢量圖斑(圖2),從農(nóng)業(yè)部耕地質(zhì)量調(diào)查與評(píng)價(jià)項(xiàng)目獲取的建甌市1∶50 000 土地利用現(xiàn)狀數(shù)據(jù)庫(kù)中提取。

圖2 建甌市行政區(qū)劃及耕地管理單元Fig.2 Administrative division and cultivated land management unit of Jian’ou City,Fujian Province,China

本研究選取的預(yù)測(cè)變量包括土壤pH、氣象變量、地形變量和遙感變量。氣象變量包括研究區(qū)年均降水量(MAP)和年均溫(MAT),來(lái)源于世界氣象數(shù)據(jù)庫(kù)1 km×1 km 的1970?2000年月平均降水量(1?12月)和月平均氣溫(1?12月)柵格數(shù)據(jù)。地形變量包括高程(DEM)、地形濕度指數(shù)(TWI)等,從ALOS衛(wèi)星PALSAR 空間分辨率為12.5 m 的數(shù)字高程模型中提取。遙感變量包括紅邊波段(B6)、增強(qiáng)植被指數(shù)(EVI)和第一主成分(PCA1)等,從2017年12月建甌市4 景10 m 分辨率的Sentinel-2 數(shù)據(jù)提取。Sentinel-2 數(shù)據(jù)預(yù)處理主要包括輻射定標(biāo)、大氣校正、圖像鑲嵌與裁剪、幾何精校正和主成分分析,其中輻射定標(biāo)和大氣校正在SNAP 軟件處理,將遙感圖像的DN 值轉(zhuǎn)化為地表反射率; 幾何精校正、主成分分析和相關(guān)植被指數(shù)獲取在ENVI 5.3 中完成,Sentinel-2 幾何精校正以精校正的SPOT 影像為參考影像,誤差控制在1 個(gè)像素以內(nèi)。參考數(shù)字土壤制圖相關(guān)文獻(xiàn),將獲取的地形及氣象數(shù)據(jù)通過(guò)最鄰近元法重采樣至10 m 空間分辨率,該方法相對(duì)于其他尺度轉(zhuǎn)換方法減少了數(shù)據(jù)因重采樣帶來(lái)信息損失,地形數(shù)據(jù)重采樣結(jié)果通過(guò)采樣時(shí)獲取的高程數(shù)據(jù)修正誤差。所有輔助建模數(shù)據(jù)統(tǒng)一到西安80 坐標(biāo)系,檢查與土壤樣點(diǎn)屬性空間匹配的一致性,無(wú)誤后進(jìn)行空間屬性連接。

為進(jìn)一步簡(jiǎn)化模型輸入,建模前需對(duì)輔助預(yù)測(cè)因子進(jìn)行優(yōu)選。先通過(guò)皮爾遜相關(guān)分析保留與土壤有效磷相關(guān)性顯著的因子,進(jìn)一步通過(guò)后向剔除法,依據(jù)建模時(shí)將每個(gè)因子依次排除模型后均方根誤差(RMSE)的增減對(duì)變量進(jìn)行篩選,RMSE 增加則保留因子,反之剔除,最終保留8 個(gè)參與建模的輔助因子(表1)。

表1 土壤有效磷建模輔助變量Table 1 Auxiliary variables of soil available phosphorus modeling

1.3 隨機(jī)森林模型

RF 是將多棵相互沒(méi)有關(guān)聯(lián)的決策樹(shù)組成集成決策樹(shù)的算法。每棵決策樹(shù)在構(gòu)建過(guò)程中通過(guò)bootstrap 方法隨機(jī)抽取約2/3 樣本子集,并使用這些子集組建回歸樹(shù),所有決策樹(shù)預(yù)測(cè)結(jié)果的均值為最終的預(yù)測(cè)結(jié)果。未被抽取的樣本則構(gòu)成袋外樣本,進(jìn)一步計(jì)算袋外誤差,其計(jì)算方法與交叉驗(yàn)證類似。因此,RF 模型不需要單獨(dú)進(jìn)行交叉驗(yàn)證。此外,預(yù)測(cè)變量的相對(duì)重要性通過(guò)置換變量時(shí)的平均預(yù)測(cè)精度估算。RF 建模在Python scikit-learn 庫(kù)中Random Forest Regressor 包實(shí)現(xiàn),決策樹(shù)的數(shù)量與決策樹(shù)的最大深度是兩個(gè)關(guān)鍵參數(shù),決策樹(shù)的數(shù)量設(shè)定范圍為100～1000,根據(jù)模型精度選擇最優(yōu)值,決策樹(shù)的最大深度一般設(shè)置為所選環(huán)境因子數(shù)量的1/3。

1.4 模型精度驗(yàn)證

在總樣本集中隨機(jī)選取訓(xùn)練集和驗(yàn)證集,訓(xùn)練集用于土壤有效磷的建模,驗(yàn)證集用于評(píng)價(jià)模型精度,訓(xùn)練集和驗(yàn)證集的樣本數(shù)量比例為9∶1。選取決定系數(shù)()、平均絕對(duì)誤差(MAE)和RMSE 對(duì)模型精度進(jìn)行驗(yàn)證評(píng)價(jià),越接近1,MAE 和RMSE 越小,表明模型精度越高。

2 結(jié)果與分析

2.1 土壤有效磷含量統(tǒng)計(jì)分析

表2 為建甌市表層土壤(0～20 cm)有效磷含量統(tǒng)計(jì)分析結(jié)果。其中土壤樣本總計(jì)96 個(gè),土壤有效磷含量范圍為4.42～211.94 mg?kg,平均值為54.04 mg?kg,

表2 建甌市耕地表層土壤(0～20 cm)有效磷含量統(tǒng)計(jì)分析Table 2 Statistical analysis of available phosphorus content in surface soil (0?20 cm) of cultivated land in Jian’ou City,Fujian Province,China

標(biāo)準(zhǔn)差為44.53 mg?kg,變異系數(shù)為82.39%。訓(xùn)練集與驗(yàn)證集土壤有效磷含量平均值與總樣本相差較小,分別為53.93 mg?kg和55.06 mg?kg,同時(shí)3 個(gè)樣本組的變異系數(shù)范圍為75.92%～83.57%,均為中等變異。

表3 為土壤有效磷含量與環(huán)境變量的相關(guān)性分析結(jié)果。土壤有效磷含量與年均溫(MAT,=0.357,<0.01)、紅邊波段(B6,=0.292,<0.01)、第一主成分(PCA1,=0.278,<0.01)、增強(qiáng)植被指數(shù)(EVI,=0.252,<0.05)顯著正相關(guān),與年均降水量(MAP,=?0.366,<0.01)、高程(DEM,=?0.336,<0.01)、pH(=?0.268,<0.01)和地形濕度指數(shù)(TWI,=?0.216,<0.05)顯著負(fù)相關(guān)。在篩選的所有變量中,氣象變量(MAP、MAT)與土壤有效磷的相關(guān)性最高,高程(DEM)次之,紅邊波段(B6)、第一主成分(PCA1)、pH 和增強(qiáng)型植被指數(shù)(EVI)緊隨其后,地形濕度指數(shù)(TWI)與土壤有效磷的相關(guān)性最低。這表明在南方丘陵地區(qū)縣市域尺度土壤有效磷含量與氣象和地形變量關(guān)系甚為密切,其次地表信息也是影響土壤有效磷含量的重要因素。

表3 建甌市耕地表層土壤(0～20 cm)有效磷含量與環(huán)境因子的相關(guān)性分析Table 3 Correlation analysis of available phosphorus content and environmental factors in surface soil (0?20 cm) of cultivated land in Jian’ou City,Fujian Province,China

2.2 環(huán)境變量的相對(duì)重要性

構(gòu)建基于全部環(huán)境變量的RF 模型時(shí),通過(guò)多次迭代RF 模型取平均值得到預(yù)測(cè)土壤有效磷過(guò)程中變量的相對(duì)重要性分值(圖3),氣象變量、地形變量、遙感變量和土壤pH 的相對(duì)重要性分值依次為30.64%、30.38%、22.87%和16.11%,結(jié)果表明對(duì)于南方丘陵地區(qū)耕地土壤有效磷含量的空間分布主要受氣象變量、地形變量和土壤pH 的綜合影響,而且氣象和地形變量的影響顯著高于遙感反映的地表信息。從單個(gè)因子看,年均溫(MAT)對(duì)南方丘陵地區(qū)耕地土壤有效磷預(yù)測(cè)結(jié)果的貢獻(xiàn)最大,可解釋22.17%的土壤有效磷含量的空間分布; 土壤pH、地形濕度指數(shù)(TWI)和高程(DEM)對(duì)南方丘陵地區(qū)耕地土壤有效磷預(yù)測(cè)結(jié)果的貢獻(xiàn)較大,分別可解釋16.11%、15.51%和14.87%的土壤有效磷含量的空間分布; 增強(qiáng)植被指數(shù)(EVI)、第一主成分(PCA1)和紅邊波段(B6)的貢獻(xiàn)相對(duì)較低,分別可解釋8.26%、7.73%和6.88%的土壤有效磷含量的空間分布。

圖3 全部環(huán)境變量下環(huán)境因子對(duì)預(yù)測(cè)土壤有效磷空間分布的相對(duì)重要性Fig.3 Relative importance of environmental factors for estimation of spatial distribution of soil available phosphorus content under total environmental variables

2.3 土壤有效磷含量預(yù)測(cè)模型及精度評(píng)價(jià)

為了研究不同變量對(duì)土壤有效磷的預(yù)測(cè)效果,分別對(duì)氣象變量、地形變量、遙感變量和土壤pH進(jìn)行不同的組合建立RF 模型以預(yù)測(cè)土壤有效磷含量的空間分布。模型A 為所有變量的組合,包括地形變量、氣象變量、遙感變量和土壤pH,模型B、D 和E 為均包含遙感變量的組合模型,模型C 為不包含遙感變量的土壤pH、地形變量和氣象變量組合(表4)。研究結(jié)果表明基于全變量組合(模型A)構(gòu)建的RF 模型預(yù)測(cè)精度最高(驗(yàn)證集為0.59,MAE 為19.04 mg?kg,RMSE 為25.26 mg?kg),模型可以解釋研究區(qū)59%的土壤有效磷含量空間分布差異。此外,訓(xùn)練集與驗(yàn)證集的相差很小,說(shuō)明模型沒(méi)有過(guò)擬合或者欠擬合,穩(wěn)定性較好。其次是包括氣象變量、地形變量和遙感變量的模型E,再次為模型B 和D,不包含遙感變量組合(模型C)預(yù)測(cè)精度最 低(驗(yàn) 證 集=0.36,MAE=22.66 mg?kg,RMSE=31.71 mg?kg),加入遙感變量的組合(模型A),明顯提高,MAE 和RMSE 誤差明顯下降,說(shuō)明引入Sentinel-2 遙感變量可以提高有效磷含量的預(yù)測(cè)精度,而且結(jié)果表明引入包括地形變量、氣象變量、遙感變量和土壤pH 的全部變量的RF 模型可基于稀疏樣點(diǎn)有效地預(yù)測(cè)丘陵區(qū)土壤有效磷含量。

2.4 土壤有效磷空間分布預(yù)測(cè)

根據(jù)擬合的5 種RF 預(yù)測(cè)模型得到典型南方丘陵地區(qū)福建省建甌市耕地表層土壤(0～20 cm)有效磷含量的空間分布圖(圖4)。在5 種RF 預(yù)測(cè)模型下土壤有效磷空間分布規(guī)律具有一定相似性,呈中西部含量高、東南部含量低的分布趨勢(shì),各模型在東南部及北部地區(qū)預(yù)測(cè)值較為一致,主要集中在小于30 mg?kg的區(qū)間,中西部區(qū)域預(yù)測(cè)差別較為明顯,其中模型C (圖4c)和模型D (圖4d)制圖顯示中西部區(qū)域預(yù)測(cè)值分別集中分布在60～70 mg?kg和50～60 mg?kg區(qū)間。與此相比,模型A (圖4a)、模型B (圖4b)和模型E (圖4e)預(yù)測(cè)分布圖在中西部區(qū)域分級(jí)變化更明顯,空間表達(dá)差異更顯著,但總體來(lái)看,模型A預(yù)測(cè)結(jié)果對(duì)局部刻畫(huà)最為細(xì)致,且模型精度驗(yàn)證效果最好(表4)。同時(shí)模型E 預(yù)測(cè)土壤有效磷空間分布與模型A 差異極小(圖4f),大多數(shù)區(qū)域土壤有效磷的差異值不超過(guò)5 mg?kg左右,考慮到土壤pH 等土壤屬性數(shù)據(jù)獲取難度較大,在缺乏該變量時(shí)完全可以考慮采用模型E 預(yù)測(cè)土壤有效磷含量的空間分布。

表4 5 種環(huán)境變量組合下表層土壤(0～20 cm)有效磷的隨機(jī)森林預(yù)測(cè)結(jié)果Table 4 Performance of random forest under five combinations of environmental variables for estimation of topsoil (0?20 cm) available phosphorus content

圖4 建甌市耕地土壤有效磷含量空間分布預(yù)測(cè)圖Fig.4 Prediction resultes of spatial distribution of soil available phosphorus (SAP) content of surface soil (0?20 cm) of cultivated land in Jian’ou City,Fujian Province,China

進(jìn)一步對(duì)最優(yōu)模型(模型A)預(yù)測(cè)的土壤有效磷空間分布格局進(jìn)行分析,研究區(qū)土壤有效磷含量范圍為20.10～124.14 mg?kg,平均值為44.22 mg?kg,變異系數(shù)為40.29%。從空間分布統(tǒng)計(jì)來(lái)看(表5),福建省建甌市耕地表層土壤有效磷含量高(≥70 mg?kg)的區(qū)域面積為3695.0 hm,占耕地總面積的9.13%,主要分布于建甌市中部的東游鎮(zhèn)、東峰鎮(zhèn)和西北部的徐墩鎮(zhèn)等鄉(xiāng)鎮(zhèn)。耕地表層土壤有效磷含量較低(<40 mg?kg)的區(qū)域面積為16 722.3 hm,占耕地總面積的41.32%,主要分布于建甌市南部的迪口鎮(zhèn)、玉山鎮(zhèn)和西部的房道鎮(zhèn)等鄉(xiāng)鎮(zhèn)。耕地表層土壤有效磷含量中高值(40～70 mg?kg)區(qū)域面積為20 053.4 hm,占耕地總面積的49.55%,主要分布于建甌市中部的東峰鎮(zhèn)、西部的南雅鎮(zhèn)和西北部的徐墩鎮(zhèn)等地區(qū)。相關(guān)研究表明土壤中過(guò)量的有效磷含量易增加生態(tài)環(huán)境風(fēng)險(xiǎn),研究區(qū)土壤有效磷含量高于40 mg?kg的區(qū)域面積為23 748.4 hm,占耕地總面積的58.68%,主要分布于中西部東峰鎮(zhèn)、徐墩鎮(zhèn)、南雅鎮(zhèn)和東游鎮(zhèn)等鄉(xiāng)鎮(zhèn)?？傮w來(lái)看,福建省建甌市耕地表層土壤有效磷含量較為豐富,中西部各鄉(xiāng)鎮(zhèn)含量相對(duì)較高,東南部區(qū)域含量相對(duì)較低?？臻g分布規(guī)律與高程、年均降水量變化規(guī)律相反,與年均溫變化一致(表6)。地勢(shì)相對(duì)較低、年均降水量較低、年均溫較高的中西部河谷平原和盆地區(qū)土壤有效磷含量高; 地勢(shì)相對(duì)較高、年均降水量較高、年均溫較低的東南部山地區(qū)土壤有效磷含量低。

表5 建甌市耕地土壤有效磷含量空間分布面積及其比例統(tǒng)計(jì)Table 5 Statistics of spatial distribution area and proportion of soil available phosphorus content of cultivated land in Jian’ou City,Fujian Province,China

表6 按照建甌市土壤有效磷含量的氣象和地形變量統(tǒng)計(jì)Table 6 Descriptive characteristics of topography and climate conditions described by soil available phosphorus content in Jian’ou City,Fujian Province,China

3 討論

3.1 基于稀疏樣點(diǎn)的土壤有效磷含量預(yù)測(cè)方法

近年來(lái),很多研究利用遙感變量、氣象變量和地形變量預(yù)測(cè)土壤屬性,然而很少利用優(yōu)于30 m的高空間分辨率的遙感數(shù)據(jù),利用紅邊波段作為預(yù)測(cè)因子的研究則更少,主要是因?yàn)檩^難獲取高空間分辨率的遙感數(shù)據(jù),而且具有紅邊波段的遙感數(shù)據(jù)較少。Henderson 等通過(guò)Landsat TM (30 m)、氣象、地形等預(yù)測(cè)土壤有效磷為0.35。本研究引入空間分辨率為10 m 的Sentinel-2 多光譜變量,結(jié)果表明加入Sentinel-2 遙感變量后土壤有效磷的預(yù)測(cè)精度顯著提高(模型A:=0.59),而且紅邊波段是有效預(yù)測(cè)土壤有效磷含量空間分布的重要變量之一。相關(guān)研究也證明引入高空間分辨率的衛(wèi)星遙感數(shù)據(jù)可提高土壤屬性預(yù)測(cè)精度。Gholizadeh 等驗(yàn)證了Sentinel-2 遙感數(shù)據(jù)提高耕地土壤屬性的制圖精度,紅邊波段也被證實(shí)可以提高土壤屬性預(yù)測(cè)精度。在未來(lái)的數(shù)字土壤制圖研究中,可以考慮使用更多的紅邊波段和更高空間分辨率的遙感數(shù)據(jù)。

在土壤樣點(diǎn)有限的情況下,利用Sentinel-2 數(shù)據(jù)結(jié)合其他環(huán)境數(shù)據(jù)作為輔助變量建立預(yù)測(cè)福建省建甌市土壤有效磷含量的RF 模型,模型的穩(wěn)定好,預(yù)測(cè)精度較高,說(shuō)明RF 在土壤屬性制圖研究中受采樣點(diǎn)數(shù)量的影響較小,這與劉明杰等在稀疏樣本下基于RF 方法實(shí)現(xiàn)土壤有機(jī)質(zhì)預(yù)測(cè)的結(jié)論一致。所有模型訓(xùn)練集與驗(yàn)證集基本一致,進(jìn)一步驗(yàn)證了隨機(jī)森林模型不易過(guò)擬合的優(yōu)勢(shì)。在土壤有效磷含量的空間分布預(yù)測(cè)中也取得較好的效果,表明RF模型可用于稀疏采樣的土壤有效磷含量的預(yù)測(cè)制圖。

本研究基于南方典型丘陵地區(qū)稀疏的土壤樣點(diǎn)數(shù)據(jù),結(jié)合全部變量(包括Sentinel-2 遙感數(shù)據(jù)、氣象變量、地形變量和土壤pH)建立的RF 模型A 在所有變量組合模型中取得最佳的預(yù)測(cè)效果。Wang等也證明基于全部環(huán)境變量的RF 模型預(yù)測(cè)土壤有機(jī)碳比其他組合模型效果好。徐劍波等預(yù)測(cè)南方丘陵縣域(廣東省高州市)土壤有效磷含量,采樣點(diǎn)為100 時(shí),MAE 為22.66 mg?kg,RMSE 為30.97 mg?kg。與其相比本研究加入輔助環(huán)境變量后基于稀疏樣點(diǎn)土壤有效磷預(yù)測(cè)精度提升明顯(MAE=19.04 mg?kg,RMSE=25.26 mg?kg)。不包含遙感變量的模型C 預(yù)測(cè)精度最低(驗(yàn)證集=0.36,MAE=22.66 mg?kg,RMSE=31.71 mg?kg),加入遙感變量后模型A 預(yù)測(cè)精度明顯提高(從0.36 提升至0.59),MAE 和RMSE 明顯下降(RMSE 降低了20.34%),進(jìn)一步證明了遙感光譜信息在土壤屬性預(yù)測(cè)的重要性,這與許多學(xué)者研究結(jié)論一致。本研究結(jié)果表明Sentinel-2 遙感數(shù)據(jù)結(jié)合其他環(huán)境變量作為二次數(shù)據(jù)源可以實(shí)現(xiàn)稀疏樣點(diǎn)的土壤屬性制圖,這與de Carvalho 等研究結(jié)論一致。

3.2 變量相對(duì)重要性討論

篩選的最佳RF 模型(模型A)變量相對(duì)重要性結(jié)果顯示年均溫、土壤pH、地形濕度指數(shù)和高程是影響南方丘陵地區(qū)福建省建甌市稀疏樣點(diǎn)土壤有效磷含量空間分布的最重要因子,它們的重要性分值累計(jì)達(dá)68.66%。研究區(qū)土壤有效磷含量空間分布基本上可以用遙感變量(22.87%)、氣象變量(30.64%)和地形變量(30.38%)解釋。遙感波段反射率和光譜信息可以反映地表的土壤和植被生長(zhǎng)狀況,從而可以直接或間接地反映土壤物理化學(xué)性質(zhì); 作為二次數(shù)據(jù)源使用還可以彌補(bǔ)樣點(diǎn)數(shù)量的不足,提高制圖精度。Shen 等研究表明在平坦地區(qū)遙感影像提取的亮度指數(shù)為土壤全磷的最優(yōu)預(yù)測(cè)變量之一,而本研究中3 個(gè)遙感變量相對(duì)重要性均偏低,這可能是因?yàn)槟戏角鹆甑貐^(qū)耕地田塊小,圖斑破碎,一個(gè)土壤樣點(diǎn)對(duì)應(yīng)的光譜信息為不同田塊的混合光譜信息,造成研究誤差進(jìn)而影響建模中的重要性。雖然遙感變量增強(qiáng)植被指數(shù)、第一主成分和紅邊波段單因子重要性較低,但累計(jì)分值達(dá)22.87%,仍是影響土壤有效磷含量空間分布的重要變量。3 個(gè)因子中,增強(qiáng)植被指數(shù)和第一主成分相對(duì)重要性略高于B6 紅邊波段,這是因?yàn)樵鰪?qiáng)植被指數(shù)等相對(duì)于單波段一定程度可以減少大氣等外部干擾因素,提高土壤屬性的預(yù)測(cè)精度。此外,紅邊波段是綠色植物生長(zhǎng)狀況的敏感性波段,可間接地反映土壤養(yǎng)分情況,研究表明紅邊波段是土壤有機(jī)質(zhì)的最佳預(yù)測(cè)波段之一,本研究中B6 紅邊波段重要性分值僅有6.88%,未達(dá)到預(yù)期效果,但其制圖潛力值得進(jìn)一步挖掘。

氣象變量和地形變量是南方丘陵地區(qū)土壤有效磷含量空間分布最重要的影響因素,這與楊之江等研究結(jié)論一致。氣象變量中年均溫(22.17%)是影響建甌市土壤有效磷含量空間分布的最重要因子,詹秋麗等也認(rèn)為年均溫與福建省土壤有效磷含量分布顯著相關(guān)。因?yàn)榍鹆甑貐^(qū)的地形起伏大,受水平及垂直地帶差異影響年均溫差異明顯(約7.68 ℃),年均溫可能是主要通過(guò)影響成土過(guò)程中含磷礦物風(fēng)化速率及土壤中微生物活性等間接影響土壤磷素的有效性,從而影響土壤有效磷的含量。而年均降水量(8.47%)主要通過(guò)土壤淋溶作用影響土壤有效磷含量的空間分布,年均降水量越多,土壤淋溶作用越強(qiáng),導(dǎo)致土壤中有效磷的流失。地形變量主要通過(guò)支配地表徑流影響土壤磷流失及通過(guò)影響水熱再分配兩方面對(duì)土壤有效磷含量的空間分布產(chǎn)生影響。地形變量中地形濕度指數(shù)反映了地表徑流和土壤水分含量狀況,地勢(shì)低的地方,地表徑流流失的有效磷含量越少,土壤磷素更容易積累,研究表明地形濕度指數(shù)為土壤有效磷預(yù)測(cè)模型中的重要影響因子(15.51%),與前人研究結(jié)論一致。

僅考慮單一因子,土壤pH (16.11%)是除了年均溫外最重要的因子,與王燕等在江蘇省金壇區(qū)土壤有效磷含量的預(yù)測(cè)研究結(jié)論相吻合。研究表明土壤pH 是通過(guò)影響土壤磷素吸附固定作用從而影響土壤磷素有效性。土壤pH 越低,酸性越強(qiáng),磷素被鐵、鋁氧化物吸附固定作用越大,以磷酸鹽的形態(tài)存在,此外農(nóng)業(yè)管理中使用的磷肥也易被大量吸附固定,從而增加土壤中有效磷的含量; 反之,pH 越大,土壤呈堿性,磷素易被土壤中的鈣吸附固定而大多數(shù)以難溶解的無(wú)機(jī)磷的形態(tài)存在,從而導(dǎo)致有效磷的含量減少。

3.3 土壤有效磷含量的空間分布

典型南方丘陵地區(qū)福建省建甌市土壤有效磷含量總體較為豐富,呈現(xiàn)中西高、東南低的分布趨勢(shì)。高值區(qū)主要分布在年均降水量低,且年均溫高的中西部低海拔區(qū)域,與其他亞熱帶丘陵地區(qū)土壤有效磷含量空間分布研究規(guī)律一致。中西部河谷平原與盆地地區(qū),年均溫較高加快了礦物風(fēng)化速率和脫硅富鐵鋁作用,且地形相對(duì)平坦,年均降水量較少,施肥中磷素不容易淋失,故有效磷含量較高; 東南部地貌以山地為主,降水量較多,通過(guò)地表徑流加劇了土壤磷素的流失,而且年均溫較低限制土壤母質(zhì)含磷礦物的風(fēng)化速率導(dǎo)致東南部山地土壤的有效磷含量較低。此外,超過(guò)一半以上的地區(qū)土壤有效磷含量高于40 mg?kg,超過(guò)該閾值則易增加從土壤中淋失到地下水并進(jìn)一步遷移至地表水體中的磷,不利于農(nóng)業(yè)可持續(xù)發(fā)展與環(huán)境保護(hù)。因此,在未來(lái)的農(nóng)業(yè)施肥管理中應(yīng)根據(jù)不同區(qū)域精準(zhǔn)制定合理的磷肥施用量。

4 結(jié)論

本研究基于實(shí)測(cè)稀疏土壤樣點(diǎn),以Sentinel-2 遙感數(shù)據(jù)、氣象、地形和土壤pH 作為輸入變量,利用皮爾遜相關(guān)性分析及后向剔除法篩選參與建模的環(huán)境因子,比較不同環(huán)境變量組合的RF 模型預(yù)測(cè)效果,篩選最佳的環(huán)境變量組合預(yù)測(cè)土壤有效磷含量(驗(yàn)證集=0.59,MAE=19.04 mg?kg,RMSE=25.26 mg?kg),繪制典型南方丘陵地區(qū)福建省建甌市耕地表層(0～20 cm)的土壤有效磷含量空間分布圖。結(jié)果顯示增加Sentinel-2 遙感變量能顯著提高模型的擬合效果和預(yù)測(cè)精度,利用Sentinel-2 遙感數(shù)據(jù)聯(lián)合氣象、地形和土壤pH 的全部變量的RF 模型預(yù)測(cè)精度最高,可作為稀疏采樣點(diǎn)下南方丘陵地區(qū)耕地表層土壤有效磷含量的有效預(yù)測(cè)方法。氣象變量和地形變量是解釋典型南方丘陵地區(qū)福建省建甌市耕地土壤有效磷空間分布的主要變量,遙感變量是影響其土壤有效磷分布的重要變量。本研究可為快速而準(zhǔn)確地繪制其他縣市尺度的表層土壤屬性空間分布圖提供借鑒。本研究為基于稀疏樣點(diǎn)的南方丘陵地區(qū)土壤有效磷制圖提供了方法上的嘗試,但也存在一定的應(yīng)用局限。在建模因子的選擇方面,受數(shù)據(jù)的可獲取性限制本研究探討的輔助變量有限,在未來(lái)研究中應(yīng)該考慮在模型中加入施肥等人為管理措施、成土母質(zhì)等因素,進(jìn)一步提高制圖精度。此外,在數(shù)據(jù)預(yù)處理方面,為了統(tǒng)一輔助數(shù)據(jù)的空間分辨率與土壤有效磷樣點(diǎn)數(shù)據(jù)進(jìn)行尺度匹配,將氣象數(shù)據(jù)和地形數(shù)據(jù)重采樣至10 m 空間分辨率,但這可能會(huì)帶來(lái)一定的數(shù)據(jù)誤差,尤其是目前可獲取的氣象數(shù)據(jù)相對(duì)粗糙,如何進(jìn)行輔助數(shù)據(jù)降尺度以更好地滿足空間匹配要求值得進(jìn)一步探究。