土壤有机质(soil organic matter，SOM)是存在于土壤中动植物的残体、微生物体及其分解和合成物质，是土壤固相部分的重要组成成分，也是陆地表层重要的碳库，其含量不仅是土壤肥力的重要评判标准，也是土壤退化状态的重要指标^[1-3]。目前，有机质含量测定方法，主要有三大类。传统的测定有机质的方法主要有容量分析法、重铬酸钾比色法、重铬酸钾氧化外加热法等^[4]。这些基于化学分析的测定方法，具有操作方便、价格便宜等不可比拟的优势^[5]，虽然精度较高，但实验过程相对繁杂、费力多、周期长^[6-7]。近年来，GIS技术和地统计学相继用于土壤有机质研究中，借助空间插值比如克里格插值、反距离插值法等，实现了中大尺度范围内的推演。但空间插值方法由于采样点的数量有限且分布不均等原因，难以获得较为准确的空间预测结果^[8-9]。如张枝枝等^[10]结合GIS技术分析了渭河两岸缓冲带的土壤有机质含量的分布特征及其影响因子；柳富坤等^[11]基于GIS技术对武功山山地草甸土壤有机质的分布特征以及空间异质性进行了研究分析；吴才武等^[12]基于地统计与遥感反演相结合的方法得到了吉林省黑土区土壤有机质空间分布。20世纪70年代以来，遥感技术作为采集地球信息及其变化的一种新型手段，因其具有宏观性、动态性、经济性被广泛运用于土壤有机质估算等方面的研究中。其原理是不同类型的土壤具有不同的反射曲线，通过遥感模型，借助尺度推绎方法^[13]，实现像元尺度上的土壤有机质空间格局模拟^[14]。遥感模型估算土壤有机质的首要前提和基础是确定有机质的遥感响应波段，即土壤有机质含量的遥感信息特征波段。由于土壤是高度不均一的历史自然综合体，不同土壤类型其有机质含量本身的高低，以及同一种类土壤其成土母质、土壤水分含量、植被状况等的差异，综合地影响土壤有机质的遥感响应，导致了遥感响应波段的差异。Krishnan等^[15]发现土壤光谱在近红波段与有机质含量的相关性弱，可见光波段与有机质含量有较强的相关性；Gunsaulis等^[16]研究表明，土壤有机质含量与红光波段具有较强的相关性；刘焕军等^[17]认为620 ~ 810 nm波段范围是黑土有机质的主要响应波段，土壤有机质响应波段集中在可见光和近红外范围内；卢艳丽等^[18]认为在480 ~ 740 nm波段，黑土的有机质含量与光谱呈极显著负相关，而在816 ~ 1 415 nm波段范围内，其有机质含量与光谱呈极显著正相关；刘炜等^[19]认为530 ~ 580 nm波段范围是褐土有机质的主要响应波段；刘磊等^[20]认为红壤土有机质的响应波段在400 ~ 1 000 nm，而沈润平等^[21]则认为红壤土有机质响应波段在540 ~ 1 000 nm，即可见光与近红外区域。

我国草地面积广阔，约占全国陆地总面积的41%，青藏高原广泛分布着高寒草甸、高寒草原和温性草原等草地生态系统类型^[22-23]，其面积约为1.28×10⁸km²，其草地土壤有机碳量达到C 33.5×10⁹ t，占全国土壤有机碳量的23.44%，占全球土壤碳库的2.4%^[24]。青海湖流域是青藏高原东北部重要的牧区，流域内有国家级自然保护区，生态环境极其脆弱。在人类不合理活动和全球气候变化的双重影响下，草地退化严重，裸土、黑土滩随处可见，土壤肥力下降。因此，快速测定大尺度的土壤有机质成为青海湖流域亟需解决的科学与现实问题。

基于以上分析，本文以Landsat8-OLI多光谱遥感数据和野外实测有机质数据，以青海湖流域为实验区，利用RS、GIS技术和多元回归方法，探寻适合青藏高原有机质反演的响应波段及参数并建立模型，从而实现区域像元尺度的土壤表层有机质遥感估算，以为青藏高原草地生态系统的草地持续利用与生态修复提供一定的科学依据，为地方层面制定草地管理策略，为国家层面制定与实施生态补偿政策提供现实指导。

1 研究区概况与研究方法 1.1 研究区概况

青海湖流域地处青藏高原东北部，流域总面积为29 661 km²，地理位置36°15′ ~ 38°20′ N，97°50′ ~ 101°20′E，地形复杂，地势起伏大(图 1)。该区属典型的高寒干旱大陆性气候，以干旱、寒冷、多风为主要特征，多年平均气温-1.4 ~ 1. 7℃，平均年降水量为340 mm，蒸发量800 ~ 1 000 mm^[25]。主要植被类型为高寒草甸、高寒草原、高寒流石坡稀疏植被、沙生植被、盐生草甸、寒漠草原和沼泽草甸等。土壤主要有沼泽土、草毡土、栗钙土、寒钙土、黑毡土等。其中，草毡土、薄黑毡土和栗钙土分别占整个流域面积的31.1%、15.1%、11.6%。

1.2 研究方法 1.2.1 数据来源与处理

土样采集采用网格采样法，将研究区划分为10 km×10 km的网格。2016年9月27至10月4日，通过手持GPS定位网格上的预采样点，然后用环刀法采样。在土壤剖面每个土壤层中部平稳打入环刀，待环刀全部进入土壤后，挖去环刀周围的土壤，取出环刀，小心脱出环刀上端的环刀托，然后用削土刀削平环刀两端的土壤，用米尺量取0 ~ 10、10 ~ 20、20 ~ 30 cm土层土壤，分别取样装入样品袋中。将样品带回实验室后，称取每个采样点0 ~ 10、10 ~ 20 cm土层的土样各50 g，混匀后装入新样品袋，土壤自然风干后，然后取少量土壤样品，进行除杂，研磨，过100目筛。土壤有机质含量用重铬酸钾氧化外加热法测定。本文共得到41个采样点数据，剔除实验失败数据与异常数据，共得到29个采样点有效数据(图 1)。

1.2.2 遥感数据获取

考虑影像时间应与采样时间相近，以提高估算的准确性，本研究从地理空间数据云(http://www.gscloud.cn/)获取了Landsat8-OLI多光谱卫星遥感数据(表 1)，即2016年9月28日1景、10月8日1景和10月17日2景，共4景影像。在ENVI5.1软件中，将影像进行辐射定标、大气校正和影像拼接等预处理。用定标工具(radiometric calibration)将图像的灰度值转换成辐射亮度值，完成辐射定标；利用FLAASH进行大气校正，校正过程中主要参数设置：大气模型设置为Mid-Latitude Winter，气溶胶设置为Rural，气溶胶反演模型为2-Band(K- T)，能见度设置为默认值(40)，输入图像采集时间、研究区平均海拔为3.6 km等；用Mosaicking下的Seamless Mosaic工具进行影像拼接^[26]。为了去除微地形、影像处理与扫描过程中的噪声，利用3×3模板的均值滤波处理影像。同时，为了计算方便，将上述处理后的各波段除以1 000，然后用于相关计算与模型建立。

1.2.3 模型构建与检验

本文将采样点分为建模样本和验证样本2部分，分别为17和12个采样点数据。运用回归分析法构建有机质遥感反演模型。模型的稳定性用R²检验，R²越大，模型越稳定；模型的精度用均方根误差RMSE(式1)检验，RMSE越小，模型精度越高、预测能力越强；预测模型估计趋势用平均误差ME(式2)验证，平均误差ME < 0表明估计趋势偏高，反之ME > 0则估计趋势偏低^[28]。计算公式如下：

$ {\rm{RMSE}} = \sqrt {\frac{1}{n}\sum\limits_{i = 1}^n {{{\left( {{x_{\rm{a}}} - {x_{\rm{p}}}} \right)}^2}} } $

(1)

${\rm{ME}} = \frac{1}{n}\sum\limits_{i = 1}^n {({x_{\rm{a}}} - {x_{{\rm{p}})}}} $

(2)

式中：x_a和x_p分别为表层土壤有机质含量的实测值与预测值；n为土壤样本数。

2 结果与分析 2.1 土壤有机质遥感信息特征波段筛选

对实测土壤有机质含量和处理后Landsat8-OLI各波段影像值进行相关性分析(表 2)，结果表明，土壤有机质含量与波段1、波段2、波段3和波段4的相关性不高，与波段5、波段6和波段7的相关性较高，相关系数分别为0.736、0.752和0.584。因此，LandSat-8 OLI的5、6和7波段含有土壤有机质遥感信息，是青海湖流域土壤有机质含量遥感反演的信息特征波段。

2.2 土壤有机质遥感反演模型构建

将Landsat8-OLI的波段5、波段6和波段7这3个波段影像值进行倒数变换后，利用SPSS20.0软件中的回归分析模块构建实测土壤有机质含量与敏感波段的模型(表 3)。通过综合比较7个回归模型的预测结果及其验证指标，双波段与多波段模型的决定系数(R²=0.704)最大，表明模型较好。但双波段模型预测结果不如多波段(ME=0.51，RMSE=8.66)模拟精度高，表明多波段模型比其他双波段和单波段模型更能准确地反演土壤有机质含量，其模型为：y=205.772+ 1053.367x₃-3175.759x₂-101.319x₁(R²=0.704, P < 0.001)。据此，可得到研究区有机质含量反演的空间分布格局图(图 2)。

2.3 土壤有机质遥感反演模型评价 2.3.1 基于实测值的评价

本研究将12个验证点的有机质含量实测值与多波段模型模拟值进行比对，评价反演模型的准确性与稳定性(表 4)。通过计算得出，本研究绝对误差ME为0.51，为正数，表明模拟结果有稍低的预测趋势，由于ME较小，不超过1，对于区域尺度的土壤有机质估算，其误差完全可以接受。另外，本研究的均方根误差RMSE为8.66，与文献[29](RMSE=5.70)和文献[30](RMSE=6.93)的研究结果接近，表明本研究构建的模型能够较好地预测流域内土壤有机质含量。

2.3.2 基于文献结果的评价

为了验证土壤类别间的模拟误差，本研究结合土壤类型图，经空间统计分析，得到相应各土壤类型多波段模型模拟值的均值，然后与其他文献中各土壤类型的有机质含量平均值进行比较，得出两者间的RMSE为8.85，进一步表明模型具有很好的模拟精度(表 5)。

2.4 各类土壤有机质含量特征分析

统计分析表明，青海湖流域有机质含量63.57 ~ 85.52 g/kg土壤所占的比例最大，占整个青海湖流域总面积的36.7%，为10 897.89 km²；其次有机质含量41.63 ~ 63.57 g/kg土壤所占的面积为7 734 km²，而有机质含量0 ~ 19.69 g/kg土壤所占比例最小，占整个研究区面积的3.2%，为942.12 km²。主要土壤类型的平均有机质含量以黑钙土最高，为65.72 g/kg；暗栗钙土次之，为65.49 g/kg；石质土最低，为8.69 g/kg。

青海湖流域土壤有机质含量有随海拔上升而变大的趋势，这与李欣宇等^[35]的研究结果一致。一般气温低于0℃，有机质的分解速率小，不利于胡敏酸的形成。Jenkinson和Rayner^[36]研究表明，年平均气温下降1℃，土壤有机质含量增加2倍~ 3倍，可能是土壤温度随海拔上升而下降，气候干旱寒冷，土壤微生物活动微弱，导致土壤表层有机质含量变高。这种现象在青海湖流域分布很明显，如图 3中的区域A、B和C。土壤有机质含量与土壤含水量的关系，特别是在河谷地区，离河流越远，有机质含量呈现变大趋势。这与已有研究不符，如董智君和赵英^[37]及张帅普和邵明安^[38]研究表明，土壤水分与有机质含量呈正相关。水分越丰富，有机质分解速度越慢。本研究区出现的这种情况，可能与过度放牧有关。河谷地区一般植被覆盖度较高，是牛羊喜食地区，随着牲畜不断大量啃食、踩踏引起了局部草地退化，从而导致大量有机质流失。沙地有机质含量少，但青海湖区东面沙漠地区模拟的有机质含量值较大，不符合沙地实际有机质含量。原因可能是，沙地反射率高，导致遥感图像值大，因而模拟结果偏高。

由于研究区地形复杂，土壤类型多样，加之交通不便，本研究对部分土壤类型未采样(如沼泽土、泥炭沼泽土等)，导致模型不能较好地模拟这些土壤的有机质含量。另外，在遥感影像处理中，对云处理不够，加之地形的影响，导致同物异谱现象，对于海拔较高的土壤类型模拟结果有较大误差，如黑毡土、石灰性草甸土。未来可通过野外地物光谱仪获取实测土壤光谱，以深入分析土壤有机质含量的遥感影响因素。国内外关于土壤有机质遥感反演的研究结果差异较大，但大多数研究表明土壤有机质与波段DN值(digital number，指遥感影像像元亮度值)或反射率呈负相关，如张法升等^[39]和Chen等^[40]研究结果。但本文研究表明，青海湖流域的土壤有机质含量与多光谱遥感波段DN值呈正相关，这可能与高寒地区的区域特征有关，如高海拔、复杂多样的地形、类型多样的土壤等。

3 结论

Landsat8-OLI遥感影像反演青海湖流域土壤有机质含量的特征波段是5、6和7波段，即波长较长的红外波段，这对其他多光谱遥感平台，如MODIS、ASTER等卫星遥感影像反演青藏高原高寒草地生态系统的土壤有机质含量具有一定的借鉴意义。

本研究通过比较，确定了基于Landsat8-OLI影像3个红外波段的，精度最高、稳定性强、预测趋势平稳的土壤有机质遥感反演模型：y=205.772+ 1053.367x₃-3175.759x₂-101.319x₁。通过实测点验证与相关文献研究结果验证，该模型可以估算青海湖流域的土壤有机质含量。这为面积广大的青藏高原地区土壤有机质含量遥感估算提供了参考，即运用多光谱卫星遥感数据的红外波段估算青藏高原的土壤有机质含量具有一定的可行性。

青海湖流域内主要土壤类型平均有机质含量为：黑钙土最高，为65.72 g/kg；暗栗钙土次之，为65.49 g/kg；石质土最低，为8.69 g/kg。