土壤颜色是土壤的重要属性，与土壤组成有关。Stoner和Baumgardner^[1]对全美485个土壤研究发现，土壤颜色主要受土壤有机质和氧化铁不同程度的影响。章明奎等^[2]指出，碳酸盐岩发育的土壤其颜色差异性主要决定于土壤氧化铁矿物类型及有机质含量。张蕾等^[3]发现氧化铁随着黏粒在土壤孔隙中的迁移，淀积在土壤结构体表面上，导致某些黏化层结构体表面较其内部土壤基质的颜色更红。当赤铁矿含量较高时土壤呈现红色，而针铁矿含量较高时呈现黄色^[4]。同时，土壤颜色还与淋溶作用有关，当淋溶作用较强时，土壤颜色通常较亮^[5-6]。此外，土壤颜色可以指示土壤肥力状况。

土壤颜色的获取方式主要有Munsell色卡目视比色法和CIE三刺激值法^[7-9]。目视比色法通过人眼匹配土壤样品和标准色卡来确定土壤颜色，简单便捷成本低，但是存在主观性比较强、系统性偏差比较大的问题，难以满足土壤颜色定量化表达的需求^[10-11]。三刺激值法利用土壤可见光光谱，结合标准观察者配色函数和光源相对光谱功率分布计算土壤颜色(CIE XYZ)^[12-14]，具有精度高、抗人为因素干扰强的特点，方便进行定量化表达。虽然目前土壤颜色表达主要采用Munsell颜色空间，但是由于该系统将颜色划分为一系列不连续的色块，不适合进行数值统计分析^[15]。因此，Melville和Atkinson^[16]对比土壤颜色表达的不同方法后，建议使用CIELAB颜色系统。Rossel等^[17]建议将CIE XYZ颜色空间作为不同色彩之间转换的媒介。Fernandez和Schulze^[7]使用Cary 17D光谱分光仪获取土壤反射光谱，并采用分光积分方式计算土壤CIE XYZ颜色值，证明相比于目视匹配颜色的方式，这种积分计算土壤颜色的方式更能量化土壤颜色的细微差异。

由于土壤光谱及土壤遥感数据直接包含红、绿、蓝等不同波段的地面光谱反射信息，因此，常用于土壤颜色预测。李怡春等^[18]基于土壤高光谱反射率和偏最小二乘回归(partial least squares regression，PLSR)方法进行土壤颜色预测，并与色系转换法进行对比，证明了利用PLSR方法预测土壤Munsell颜色较色系转换法更优。Mattikalli^[19]利用Landsat MSS传感器获取土壤反射率数据并与实测土壤Munsell颜色进行建模拟合，结果具有很好的相关性。Escadafal^[20]讨论了遥感数据与土壤颜色的相关性和相关研究现状，认为利用遥感数据进行土壤颜色的建模预测是可行的。

虽然近年来使用可见–近红外光谱直接测量土壤颜色的方法研究越来越多，但是利用遥感数据进行区域土壤颜色制图还需要大量研究。如何挑选最能代表裸露地表遥感信息的反射光谱，并基于该反射光谱，进行CIE XYZ色系的土壤颜色反演是需要解决的关键问题。因此，本研究基于30 m空间分辨率Landsat-8 OLI原始反射率数据，提出遥感指数的筛选条件，并进一步采用PLSR方法进行土壤CIE XYZ颜色各分量的预测研究，为遥感获取大范围土壤颜色提供新途径。

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

研究区位于黑龙江省，该地区地势大致是西北、北部和东南部高，东北、西南部低，西部属松嫩平原，东北部为三江平原，是我国黑土分布的核心地带。黑龙江省属于寒温带与温带大陆性季风气候，年平均降水量为526.64 mm，年平均气温为3℃。该区分布着多种地貌类型，包括平原、台地、丘陵和起伏山地；土地利用类型以耕地、林地、草地为主，且耕地以旱地为主；成土母质主要为冲积–洪积物、中性冲积物、浅色结晶岩风化物、砂页岩风化物、石灰岩风化物以及黏质黄土；土壤类型包括草甸土、黑土、暗棕壤、水稻土、黑钙土、白浆土、沼泽土、风沙土和棕色针叶林土。

1.2 土壤样品采集与制备

根据地形、母质、土地利用类型、土壤类型、交通条件等在研究区布置采样点292个(图 1)，于2016年采用五点法采集土壤样品，采样时在8 m × 8 m的方形范围内4个顶点和中心位置各采1个土样，混合后装袋，取样深度为0 ~ 15 cm，每个样品1 kg。土壤样品自然风干后去除植物残体，将大土块敲碎，并取出适量样品研磨过100目(0.149 mm)尼龙筛后备用。

1.3 土壤光谱与有机碳测定

土壤样品高光谱采集采用美国ASD FieldSpec4 Hi-Res地物光谱仪，其波段范围为350 ~ 2 500 nm。土壤样品进行室内高光谱测定时，尽可能均匀地将土样填充于涂黑的样品盘中。位于样品两侧等距位置45°角的斜上方安置两个50 W的卤素灯，保证最大限度消除土壤颗粒阴影带来的测量误差。ASD探头视场角为25°，探头放置于样品正上方固定高度处，保证入射视场合理地覆盖于土壤样品盘中心。进行光谱测量前，提前30 min预热光谱仪并进行暗电流校正和白板校正，光谱测量时每隔10 min进行一次白板校正，白板校正使用的是25 cm×25 cm的漫反射标准参考板。对于每个土壤样品，为保证测量结果不受土壤表面状态的干扰，在依次旋转90°的4个方向上各测量5条光谱曲线，共20条光谱曲线，以这20条光谱曲线的平均值作为样品的测定值^[1]。通过对光谱曲线进行拼接校正，获得土壤样品在不同波段的光谱反射率值，本研究选取的波段范围位于可见光部分380 ~ 780 nm。

土壤样品有机碳含量基本信息见表 1，测定方法参照《土壤农业化学分析方法》^[21]。

1.4 土壤颜色获取

土壤颜色利用CIE XYZ颜色空间模型来表示，CIE XYZ颜色空间是由国际照明委员会在1931年提出的，其中Y代表亮度，X和Z是主光谱的虚拟组成部分。三刺激值XYZ计算公式为^[22]：

$ \left. \begin{gathered} X = K\int_{380}^{780} {S(\lambda )\bar x(\lambda )R(\lambda )d\lambda } \hfill \\ Y = K\int_{380}^{780} {S(\lambda )\bar y(\lambda )R(\lambda )d\lambda } \hfill \\ Z = K\int_{380}^{780} {S(\lambda )\bar z(\lambda )R(\lambda )d\lambda } \hfill \\ K = \frac{{100}}{{\int_{380}^{780} {S(\lambda )\bar y(\lambda )d\lambda } }} \hfill \\ \end{gathered} \right\} $

(1)

式中：$ S(\lambda ) $为照明光源的相对光谱功率分布；$ \bar x(\lambda ) $、$ \bar y(\lambda ) $、$ \bar z(\lambda ) $为CIE 2°标准观察者配色函数；$ R(\lambda ) $为物体的光谱反射比。其中，$ S(\lambda ) $、$ \bar x(\lambda ) $、$ \bar y(\lambda ) $、$ \bar z(\lambda ) $都是已知的，本文涉及的计算中$ S(\lambda ) $为C光源的相对光谱功率分布。前人认为，对大多数实用的计算，波长范围取380 ~ 780 nm，波长间隔取$ \Delta \lambda $=5 nm，不会产生明显的计算误差^[23]，实际使用中，由于光谱数据为离散型数据，通常以求和的方法替换上式中的积分方式来简便计算^[24]。因此，本文采用5 nm波长间隔，以求和方式完成土壤颜色的积分计算。

1.5 遥感数据的获取与处理

Landsat-8遥感影像数据集来自于Google Earth Engine(GEE)，采用经过大气校正和正射校正的Landsat-8表面反射率数据产品(Surface Reflectance Tier 1)^[25]。本研究采用Landsat-8 OLI遥感影像数据集中的可见–近红外波段以及短波红外波段(共7个波段)进行土壤颜色预测。考虑到遥感影像由于云、雪等对下垫面遮挡造成的遥感数据缺失以及土壤颜色在短期内变化幅度很小两个因素，对使用遥感光谱数据的时间进行适当放宽，在土壤样品采集年份(2016年)前后各延伸1年，最终选取了黑龙江省2015—2017年遥感影像，并通过质量评估波段去除包含云、云阴影、雪等干扰的影像，获得了所有土壤样点在1 248景遥感影像上的共9 885条光谱数据。

本研究首先提出筛选建模光谱的遥感指数条件，然后通过遥感指数构建建模光谱数据的筛选条件。选用的遥感指数分别为归一化差值植被指数(normalized difference vegetation index, NDVI)、归一化差值水体指数(normalized difference water index, NDWI)、归一化差值湿度指数(normalized difference moisture index, NDMI)、归一化差值不透水面指数(normalized difference impervious surface index, NDISI)^[26-27]。各指数具体计算公式为：

$ \left. {\begin{array}{*{20}{c}} {NDVI = ({B_5} - {B_4})/({B_5} + {B_4})} \\ {NDWI = ({B_3} - {B_5})/({B_3} + {B_5})} \\ {NDMI = ({B_5} - {B_6})/({B_5} + {B_6})} \\ {NDISI = \frac{{({B_{10}} - ({B_3} + {B_5} + {B_6})/3)}}{{({B_{10}} + ({B_3} + {B_5} + {B_6})/3)}}} \end{array}} \right\} $

(2)

式中：B₃、B₄、B₅、B₆、B₁₀分别与Landsat-8 OLI遥感影像数据集波段相对应，代表绿色、红色、近红外、短波红外、热红外波段。通过对遥感影像光谱反射率数据和土壤样品实测反射率数据的对比分析，确定了筛选光谱约束条件为：

$ \left. {\begin{array}{*{20}{c}} {0.1 \leqslant NDVI \leqslant 0.4} \\ {NDWI \leqslant - 0.3} \\ {NDMI \leqslant 0} \\ {NDISI \leqslant 0.4} \end{array}} \right\} $

(3)

对于每个土壤样点，满足所有筛选条件的遥感影像光谱通常不只一条，还需进一步选取。先计算各波段光谱平均值，然后计算筛选出来的各条光谱与该平均光谱的欧氏距离，并选择欧氏距离最小的一条光谱曲线作为对应点用来进行颜色预测的光谱数据。欧氏距离计算公式为：

$ D = \sqrt {\sum\nolimits_{i = 1}^n {{{({B_i} - {{\bar B}_i})}^2}} } $

(4)

式中：D为两条光谱曲线之间的欧氏距离；i、n分别为波段序号和波段总数；$ \rm{B_i} $为波段i处光谱反射率；$ {\bar B_i} $为波段i处平均光谱反射率。

1.6 建模与评价方法

本研究采用PLSR模型对土壤颜色进行建模预测。PLSR算法是一种多元线性回归建模方法，在建模过程中集成了主成分分析、典型相关分析和线性回归分析方法的特点，能够有效处理变量间存在多重相关性而样本量较少的问题^[28]。主成分数选取对PLSR建模尤为重要，主成分数过少，没有足够信息进入模型拟合，会导致拟合不充分，降低模型的预测精度；主成分数过多，很多无用信息会进入模型，容易出现过拟合，同样会降低模型预测精度。本文利用赤池信息量准则值(akaike information criterion, AIC)来选择模型的最佳主成分数：

$ \text{AIC}=n \text{ lnRMSE}+2p $

(5)

式中：n是建模集样本数量；RMSE为验证均方根误差；p为对应的主成分数，以AIC最小时的主成分数作为PLSR模型的最优主成分数^[29]。评价建模精度时使用验证决定系数R²、RMSE和RPD(测定值标准偏差与标准预测误差的比值)：

$ \rm RPD = \frac{{SD}}{{RMSE}} $

(6)

式中：SD表示数据集实测值标准差。已有研究结果表明，当RPD≥2.50时，模型可以精确预测；当2.50 > RPD≥2.00时，模型可以很好地预测；当2.00 > RPD≥1.80时，模型可以较好地预测；当1.80 > RPD≥1.40时，模型可以近似地预测；当1.40 > RPD≥1.00时，模型可以区分变量值的大小；当RPD < 1.00时，模型不适合预测^[29]。

建模过程中数据集被随机划分成70% 的建模集(204个土样)和30% 的验证集(88个土样)，并选择Landsat-8遥感产品数据集可见–近红外、短波红外共7个波段作为建模变量，因变量分别为CIE XYZ三种颜色分量。根据AIC标准，PLSR方法预测CIE XYZ颜色各分量的最佳主成分数参数设定为5、5、6，其余参数为默认值。

文中涉及的光谱值处理、色系值等计算在Python 3.6中完成，PLSR建模采用基于Python语言编制的scikit-learn机器学习库^[30]完成，图表由Excel 2013、ArcGIS 10.2制作完成。

2 结果 2.1 土壤颜色实测结果统计特征

实测土壤样品CIE XYZ颜色结果如图 2所示，各颜色分量呈近似正态分布。其中，X值可以指示人眼对红色色光的感应强弱，其均值为9.11，标准差为2.52，最大值为14.63，最小值为4.57，总体跨度为10.06；Y值可以指示人眼对绿色色光的感应强弱，同时反映着土壤颜色亮度大小，其均值为9.01，标准差为2.47，最大值为14.48，最小值为4.59，总体跨度为9.89，表明土壤颜色总体亮度偏低；Z值可以指示人眼对蓝色色光的感应强弱，其均值为7.33，标准差为1.73，最大值为11.57，最小值为4.33，总体跨度为7.24。总体来看，CIE XYZ颜色各分量跨度不大，均值偏低，这是因为土壤样品本身光谱反射率较低，加之样品中有机碳含量丰富，导致土壤颜色偏暗。

2.2 土壤颜色预测结果 2.2.1 光谱筛选对土壤颜色预测的影响

光谱筛选前后土壤颜色预测结果对比如图 3所示。基于筛选后的建模光谱进行土壤颜色的建模预测，结果如图 3A所示。在X颜色分量上，PLSR预测结果R²、RMSE、RPD分别为0.76、0.94和1.74，说明该模型可以近似地预测X颜色分量；在Y颜色分量上，PLSR预测结果R²、RMSE、RPD分别为0.76、0.91和1.76，说明该模型可以近似地预测Y颜色分量；在Z颜色分量上，PLSR预测结果R²、RMSE、RPD分别为0.69、0.66和1.68，说明该模型也可以近似地预测Z颜色分量。总体来看，基于PLSR方法进行建模预测3种土壤颜色分量是可行的，均可实现近似预测。为了比较光谱筛选前后土壤颜色预测精度，本文对光谱筛选前遥感数据进行建模预测。图 3B显示了进行光谱筛选前的土壤颜色的建模预测结果，可见土壤颜色CIE XYZ各分量PLSR预测结果RPD分别为1.26、1.27和1.26，均低于光谱筛选后的建模预测结果，说明进行本文光谱筛选能够有效提升土壤颜色的建模预测精度。

2.2.2 土地利用对土壤颜色预测的影响

为了探究土地利用类型对土壤颜色预测精度的影响，按照不同土地利用类型分别对预测结果进行拟合分析，结果如图 4所示。在土壤颜色CIE XYZ各颜色分量上，耕地类型的预测结果均优于草地和林地类型的预测结果。其中，耕地类型土壤颜色预测RPD分别为2.42、2.45和2.21，这表明模型能够对耕地土壤颜色进行很好的预测；林地类型土壤颜色预测RPD分别为1.69、1.69和1.42，这表明模型能够近似预测林地土壤颜色；草地类型土壤颜色预测RPD分别为0.98、1.01和1.28，表明模型对于草地土壤颜色的预测效果较差。可见，在不同土地利用类型下土壤颜色的预测结果具有明显的差异性。这是由于耕地土壤反射光谱较为一致，而林地、草地由于地表覆盖、腐殖质积累等差异较大，导致反射光谱中混杂有较多的噪声信息。此外，本研究中林地和草地的样本数量相对偏少，这也会导致土壤颜色的预测结果出现这种差异性。

2.2.3 有机碳含量对土壤颜色预测的影响

为了探究不同有机碳含量对土壤颜色建模预测结果的影响，进一步按照不同有机碳含量将验证集土壤样本划分成高值区、中值区和低值区3个子集，并分别对3个子集进行建模预测，结果如图 5所示。其中，高值区子集有机碳含量范围为28.25 ~ 88.92 g/kg，共30个土壤样品；中值区子集有机碳含量范围为21.00 ~ 28.19 g/kg，共29个土壤样品；低值区子集有机碳含量范围为7.42 ~ 20.82 g/kg，共29个土壤样品。分别对有机碳含量高、中和低值区3个子集的土壤颜色建模预测结果进行验证，结果显示，在CIE X颜色分量上，预测土壤颜色验证RPD分别为1.39、1.50和1.77；在CIE Y颜色分量上，预测土壤颜色验证RPD分别为1.35、1.55和1.82；在CIE Z颜色分量上，预测土壤颜色验证RPD分别为1.02、1.97和2.07。总体来说，当土壤有机碳含量较低时，各颜色分量验证RPD较高，预测效果较好。

3 讨论

本文基于土壤反射光谱，依据三刺激值法计算土壤CIE XYZ色彩空间各颜色分量值，并结合遥感数据进行土壤颜色建模预测。首先定义了光谱数据筛选方法对遥感光谱数据进行筛选，然后采用PLSR方法对土壤颜色各颜色分量进行建模预测，结果显示，CIE XYZ颜色各分量的验证RPD分别为1.74、1.76和1.68，表明本方法建立的模型可以对土壤颜色进行可近似预测。土壤颜色的预测受多种因素影响，本文从光谱筛选、土地利用类型和有机碳含量3个方面进行讨论。

从光谱筛选过程来看，本研究通过确定遥感指数NDVI、NDWI、NDMI、NDISI的阈值，进行建模光谱数据的筛选。其中，NDVI阈值主要消除植被信息，NDWI、NDMI阈值主要消除水分信息，NDISI阈值主要消除不透水地面信息。对比光谱筛选前后的预测结果，可以看见土壤颜色预测结果RPD由光谱筛选前的1.26、1.27和1.26分别提升至1.74、1.76和1.68，这说明本文提出的光谱筛选方法能够有效提升土壤颜色的预测精度，通过筛选获得的光谱，可以在最大程度上排除植被等外部参数的干扰。实际上，卫星过境与地面纯裸露时间不一定正好匹配，即便是裸露地表也或多或少存在诸如腐烂秸秆、稀疏植被、动植物残体等的覆盖，不可能获得像在实验室条件下制作非常干净的裸土光谱，因此，进一步消除植被等外部影响的干扰，提升遥感数据预测土壤颜色的精度是有可能的。

从土地利用类型影响来看，对比不同土地利用类型下土壤颜色的预测结果，发现存在差异，耕地预测效果要优于林地、草地。耕地类型土壤颜色CIE XYZ各分量值预测RPD分别为2.42、2.45和2.21，这表明模型能够对耕地土壤颜色进行很好的预测。这是因为耕地土壤反射光谱较为一致，而林地、草地由于地表覆盖、腐殖质积累等差异较大，反射光谱中混杂有较多的噪声信息。此外，由于本研究中林地和草地的样本数量相对偏少，在土壤颜色建模预测的过程中偏差相对较大，这也会导致土壤颜色的预测结果较差。

从有机碳含量影响来看，Stoner和Baumgar­dner^[1]研究表明，土壤颜色主要受土壤有机碳和氧化铁组分影响。本研究采集的土壤样品分布于黑龙江省全境，土壤类型以草甸土、黑土、暗棕壤和水稻土居多，有机碳含量丰富。土壤样品颜色主要受有机碳含量影响。本研究中，不同有机碳含量土壤颜色预测结果表明，有机碳含量较低时，土壤颜色预测结果较好。出现这种结果可能与土壤本身反射光谱的能力有关，当有机碳含量较低时，土壤颜色较亮，反射光谱能力较强，反射光谱中包含更多的有效信息，预测结果较好。

本研究中土壤样品全部来自于黑龙江省，具有一定的区域性。土样采集涵盖当地典型土地利用类型和土壤类型，土壤颜色主要受有机碳影响。对于其他受不同土壤颜色主导因素影响的土壤如红壤、紫色土等是否具有相同结果还需进一步验证。

4 结论

本研究基于土壤实测反射光谱和Landsat-8遥感数据进行土壤颜色建模预测，并提出了建模光谱数据筛选方法。采用PLSR建模预测CIE XYZ颜色各分量，验证结果RPD分别为1.74、1.76和1.68，说明按本方法建立的模型可以对X、Y、Z各分量进行近似预测，且耕地土壤颜色预测效果要优于林地和草地。分别对有机碳含量高、中和低值区土壤颜色进行预测，各颜色分量的验证RPD随着有机碳含量的降低而升高，说明当有机碳含量较低时，土壤颜色预测结果较好。