基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演.pdf

基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演聂小军1，赵星辉2，AMMARAGill1，洪雯雯1，张合兵1(1.河南理工大学测绘与国土信息工程学院,河南焦作454003；2.河南农业大学林学院,河南郑州450002)R2m摘要：煤颗粒作为一种高碳有机质，即使少量扩散至土壤中便会以煤源碳的干扰形式导致土壤有机碳(植物源碳)质量分数的明显高估，从而增加土壤碳固存评估的不确定性。然而，目前缺少土壤煤源碳定量评估的有效方法。本研究以具有 100 多年无烟煤开采历史、土壤中煤颗粒累积(煤累积土壤)普遍的河南焦作矿区为研究区，采集当地的煤粉与不含煤颗粒的土壤，配置煤源碳质量分数不同的土壤样品，利用高光谱遥感技术分析了矿区土壤的光谱特征。综合 8 种光谱数学变换方法和 2 种光谱特征筛选方法，构建了 6 种矿区土壤煤源碳质量分数反演模型，包括 3 种光谱特征指数(弓曲差、差值指数、比值指数)、3 种机器学习(偏最小二乘回归(PLSR)、支持向量机(SVM)与随机森林(RF)模型。本研究也对最优反演模型的适用性进行了检验。研究发现：波长 3502500nm内，煤的光谱曲线特征明显不同于植物源有机质与不含煤颗粒的土壤，煤累积土壤的光谱反射率(R)随煤源碳质量分数的增加而减小，这为高光谱遥感技术定量反演土壤煤源碳提供了理论基础。在光谱特征筛选方面，综合竞争性自适应重加权采样法(CARS)提取出的煤源碳特征波段数远高于光谱特征指数相关系数法，而且特征波段在波长 3502500nm 内分布均匀。原始光谱 R 经光谱数学变换后，构建的弓曲差、差值、比值光谱特征指数反演模型对土壤煤源碳质量分数的估测精度明显提升，其中，基于倒数 1/R 变换的差值指数模型反演效果最好。相较于光谱特征指数模型，结合 CARS 的机器学习模型对煤源碳质量分数的估测精度进一步提升。3 种机器学习模型中，1/R-CARS-RF 煤源碳质量分数反演模型的估测精度最高，验证集=0.998、RMSE=0.348、RPD=29.943。适 用 性 检 验 表 明，1/R-CARS-RF 煤 源 碳 质 量 分 数 反 演 模 型 的 适 用 性 良 好(RMSE=1.88%、RPD=4.97)，可以较准确地估测焦作矿区土壤中的煤源碳质量分数。本研究可为矿区土壤煤源与植物源碳区分、土壤碳固存精确评估提供方法支撑。关键词：煤源碳；土壤有机碳；高光谱反演；碳固存；煤矿区中图分类号：TD88文献标志码：A文章编号：02539993(2023)07286912Inversion of coal-derived carbon content in mine soils based on hyperspectralindex and machine learningNIEXiaojun1,ZHAOXinghui2,AMMARAGill1,HONGWenwen1,ZHANGHebing1(1.School of Surveying and Land lnformation Engineering,Henan Polytechnic University,Jiaozuo454003,China;2.College of Forestry,Henan Agricul-tural University,Zhengzhou450002,China)收稿日期：20230315修回日期：20230428责任编辑：王凡DOI：10.13225/ki.jccs.CN23.0329基金项目：国家自然科学基金资助项目(41977284)；河南省科技攻关资助项目(222102320032)作者简介：聂小军(1977)，男，山西曲沃人，教授，博士生导师。E-mail：通讯作者：张合兵(1975)，男，河南滑县人，教授，博士生导师。E-mail：引用格式：聂小军，赵星辉，AMMARAGill，等.基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演J.煤炭学报，2023，48(7)：28692880.NIEXiaojun，ZHAOXinghui，AMMARAGill，etal.Inversionofcoal-derivedcarboncontentinminesoilsbasedonhyperspectralindexandmachinelearningJ.JournalofChinaCoalSociety，2023，48(7)：28692880.第48卷第7期煤炭学报Vol.48No.72023年7月JOURNALOFCHINACOALSOCIETYJul.2023R2mAbstract:Coalparticles,atypicaloforganicmatterwithhighcarboncontent,widelydiffusedinsoilenvironmentduetomassiveconsumptionofcoalenergyforcenturies.Eveninsmallamountspresentedinsoils,coalparticlescanproduceob-viousoverestimationofsoilorganiccarbon(SOC)andthusincreaseanuncertaintyofsoilCsequestrationassessment.However,thereisalackofdeterminationmethodsfordistinguishingcoal-derivedCfromSOC.ThisstudytakestheJiaozuominingareaasthestudyarea,whichhasahistoryofanthraciteminingformorethan100yearsandwherecoal-contaminatedsoilsarewidespread.Thehyperspectralcharacteristicsoftheminesoilswereanalyzedbycollectingcoalparticlesandcoal-freesoilsandthenmixingknownquantitiesofthetwosamplesmanually.Basedoneightspectralmath-ematicaltransformationandtwospectralfeaturescreeningmethods(i.e.,traditionalcorrelationcoefficiencies,andcom-prehensivecompetitiveadaptivereweightedsampling(CARS),sixinversionmodelsincludingspectralfeatureindexmodels(i.e.,deviationofarch(IDOA),differenceindex(ID),ratioindex(IR),andthreemachinelearningmodels(i.e.,par-tialleastsquaresregression(PLSR),supportvectormachine(SVM)andrandomforest(RF),wereestablishedtoinvers-ingcoal-derivedCcontent.Theapplicabilityoftheestablishedoptimalinversionmodelwasalsoexamined.Itwasfoundthatinthewavelengthrangefrom350nmto2500nm,thespectralcurvesofcoalparticlesareobviouslydifferentfromthoseofplant-derivedorganicmatterandcoal-freesoils.Moreover,thespectralreflectance(R)ofcoal-contaminatedsoilsdecreaseswithincreasingcoal-derivedCcontent.Thesefindingsprovideatheoreticalbasisfortheapplicationofhyper-spectralremotesensingtechnologytoquantitativelyinversecoal-derivedCofminesoils.TheCARS-extractedfeaturewavebandsofcoal-derivedCdistributedevenlyinthewavelengthrangebetween350nmand2500nm,andthenumberoffeaturewavebandextractedbytheCARSwasfarhigherthanthatbytraditionalcorrelationcoefficiencymethods.AftermathematicaltransformationofRdata,theestimatedaccuraciesofcoal-derivedCcontentproducedbytheIDOA,IDandIRinversionmodelsweresignificantlyenhanced,amongwhich,theDImodelof1/Rspectraltransformationexhibitedthehighestestimationaccuracy.Comparedwiththetraditionalindexmodels,thethreemachinelearningmodelscombinedwithCARSfurtherenhancedestimationaccuracyofcoal-derivedCcontent.Amongthethreemachinelearningmodels,the 1/R-CARS-RF inversion model produced the highest estimation accuracy,showing=0.998,RMSE=0.348,RPD=29.943foritsvalidationset.Theapplicabilitytestshowedthatthe1/R-CARS-RFmodelexhibitedagoodapplicabil-ityfordifferentcoal-contaminatedsoilsintheJiaozuominingareaasobserved1.88%forRMSEand4.97forRPD.There-fore,itcanbeexpectedthathyperspectralremotesensingtechnologyhasapromisingapplicationprospectfordeterminingcoal-derivedCcontentinminesoils.Inaddition,thisstudycanprovideamethodologysupportfordistinguishingcoal-de-rivedCfromSOCinminesoilsandaccuratelyassessingsoilcarbonsequestration.Key words:coal-derviedcarbon；soilorganiccarbon；hyperspectralinversion；carbonsequestration；coalminingarea土壤是陆地生态系统最大的碳库，全球约 1500亿 t碳以有机碳形式储存于陆表土壤，其碳库储量约是大气碳库的 2 倍，即使轻微的土壤有机碳质量分数变化，都会导致大气圈中二氧化碳体积分数的剧烈波动1-2。双碳背景下，土壤有机碳的高效固存是21世纪我国耕地质量提升和生态环境保护的重大需求，也是国际社会需要长期面对的重大科学命题3。准确估算土壤有机碳储量，对全球气候变化和农业可持续发展等具有重要意义。全球煤炭资源分布广泛，近 250 多年的煤炭规模使用历史导致煤颗粒通过不同方式广泛沉积至陆地土壤中，形成了大面积的黑色土壤。但是，这种黑色土壤并不是富含植物源碳(plant-derivedOC)的肥沃土壤，而是富积煤源碳(化石燃料碳 fossilfuelC 或地质成因有机碳geogenicOC)的特殊黑土。千年尺度范围内的现代植物被认为是土壤有机碳的主要来源4-5。不同于植物源有机碳(植物源碳)，煤源碳成因于久远(千万年尺度以上)的地质过程，化学结构以高分子芳香环为主，元素组成中缺少植物与土壤微生物所需的氮磷钾等矿质养分，稳定性强，生物极难分解利用。因此，土壤中的煤源碳并不是真正的土壤有机碳6。研究报道：煤源碳的存在不仅会造成煤矿区土壤有机碳的严重高估(0.27.4 倍)，而且还会导致植物源有机碳的快速分解，结果使得煤矿区土壤成为一个明显的碳源7-8。我国是煤炭资源大国，大多煤矿区煤粮复合特征明显，农业土壤中的煤粒累积极为普遍。土壤中煤粒的存在会干扰土壤有机碳的识别。如何精准评估煤矿区土壤有机碳质量分数是矿区土壤碳固存与实现“碳达峰、碳中和”双碳目标亟需解决的一个科学问题。2870煤炭学报2023年第48卷目前，世界各国通用的土壤有机碳测定方法包括:重铬酸钾氧化法(WB135法、WB137法)与干烧法9。但是，这些方法都属于化学方法，在测定煤累积土壤时，它们均会不同程度地氧化煤源碳，从而高估土壤有机碳质量分数，且检测过程复杂、污染环境，因此并不适合用来评估煤矿区土壤有机碳的真实固存7。随着同位素技术的发展，稳定性(13C)与放射性(14C)碳同位素示踪技术已经被成功地应用于德国、美国一些煤矿区土壤的有机碳的源识别10-11。但是，13C 示踪技术的应用必须基于一个前提：土壤中的煤源碳与植物源碳的 13C 存在显著的差别，是可以被区分的，然而事实上很多矿区土壤中这 2 种源的有机碳 13C 值是很接近的，难以被分清，因此，该技术的应用具有很大的局限性10,12。14C 示踪技术分析使用的仪器昂贵、国内外仅有极少数高校与科研院所拥有该仪器且测试费用极高，普及性低。综上，目前常规的土壤有机碳测定方法与碳同位素示踪技术从适用性或普及性方面都不能很好地服务于煤矿区土壤有机碳准确定量与土壤碳固存评估。因此，亟需开展科学有效的方法研究来定量识别煤矿区土壤中的煤源碳。高光谱遥感技术以其快速经济、安全、实时、光谱信息丰富等优点成为反演土壤属性信息的一种新方法13。目前，许多学者利用高光谱遥感技术开展了土壤有机碳的反演估算研究，展现了良好的应用前景14-16。有关煤矿区土壤有机碳的高光谱反演研究也被报道。南锋等17利用偏最小二乘回归构建了黄土高原煤矿区有机碳质量分数估测模型。陈玮等18估测了煤炭矿区周边不同沉陷阶段的耕地土壤有机碳。ZHU 等19使用线性模型和随机森林对新疆省准东煤田矿区土壤有机碳进行了反演。WEI 等20分别构建了中国西北的农田区和矿区土壤有机碳质量分数反演模型，并发现不同土壤样品对建模效果有显著影响。但是，目前的矿区土壤有机碳光谱反演研究所得出的反演精度较低，这可能是由于忽略了土壤中煤源碳的存在，煤源碳干扰了土壤有机碳的光谱信息。因此，如何利用高光谱遥感技术实现矿区土壤煤源碳的定量识别值得关注。基于以上背景，笔者以河南省焦作矿区为研究区，利用高光谱遥感技术，通过配置煤源碳质量分数不同的土壤样品，室内采集样品的光谱数据，分析煤累积土壤的高光谱特征，构建基于光谱特征指数与机器学习方法的煤源碳反演模型并进行验证，为定量识别矿区土壤煤源碳与准确评估土壤碳固存本底提供技术支撑。1材料与方法1.1研究区概况焦作矿区(39.57N，113.70E)地处河南省太行山山前平原冲洪积扇地带，属温带大陆性季风气候，平均气温 15，年降水量 600mm，大约 70%的降水发生在 69 月21。土壤类型为石灰性褐土(pH 值为8.1，CaCO3质量分数为 14%)，母质为第四系黄土和灰岩、砂岩风化物，质地为粉质壤土(砂粒粉粒黏粒为 23%67%10%)22。矿区光热资源丰富，煤粮复合特征典型，是中国粮食高产区之一。农田种植模式主要为夏玉米冬小麦。矿区 110 多年的煤矿开采、分选、运输、矸石堆积、矿井水灌溉等活动，导致煤粒广泛沉积在该地区的土壤中，煤累积土壤普遍。1.2样品采集与配置由于目前没有可靠的方法获取研究区土壤中的煤源碳质量分数，因此，笔者采用室内样品配置方法来开展土壤煤源碳质量分数的光谱反演研究。通过田间调查、土壤颜色(黑度)识别(卡兰德 EDA-100X，中国)与便携式光学显微镜与场发射扫描电镜(Mer-linCompact，德国)检测，采集矿区范围内不含煤粒的耕地表土(厚度 020cm)样品 175 个(图 1)。随机挑选 100 个样品等量均匀混合，形成混合样品 50kg，以此用来光谱建模研究。其余 75 个表土样品不进行混合，用来检验模型的适用性效果。实地调查也发现，矿区存在个别地块利用南水北调中线工程挖掘的第四系深层古红土来进行复垦，为此，增加采集该类型土壤样品 25 个(图 1)，用来检验模型的适用性效果。此外，根据实地调查与测试分析发现，研究区煤累积土壤中煤粒的粒径0.25mm、有机碳质量分数60%(600g/kg)，这与当地煤矿煤炭分选过程中产生的煤粉特征一致。为此，采集当地煤矿煤炭分选过程中产生的煤粉样品 20kg，以此来配置煤源碳质量分数不Nbbbbbbbb!P修武县焦作048 kmb!耕地样点!复垦红土样点矿区边界煤矿耕地图1采样点分布Fig.1Distributionofsamplingpoints第7期聂小军等：基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演2871同的土壤样品。同步采集当地农田土壤有机碳的源样品(即植物源有机质)：猪厩肥、牛厩肥、鸡厩肥、羊厩肥 4 种农家有机肥，2 种不同腐化程度的农田秸秆。此外，为更好地阐明煤源碳的光谱特征，也采集了当地煤矿开采的新鲜无烟煤。所有采集的样品带回实验室进行风干研磨过筛(0.25mm)烘干(65，96h)系列处理。为配置煤源碳质量分数不同的土壤样品，在实验室内对采集的煤粉、新鲜无烟煤以及不含煤的土壤样品进行酸化处理，然后利用 VarioTOC 仪(Elemental德国)测定其有机碳质量分数。经测定，煤粉与不含煤土壤的有机碳质量分数分别为 57.57%(575.7g/kg)、1.43%(14.3g/kg)；新鲜无烟煤煤块制成粉末样品后测定的有机碳质量分数为 88.48%(884.8g/kg)。前期研究发现，研究区大多数农田煤积累土壤的总有机碳质量分数不超过 14%，极少数高达 27%，而不含煤颗粒土壤的有机碳质量分数不超过 2%7。由此估算，研究区煤积累土壤的煤源碳质量分数主要在12%以下，最高不超过 25%。为此，根据供试煤粉的有机碳质量分数，按事先计算好的(煤+土)混合比例设置煤源碳质量分数范围：032%，总计 83 个，以此用于高光谱建模。为保证建模过程中光谱测试的准确性，每个质量分数样品配置 3 个重复，总计配置了249 个训练样品(训练集+验证集)。用于模型适用性检验的 100 个(75+25，见上文)土壤样品同样按照032%煤源碳质量分数来配置，以此作为测试集样品。所有样品都进行光谱测定，随后 249 个建模样品依次进行数据变换建模变量优选反演模型构建，100 个测试集样品用来检验反演模型的适用性。1.3光谱测定每个样品取 40g 放置于密闭的盛样器(直径5cm、深 4cm)中摇匀，然后待测光谱反射率。光谱反射率测定采用美国 ASD 公司生产的 FieldSpec4 便携式地物光谱仪。该仪器的光谱分辨率为 3nm，光谱检测范围 3502500nm，在 3501000nm 采样间隔为 1.4nm，在 10002500nm 采样间隔 2nm，重采样间隔 1nm23。该仪器提供了一种高密度接触式反射探头，可以直接接触样品表面进行光谱采集，不仅可以解决天气状况对光谱测试的影响，也可以随时创造暗室条件，有效避免杂光影响24。试验在暗室中进行，刮平盛样器内的土壤表面，利用高密度反射探头获取光谱反射信息。每个样品重复测 10 次，取其算术平均值作为该样品的反射率光谱值。测试之前进行白板校正，每隔 15min 重复 1 次白板校正。1.4光谱数据处理由于在 350400nm 波段与 24502500nm 波段反射率的噪声较大，因此选取 4002450nm 光谱数据进行数据分析25。为了得到平稳的光谱反射率，笔者采用 S-G平滑法(多项式平滑法)进行平滑波谱处理，去除波段中的少量噪声。采用的光谱数据有光谱反射率(R)、倒数(1/R)、去包络线(RC)、平方(R2)、开根号(R0.5)、二阶微分(R)、一阶微分(R)、对数(lnR)、倒数对数(ln(1/R)、对数倒数(1/lnR)。1.5建模变量优选光谱特征指数相关系数法是对光谱波段进行弓曲差、差值、比值等计算后，进行相关分析来优选建模变量的一种方法26。由于光谱仪每次光谱采集存在误差，使用光谱特征指数相关系数方法提取光谱特征，比单一波段所构建的模型具有更高的泛化力，应用性更强27。本研究通过检验两波段弓曲差(IDOA)、差值指数(ID)、比值指数(IR)与煤源碳质量分数之间的相关性来筛选最大相关系数，进而获取最大相关系数对应的波段组合来实现光谱特征波段的提取。3 种光谱特征指数的计算公式为IDOA=R(i+j)/2(Ri+Rj)/2(1)ID=RiRj(2)IR=RiRj(3)R(i+j)/2其中，Ri和 Rj分别为 i 和 j 波段的土壤光谱反射率；为 i 和 j 两个波段的光谱反射率的平均值。可见光范围内的弓曲差已经广泛被用来优选土壤有机碳光谱特征28-29。笔者通过计算 450800nm 波段范围内弓曲差与土壤煤源碳质量分数间的相关性，来优选基于弓曲差光谱特征指数的建模变量。竞争性自适应重加权采样法(CARS)是一种结合蒙特卡洛采样与偏最小二乘回归系数的特征变量选择方法，已经广泛用于光谱特征筛选。该方法借鉴了自然界生物遗传机制，模仿达尔文理论中的“适者生存”的原则，有助于建立对外界因素的鲁棒性分析模型，减少谱矩阵的信息冗余和多重共线性问题30-31。波长变量优选过程中，首先通过蒙特卡罗采样，建立偏最小二乘回归(PLSR)模型，选取回归系数较大的波长变量，剔除回归系数较小的波长变量，从而获得多个优选变量。再通过 K 折交叉验证的方式计算模型 PLSR 中的交叉验证均方根误差(RMSECV)最小的优选变量，即为最优变量。笔者使用 10 折交叉验证，设置蒙特卡罗采样次数为 50次，对光谱数据循环迭代 1000 次，进行光谱特征波段筛选。2872煤炭学报2023年第48卷1.6模型构建与精度检验本研究通过随机划分法，将 249 个土壤样品的70%作为训练集用于建模，30%用于模型验证。采用3 种光谱特征指数(IDOA、ID、IR)作为建模变量，基于多元线性回归进行煤源碳质量分数建模。采用CARS 算法筛选出特征波段后，基于偏最小二乘回归(PLSR)、支持向量机(SVM)与随机森林(RF)3 种常用机器学习方法来构建煤源碳质量分数反演模型。PLSR 是基于因子提取的多元回归分析方法，其提供一种多对多线性回归建模的方法，特别当 2 组变量的个数很多，且都存在多重相关性，而观测数据的数量较少时，用该模型优于传统的经典回归分析等方法32。SVM 是一种基于统计学 VC 维理论和结构风险最小原理，针对二分类问题提出的一种机器学习方法，该方法可以进行线性回归，也可以通过核函数的作用解决非线性回归33。RF 是通过集成学习的思想将多棵树集成的一种算法，它能够处理高维度数据，并且不用做特征选择，对数据集的适应能力强，既能处理离散型数据，也能处理连续型数据，训练速度快，实现简单，精度高、抗过拟合能力强34。本研究中的SVM 模型，通过多次调试，使用高斯核函数(rbf)时建模效果优于其他核函数，RF 使用 K 折交叉验证经过多次调试选取最优模型。以上模型的构建均在Anaconda软件中使用 python3.7 语言编写，通过调用sklearn0.21.3 库实现。R2mR2m模型精度通过决定系数()、均方根误差(RMSE)以及相对分析误差(RPD)3 个指标进行评价。越接近于 1，表明模型拟合度越高。RMSE越小以及 RPD越大，说明模型的估测效果和稳定性越好。当 RPD3 时，表明模型估测能力极好35。2结果与分析2.1光谱曲线特征由图 2(a)可以看出，煤与植物源有机质、不含煤粒土壤的光谱曲线特征明显不同。煤样品(煤粉与新鲜无烟煤)的光谱曲线在全波段范围内(3502500nm)都比较平缓，反射率低(0.72)；其中，600nm 与 750nm两处的相关性最高(|r|0.90)，表现出极好的负相关性(表 1)。图 3 与图 4 为光谱特征指数与煤源碳质量分数的相关系数等势图，不同颜色表示两波段组合计算的ID、IR光谱特征指数与煤源碳质量分数的相关系数绝对值。相关系数等势图清晰展示出相关系数较高的特征指数构成波段的分布范围。原始光谱数据中，煤源碳质量分数与 ID、IR的高相关组合波段分布在第7期聂小军等：基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演2873400500nm，其中，差值光谱特征指数相关系数法提取的最优组合波段为 406nm&454nm，比值光谱特征指数相关系数法提取的最优组合波段为 402nm&1003nm，对应的|r|分别为 0.952 与 0.974。原始光谱数据经过不同的数学变换处理后(除微分变换外)，增强了光谱特征指数与煤源碳质量分数的相关性。光谱变换下的 ID、IR与煤源碳质量分数高相关组合波段为422nm&454nm、423nm&454nm、402nm&1003nm、表 1 可见光部分弓曲差与煤源碳质量分数的相关性(n=249)Table 1 Correlation analysis between bow difference and coal source carbon content in visible light waveband(n=249)弓曲差/nm450500550600650700750800相关系数0.0970.8730.5170.9050.7270.7670.9040.824注：n为样品的数量。2 4001 2001 6002 000800波长/nm1.00.80.60.40.20RC1/Rln RR2 4004001 2001 6002 0008004008001 6002 4002 4002 4002 4002 4008001 600ln(1/R)8001 6001/ln RRR8001 600波长/nm相关系数|r|R28001 600R0.5图3煤源碳质量分数与 ID相关系数等势图Fig.3Equipotentialdiagramofcorrelationcoefficientbetweencoal-derivedcarboncontentandIDvalue2 4001 2001 6002 000800400波长/nm2 4001 2001 6002 0008004008001 6002 400 8001 6002 400 8001 6002 400波长/nm相关系数|r|8001 6002 400 8001 6002 4001.00.80.60.40.20RC1/Rln RRln(1/R)1/ln RRRR2R0.5图4煤源碳质量分数与 IR相关系数等势图Fig.4Equipotentialdiagramofcorrelationcoefficientbetweencoal-derivedcarboncontentandIR2874煤炭学报2023年第48卷2184nm&2318nm、2127nm&2323nm 等。其中，最优数学变换为 1/R 变换，该变换下 ID、IR与煤源碳质量分数的|r|分别为 0.991 与 0.988，对应的组合波段为 2127nm&2323nm 与 402nm&1003nm，高于原始光谱下的|r|最大值(0.974)。CARS 算法对煤源碳特征波段的提取如图 5 所示。随着采样数的增加，选出的特征波段数逐渐减少。在采样时，RMSECV不断减小，表明与煤源碳无关的变量逐渐被剔除；当采样次数达到 16 次以后，RMSECV上升，这可能是剔除了与煤源碳有关的变量而导致了 RMSECV增大(图 5(a)。原始光谱与 9 种光谱变换数据中，CARS 筛选出的煤源碳特征波段在全波段范围内分布较为均匀，特征波段为 53214 个(图 5(b)，远高于 3 种光谱特征指数提取的高相关组合波段数。4008001 2001 600波长/nm2 0002 400504030201000.010.020.030.041 0001 000变量数量/个RMSECVPLSR 回归系数路径采样次数(b)(a)1/Rln(1/R)R0.5R2ln R1/(ln R)RCRRR01 0002 000图5CARS 算法提取煤源碳特征波段Fig.5Featurewavebandsofcoal-derivedcarbonextractedbyCARSalgorithm2.3基于光谱特征指数的煤源碳质量分数反演模型表 2 为基于 3 种光谱特征指数构建的煤源碳质量分数反演模型及其效果。由表 2 可知，基于 600、750nm 处弓曲差的反演模型都可以用一元二次回归模型来反演土壤煤源碳质量分数，反演精度较高(5.9RPD13)、lnR(RPD21)变换下具有较好的建模效果。综合 3种机器学习建模结果，基于 1/R 变换的 RF、PLSR 模型以及基于 lnR 变换的 SVM 模型精度最好，尤其是1/R-CARS-RF 模型反演煤源碳质量分数的精度最高，验证集达到 0.998，RMSE 为 0.348%，RPD 高达29.943(表 3)。2.5模型适用性测试用于模型测试的 100 个不含煤颗粒的原始土壤样品中，75 个为矿区分布最广泛的土壤，母质为第四系黄土和灰岩、砂岩风化物(矿区土壤)；其余 25 个为矿区复垦土壤，母质为第四系深层古红土(复垦红土)。2 种类型土壤的光谱曲线特性存在差异，主要表现为：复垦红土在 750nm 附近反射率迅速增强，12002500nm 光谱反射率远高于矿区土壤(图 7(a)。通过向采集的不含煤颗粒的原始土壤样品中添加已知量的煤粉(煤源碳质量分数 575.7g/kg)，配置煤源碳质量分数不同的矿区土壤、复垦红土，以此为光谱测试对象，检验精度最高的 1/R-CARS-RF 煤源碳质量分数反演模型的适用性。经测试发现：对于矿区土壤，1/R-CARS-RF 煤源碳质量分数反演模型的第7期聂小军等：基于光谱特征指数与机器学习的矿区土壤煤源碳质量分数反演2875R2mR2m=0.94、RPD=4.97；对于复垦红土，该模型的=0.64、RPD=1.70。此外，由图 7(b)可以看出，2 种类型土壤的煤源碳质量分数估测值与真实值极显著相关(r0.98，P0.001)，但是，复垦红土的 RMSE(4.79%)远大于矿区土壤(1.88%)。进一步分析发现：在 020%范围内，1/R-CARS-RF 模型会轻微低估矿区土壤的煤源碳质量分数，而在大于 22%时高估；而复垦红土在煤源碳质量分数小于 10%时，模型会产生轻微高估，大于 10%时明显高估。总体上，1/R-CARS-RF 煤源碳质量分数反演模型的适用性良好。3讨论煤颗粒在土壤中的累积会使得土壤中的有机碳源由最初的植物源变为植物源+煤源。煤源碳并非真正的土壤有机碳，而且其碳质量分数极高(无烟煤90%98%，烟煤 75%90%，褐煤 60%77%)，因此即使少量的煤源碳累积也会造成土壤有机碳质量分数的明显高估6-7。然而，目前缺少可靠的方法定量土壤中的煤源碳。本研究发现：3502500nm 内，植物源有机质与不含煤粒土壤的光谱曲线相似，光谱反射率高且波峰波谷明显，其中 470nm 附近的微弱的吸收特征是受铁、锰氧化物吸收的影响，1400nm 和 1900nm附近波谷是受水分子吸收影响，2200nm 附近波谷是受 AlOH 键和 C=C 键间分子伸缩运动的影响36-38。相较于植物源有机质与不含煤粒的土壤，煤的光谱曲线特征截然不同，展现出低平的形态，表明煤源碳与植物源碳在光谱特征方面具有明显的差异性。这种差异性可以解释目前矿区土壤有机碳光谱反演研究中报道的模型估测精度不高的问题17-20。研究也发现，新鲜煤的光谱反射率低于煤粉，这是因为新鲜煤的碳质量分数(88%)高于煤粉(58%)。随着土壤中煤源碳质量分数的增加，其光谱反射率降低。这些研究结果为利用高光谱遥感技术定量反演土壤中的煤源碳提供了理论基础。本研究中，600nm 与 750nm 两处的弓曲差均与煤源碳质量分数具有极好的负相关关系，该结果区别表 2 基于光谱特征指数的煤源碳质量分数反演模型Table 2 Inversion models of coal-derived carbon content based on spectral characteristic indices光谱指数数据形式波段组合模型训练集验证集R2mRMSE/%R2mRMSE/%RPDIDOA600nmy=1328.72x2619.4x+72.670.9681.6410.9771.5326.678750nmy=605x2368.72x+56.490.9651.7010.9711.7195.954IDR406nm&454nmy=921.01x+27.010.8913.0280.9133.0183.392RC1374nm&1974nmy=298.44x+16.880.9731.4990.9602.0345.0311/R2127nm&2323nmy=79.59x+10.030.9940.6950.9930.81112.613lnR423nm&454nmy=205.1x+31.910.9541.9700.9621.9965.126R534nm&673nmy=45680.86x+25.130.8373.7130.8863.4432.973R1905nm&2142nmy=331114.03x+21.780.7064.9850.7165.4521.8771/lnR423nm&454nmy=728.91x+29.580.9242.5280.9422.4534.173R21065nm&2341nmy=29389.91x2+2719.61x+64.10.9721.5310.9731.6666.144R0.5423nm&454nmy=34639.78x2+2250.97x+37.010.9791.3120.9871.1369.011ln(1/R)423nm&454nmy=472.26x+31.910.9541.9700.9621.9965.126IRR402nm&1003nmy=12.17x+37.390.8623.4090.8743.6302.819RC1974nm&1374nmy=1510.51x2+2606.83x1092.490.9142.6920.9362.5883.9541/R402nm&1003nmy=8.43x250.64x+76.720.9781.3510.9771.5276.704lnR422nm&454nmy=347.73x+376.620.9222.5680.9362.5773.971R2318nm&2184nmy=20.77x+4.150.9512.0230.9751.6166.332R2122nm&1906nmy=2.44x+6.920.1068.7010.1319.5391.0731/lnR406nm&454nmy=422.59x3