本研究使用两个数据集: EcoSIS实验数据集和PROSPECT模拟数据集。

1.1.1 EcoSIS数据集

EcoSIS数据集是由Kothari等^[12]、 Wang^[13]等于2017年— 2019年间在美国、 加拿大等各野外地区采用PSR-3500型全量程光谱仪测量得到。 其中包含近2 000个植被叶片样本, 灌木、 乔木等100种植物, 包括叶片的反射光谱以及各种生物化学参数, 如叶绿素a/b含量、 类胡萝卜素含量等, 满足光谱的多样性需求。 数据集可以在线获得(https://ecosts.org), 所采用标准化协议参考Asner^[14]的研究文献。

1.1.2 PROSPECT模拟数据集

采用改进后的PROSPECT-D模型生成模拟数据集, 该数据集用于与EcoSIS数据集融合构成机器学习建模所需的混合训练集, 达到增强机器学习数据训练的效果。 其中, 叶绿素和类胡萝卜素模拟参数的含量范围由EcoSIS数据集限制, 其他参数参考以往文献设定(见表1)。

表1

Table1

表1(Table1)

表1 PROSPECT-D 模型参数设置 Table1 Parameter configurations of PROSPECT-D model 参数 取值范围 变量数 叶绿素含量cab/(μ g· cm-2) 1.0~40.0 类胡萝卜素含量car/(μ g· cm-2) 0.1~6.0 花青素含量canth/(μ g· cm-2) 0.0~12.0 等效水分厚度cw/(g· cm-2) 0.004 5~0.041 5 4 998 干物质含量cm/(g· cm-2) 0.001~0.025 结构参数N 1.0~3.0 棕色素cbp 0.0

表1 PROSPECT-D 模型参数设置 Table1 Parameter configurations of PROSPECT-D model

1.1.3 混合数据集的构建

为了有效构建混合训练数据集, 使用K-S检验(Kolmogorov-Smirnov, K-S)算法, 按3∶ 1的比例将EcoSIS数据集划分为完全独立的A、 B两部分, A部分参与构建建模集, B部分作为模型的验证集。 每次随机挑选PROSPECT模拟数据, 与EcoSIS数据集A部分相混合, 按1∶ 3的混合比例构成机器学习的训练集。 从表2可见, A、 B两部分的最值、 平均值和标准差皆相近, 可以较为均匀地覆盖训练集样本和测试集样本。

表2

Table 2

表2(Table 2)

表2 样本集色素计量数据统计(μ g· cm-2) Table 2 Statistics of pigmentation in the sample set (μ g· cm-2) EcoSIS数据集 最小值 最大值 平均值 标准差 A部分 cab 2.05 38.91 10.54 6.18 car 0.19 5.30 2.05 0.43 B部分 cab 1.67 39.16 9.69 4.95 car 0.20 5.12 1.55 0.72

表2 样本集色素计量数据统计(μ g· cm-2) Table 2 Statistics of pigmentation in the sample set (μ g· cm-2)

根据作用的不同, 光谱预处理可以分为基线校正、 散射校正、 平滑处理和尺度缩放四类。 其中, 基线校正包括一阶导数(first-order derivative, FD)、 二阶导数(second-order derivative, SD)、 连续小波变换(continuous wavelet transform, CWT)等, 主要应用于消除背景基线漂移对信号的影响。

混合训练集中, DHRF光谱与BRF光谱反射率之间的偏差f会影响最终的反演精度。 偏差f随光谱区域(如可见光和近红外区域)和植物种类而变化, 假设该f值线性可加^{[11, 15]}, 同时可以把f看作是一种“ 系统背景” 或者“ 系统漂移信号” , 利用导数和小波变换来消除部分f值对BRF的影响。 即, 用小波系数谱和一阶导数谱代替原反射光谱, 作为混合模型反演时的输入值, 以降低表面反射效应对模型反演效果的干扰。

连续小波变换使用缩放和移位的母小波对反射光谱进行卷积。 高斯函数的二阶导数与植物叶片吸收特征光谱的形状相似^{[10, 16]}, 而第四尺度的缩放又在叶片色素反演上有更好的表现, 因此选择高斯函数的二阶导数为母小波函数, 并在第四尺度缩放。 CWT 调用Python3.7版本的pywave库对光谱进行处理; 一阶导数直接在Python软件实现。

BRF与DHRF的关系如式(1)和式(2)描述。

$\text{BR}{{\text{F}}_{\left( {{\theta }_{S}}, \ \ {{\varphi }_{S}}\ ; \ \ {{\theta }_{V}}, \ \ {{\varphi }_{V}}\ ; \ \ \lambda \right)}}=\pi \times \text{BRD}{{\text{F}}_{\left( {{\theta }_{S}}, \ \ {{\varphi }_{S}}\ ; \ \ {{\theta }_{V}}, \ \ {{\varphi }_{V}}\ ; \ \ \lambda \right)}}=\pi \times \frac{\ \ \ \ {{R}_{\left( {{\theta }_{S}}\ , \ {{\varphi }_{S}}\ ; \ {{\theta }_{V}}\ , \ {{\varphi }_{V}}\ ; \ \lambda \right)}}\ \ \ \ \ }{{{I}_{\left( {{\theta }_{S}}\ , \ {{\varphi }_{S}}\ ; \ \lambda \right)}}\ \ \ }\ \ \ \ \ $(1)

$\text{DHR}{{\text{F}}_{\left( {{\theta }_{S}}, {{\varphi }_{S}}; \ \lambda \right)}}=\int\limits_{0}^{\pi }{\underset{0}{\overset{\frac{ \pi }{2}}{\mathop \int }}\, }\text{BRD}{{\text{F}}_{\left(\ {{\theta }_{S}}, \ {{\varphi }_{S}}\ ; \ {{\theta }_{V}}, {{\varphi }_{V}}\ ; \ \lambda \ \right)}}\ \ \text{cos}{{\theta }_{V}}\ \text{sin}{{\theta }_{V}}\ \text{d}{{\theta }_{V}}\ \text{d}{{\varphi }_{V}}\ \ \ $(2)

Jay^[9]将二者关系简化为式(3)

BRFλ=DHRFλ+fλ(3)

其中

fλ=BRFsurf(λ)-DHRFsurf(λ)(4)

反射率由表面反射分量(surface)和漫反射分量(diffuse)构成, 前者描述了叶子的表面反射, 后者描述了叶片内部的吸收和散射过程。 式(4)反映了BRF_surf和DHRF_surf两个表面反射分量之间的差异, 也是DHRF光谱和BRF光谱之间的差异。

考虑到大部分波长与光谱特征相关, 且相邻波段间存在共线性问题, 这些会影响叶片的色素含量分析^[2]。 为此, 本文采用竞争性自适应重加权算法(competitive adaptive reweighted sampling, CARS)、 连续投影算法(successive projections algorithm, SPA)、 主成分分析(principle component analysis, PCA)三种算法^{[17, 18]}, 对原始光谱和分别经过一阶导处理、 小波变换处理的光谱, 进一步筛选出与叶绿素、 类胡萝卜素最具有相关性的特征变量, 再导入机器学习模型中。

已有研究结果表明在光谱子域检索叶片生化性状, 可能获得更好的反演结果^{[19, 20]}, DHRF_suff和BRF_suff之间的偏差f在不同的光谱区域有所差异。 通过对不同特征变量提取方法组合的结果分析(见表3), 综合考虑全光谱域(400~2 400 nm)和可见-近红外光谱子域(visible and near infra-red, VNIR)反演结果可能带来的不同影响, 依据前述FD、 CWT 2种预处理方法, 以及CARS、 SPA、 PCA 3种算法降维后的特征光谱变量, 设计12种组合, 分别对叶绿素和类胡萝卜素含量开展反演。 单独使用原始光谱和经一阶导、 小波变换处理的光谱所构模型作为对照组。 各组合的变量筛选详见后文。

表3

Table 3

表3(Table 3)

表3 不同特征变量提取方法组合 Table 3 Combination of different feature variable extraction methods 提取/处理方法 变量数 全波段 VNIR cab car cab car 原始光谱 2 001 601 FD 2 000 600 CWT 2 001 601 FD+CARS 52 68 71 74 FD+SPA 5 9 28 11 FD+PCA 7 CWT+CARS 34 59 33 26 CWT+SPA 18 15 12 15 CWT+PCA 7

表3 不同特征变量提取方法组合 Table 3 Combination of different feature variable extraction methods

PROSPECT模型约束下的机器学习反演植物叶片胡萝卜素、 叶绿素流程见图1, 其中由人工神经网络模型(artificial neural network, ANN)在混合数据集中分别构建叶绿素和类胡萝卜素含量监测模型, 对比经过单一算法和组合算法处理过的模型性能。 采用测试集决定系数(coefficient of determination of prediction, RP2)和测试集均方根误差(root mean square error of prediction, RMSEP)进行模型效果评价。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 PROSPECT模型约束下的机器学习反演植物叶片胡萝卜素、 叶绿素流程图Fig.1 Flow chart of machine learning inversion of plant leaf carotene and chlorophyll contents with PROSPECT model constraints

在EcoSIS数据集和仿真数据集中各自选取第一个样本, 其原始光谱如图2(a)所示; 分别进行一阶导和小波变换处理后的图像如图2(b), 图2(c, d)所示。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 FD与CWT处理前后的光谱图像 (a): PROSPECT模拟光谱与EcoSIS实测光谱反射率曲线; (b): 模拟与实测光谱一阶导系数谱; (c): 模拟与实测光谱400~1 000 nm小波系数谱; (d): 模拟与实测光谱1 000~2 400 nm小波系数谱Fig.2 Spectral images before and after FD and CWT preprocessings (a): PROSPECT simulated and EcoSIS measured reflectivity spectra; (b): FD spectra; (c): CWT coefficient spectra (400~1 000 nm); (d): CWT coefficient spectra (1 000~2 400 nm)

图2(a)可见, 模拟光谱DHRF与实测光谱BRF之间存在反射率的差异, 难以在振幅上将二者进行直接比较。 并且在反射率空间中, 只能观察到叶绿素在450和680 nm处和水分在1 450和1 950 nm 处较强的吸收特征。 由于叶片色素在可见蓝光、 红光区域的强吸收, 而在绿波段强反射, 因此会在420和680 nm附近出现两个吸收谷, 在约600 nm处出现一个反射峰。 对比图2(a)和(b)可见, 经一阶导处理后的光谱更显著地突出了原始光谱的两“ 谷” 一“ 峰” , 从中可见: 以拐点形式所捕捉到的叶片黄边(560~580 nm)陡然降低、 红边(680~750 nm)陡然增加, 光谱导数在增强局部光谱特征和消除微小的光谱变化方面具有一定优势, 小波变换同样可以较好地捕捉到这些差异。 如图2(c, d)所示, 第四尺度的小波系数谱展示了更为具体的窄带吸收特征。 以对植被光谱特征具有重要指示意义的两个波段, 即可见光黄波段(560~580 nm)和红波段(680~750 mm)为例, 相较于原始光谱, 小波系数谱更易于观察到关于叶片光谱黄边、 红边特征的变化细节, 也意味着可能包含了更丰富的叶片叶绿素、 类胡萝卜素的光谱信息。 这两种处理方法或能有助于减少混合训练集中实测BRF光谱和模拟DHRF光谱之间的差异, 并增强叶片色素的吸收特征, 改善植物叶片叶绿素、 类胡萝卜素的提取效果。

2.2.1 基于CARS的特征光谱变量提取

对样本的一阶导数谱和小波变换系数谱采用 CARS算法提取特征光谱变量, 设置蒙特卡罗(Monte-Carlo, MC)采样50次, 采用5折交叉验证法, 得到一阶导数谱在全光谱域筛选叶绿素含量特征波长的运行结果(见图3)。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 一阶导数谱筛选叶绿素特征波长CARS运行结果Fig.3 Screening chlorophyll characteristic wavelengths of FD spectra by CARS

图3(a)、 (b)和(c)分别表示在CARS算法运行过程中, 特征变量数、 五折交叉验证均方根误差(root mean square error of cross validation, RMSECV)和各特征变量回归系数三者随着采样次数的增加而发生的变化。 从图3(a)可以看出, 变量数经历了先“ 粗选” 后“ 精选” 2个选择过程, 在前5次采样中快速下降, 随后缓慢减少并趋于平稳。 图3(b)展示了RMSECV先逐渐减小后陡然增加的趋势, 在其达到最小值时, 采样运行14次, 各特征光谱变量的回归系数由图3(c)中的“ * ” 号垂线标示。 此后若采样次数继续增加, 光谱中某些包含叶绿素关键信息的重要变量将被剔除, RMSECV也随之增加。 因此, 选择第14次采样提取的52个变量为叶绿素含量的特征光谱变量, 现有维度相比原始维度下降了97.4%。

如图3(d)所示, 特征波段大多数分布于400~800 nm之间, 包含了可见蓝光、 绿光区域, 也是对植物光合作用活性最为敏感的区域。 叶片的“ 黄边” 特征是指在560~580 nm间的窄带, 由于其叶绿素、 类胡萝卜素从绿波段强反射过渡到黄波段强吸收, 光谱反射率曲线陡然降低的特点。 “ 红边” 特征则是指叶片光谱反射率在680~750 nm之间陡然攀升的特点, 也表征了叶片光学特性的主要决定因素从叶绿素对红波段的强吸收, 转变为叶片内部细胞结构对近红外波段的强反射。 CARS算法在黄边位置(560~580 nm区间反射率降低最快的点), 即一阶导数谱的第一个拐点附近筛选出了若干特征波长, 红边位置(680~750 nm区间反射率上升最快的点)即第二个拐点右侧也有少量波长被采集。 使用该算法提取的特征变量基本具备叶绿素的关键光谱信息。

2.2.2 基于SPA的特征光谱变量提取

对样本的一阶导数和小波变换谱采用SPA算法提取特征光谱变量, 计算在不同个数子集的特征光谱变量下的均方根误差RMSE, 选取该值最小者为最优特征变量子集。 以小波系数数谱在VNIR子域的类胡萝卜素敏感波长提取为例, 设置5~50为特征变量数的变化范围, 当提取特征变量个数为15时, RMSE取得最小值3.927。 该15个特征光谱变量占原始波段数的2.5%, 在VNIR光谱子域的分布如图4所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 基于SPA的特征变量提取 (a): RMSE迭代变化图; (b): SPA特征变量分布图Fig.4 Extracting characteristic variables by SPA (a): RMSE iterative change chart; (b): Distribution of characteristic spectral variables

2.2.3 基于PCA的主成分提取

对样本的原始高光谱、 一阶导数和小波变换光谱数据进行主成分分析, 前10个主成分(principal components, PC) 的方差贡献率如图5所示。 以原始光谱的PCA为例, 第1个PC贡献率为72.78%, 是所有PC中贡献率最大的; 前2个PC累积贡献率为86.68%, 之后各PC累积贡献率缓慢增加; 前7个PC的累积贡献率为99.89%, 仅剩0.11%的光谱信息未能表达。 以方差累积贡献率≥ 95%以上为阈值, 原始光谱、 一阶导光谱和小波系数谱均选择最佳PC个数为7, 并以该7个主成分作为特征变量进行建模, 为原始高光谱维度的0.3%。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 基于PCA的原始光谱、 一阶导数谱、 小波系数谱主成分提取Fig.5 Extractions of principal components of the original, first-order derivative, and wavelet coefficient spectra

2.3.1 ANN模型的叶绿素含量预测结果

在人工神经网络ANN模型中, 选取“ Relu” 函数作为神经网络的激活函数, “ Adam” 为权重优化器, 并设置ANN的隐含层神经元个数为50(4层), 最大迭代次数2 000。 用原始光谱及6种光谱处理方法提取的特征变量建立对应的ANN模型, 得到叶绿素含量预测结果见表4。 为避免自变量数值差异过大影响建模精度, 建模前将所有自变量作归一化处理。

表4

Table 4

表4(Table 4)

表4 不同处理方法建立的ANN模型叶绿素含量预测结果 Table 4 Prediction results of chlorophyll content of ANN models established by different treatments 处理方法 全光谱域 VNIR子域 变量数 RP2 RMSEP 变量数 RP2 RMSEP 原始光谱 2 001 0.661 4 3.398 6 601 0.711 4 2.481 7 FD 2 000 0.819 9 3.810 4 600 0.779 2 4.264 2 CWT 2 001 0.805 2 3.103 5 601 0.795 6 2.480 7 FD+CARS 52 0.806 4 2.911 4 71 0.754 6 5.843 1 FD+PCA 7 0.711 0 3.902 6 7 0.645 3 3.902 6 FD+SPA 5 0.304 2 7.646 4 28 0.675 5 4.351 6 CWT+CARS 34 0.781 1 3.888 9 33 0.540 4 5.320 9 CWT+PCA 7 0.671 5 4.162 7 7 0.712 6 4.026 9 CWT+SPA 18 0.434 1 4.396 9 12 0.727 8 4.897 6

表4 不同处理方法建立的ANN模型叶绿素含量预测结果 Table 4 Prediction results of chlorophyll content of ANN models established by different treatments

表4结果表明: 对于叶绿素反演, 无论全光谱域还是VNIR子域, 经过一阶导和小波变换处理均表现出了相似的性能, 相对于原始光谱有更高的准确率。 在全光谱域下, 测试集决定系数 RP2分别达到0.819 9、 0.805 2, 表明经过该二者处理的光谱较好地削弱了“ 系统背景” , 在利用PROSPECT模拟数据增强机器学习训练集数据质量的同时, 也减小了原始混合光谱之间较大的偏差f对反演精度的影响。 此外, 还增强了叶片色素的光谱吸收特征, 叶绿素预测模型的反演效果得以提升。 与小波变换处理相比, 一阶导处理更好地提升了 RP2值, 不过RMSEP较高, 其RMSEP值比原始光谱增大12.1%, 比CWT高22.8%, 这也许意味着反演精度受各参数之间的相关性影响较大。 FD+CARS相比于FD得到了一定的改善, RMSEP降低了23.6%, RP2, RMSEP也优于CWT+CARS、 CWT+PCA组合, 而FD+SPA、 CWT+SPA的预测效果均不佳。

2.3.2 ANN模型的类胡萝卜素含量预测结果

不同处理方法ANN模型的类胡萝卜素含量预测结果如表5。

表5

Table 5

表5(Table 5)

表5 不同处理方法建立的ANN模型类胡萝卜素含量预测结果 Table 5 Prediction results of carotenoid content of ANN models established by different treatments 处理方法 全光谱域 VNIR子域 变量数 RP2 RMSEP 变量数 RP2 RMSEP 原始光谱 2 001 0.594 3 0.404 9 601 0.625 7 0.760 1 FD 2 000 0.762 3 0.355 4 600 0.720 6 0.526 2 CWT 2 001 0.796 2 0.471 7 601 0.761 3 0.774 6 FD+CARS 68 0.613 6 0.881 8 74 0.734 5 0.411 9 FD+PCA 7 0.616 3 0.568 1 7 0.661 0 0.565 3 FD+SPA 9 0.558 9 0.341 4 11 0.677 1 0.586 5 CWT+CARS 59 0.693 7 0.442 1 26 0.797 2 0.414 1 CWT+PCA 7 0.635 2 0.737 2 7 0.643 1 0.464 1 CWT+SPA 15 0.415 0 0.396 8 15 0.698 3 0.723 8

表5 不同处理方法建立的ANN模型类胡萝卜素含量预测结果 Table 5 Prediction results of carotenoid content of ANN models established by different treatments

从表4、 表5的分析可以发现, 大部分光谱处理组合的模型对两种色素的预测表现都高于原始光谱, 叶绿素的模型预测表现整体高于类胡萝卜素。 总体上看, 叶片的黄边特征、 红边特征不受水分和类胡萝卜素的影响, 但其红边位置, 即680~750 nm区间反射率攀升最快点会随叶绿素含量的增加而红移, 并且在400~800 nm的可见光区域, 叶片的光谱特征由叶绿素的吸收主导, 因此叶绿素的模型预测精度较高。 而类胡萝卜素的反演因其吸收特征会被可见光区域中的叶绿素和短波红外(short wavelength infra-red, SWIR)区域中的水分含量所掩盖^[16], 所以类胡萝卜素反演精度低于叶绿素的反演精度。 两种色素的最佳预测结果见散点图6。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 叶绿素和类胡萝卜素最佳预测结果 (a): FD+CARS+全光谱域+叶绿素; (b): CWT+CARS+VNIR+类胡萝卜素Fig.6 The best prediction results of chlorophyll and carotenoid (a): FD+CARS+Full spectrum domain+chlorophyll; (b): CWT+CARS+VNIR+carotenoid

2.3.3 特征提取对反演结果的影响

不同组合的特征波长提取结果见图7。 经过降维处理的模型并没有更优异的预测表现。 例如表5的FD+CARS、 CWT+CARS组合, 有效剔除了全光谱域中与叶绿素、 类胡萝卜素相关性小的波段, 提高了信噪比, 在保持准确率基本不变的基础上减少了建模所需要的波长, 但模型很难有比全谱更突出的预测性能。 FD+SPA在VNIR子域所提取波长能覆盖叶绿素、 类胡萝卜素的部分敏感波段, 却使得蓝光波段出现“ 波长聚集” 现象; CWT+SPA在全光谱域提取的特征波长分布较均匀, 然而弱相关波长和噪声较多, 预测效果差。 使用PCA能够较好地解决多重共线性的问题, 各模型性能相对稳定, 但准确率不算最优。 从结果来看, 光谱特征筛选的主要作用是降维并减少模型负载, 降低模型训练所耗时间。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 不同组合所提取特征波长位置 (a): VNIR子域下的特征波长位置; (b): 全光谱域下的特征波长位置Fig.7 Characteristic wavelengths extracted with different combinations (a): Within the VNIR subdomain; (b): Within the full spectral domain

2.3.4 不同光谱域对反演结果的影响

在VNIR子域, 原始光谱和CWT处理过的光谱在反演叶绿素方面有所改善, 拟合度和稳定性得到提高, 且较好地提高了SPA算法组合的 RP2。 另外, 尽管结合特征提取算法的预测模型性能大部分逊于对照组FD、 CWT(表4、 表5), 但在VNIR子域, 算法组合模型的预测表现相比全光谱域有一定幅度的提升, 例如CWT+CARS组合在VNIR子域的类胡萝卜素反演, RP2相对提高了14.9%, RMSEP降低6.3%。 这表明模型反演受光谱子域的影响。 但研究发现, f在整个光谱域均具有高灵敏度, 即使f仅在SWIR区域发生变化, 也会影响可见光区域中色素的反演^{[9, 11]}, 表明使用光谱子域不能充分考虑到所有参数对光谱的贡献。 这可能也是VNIR子域的模型未能比全谱域模型反演结果更好的原因。

2.3.5 FD与CWT对消除混合光谱中的f的有效性和适用性

Wan等^[11]研究结果表明: f对BRF谱的贡献可达20%~100%, 视不同谱域而定; 不同物种叶片的形状、 叶片表面的粗糙程度也会对f造成影响, 例如针叶植物比阔叶植物有更复杂的多向反射^[21]; 光滑表面会使叶片产生更高的表面反射率, 影响反演精度。 在本研究中, 通过消除表面反射率改善叶绿素、 类胡萝卜素的反演结果, 可有效减少f对模型反演的影响。 然而, 本方法不能直接用于光谱模拟, 而是通过缩小DHRF和BRF之间的差异来提高叶片生化参数的反演精度。 除此之外, 叶片表面反射率影响着叶子在不同观测角度的定向反射率的大小和分布。 导数和小波变换处理直接将光谱中的表面反射率简化为混合光谱间的差异f[见式(4)], 而未考虑BRF光谱的多角度问题, 对叶片表面特性的表征不充分, 尚不能完全消除f对BRF的影响^[22]。 如何进一步提高反演精度, 还有待进一步深入研究。