磴口县荒漠绿洲区(40° 9'— 40° 57'N, 106° 9'— 107° 10'E)位于内蒙古巴彦淖尔市西南部、 乌兰布和沙漠东北缘, 内蒙古河套平原源头, 黄河中上游, 背靠狼山山脉, 研究区总面积为2 887.98 km²。 该区域属于温带大陆性季风气候, 四季分明, 年均降雨量143.4 mm, 年均蒸发量2 387.6 mm, 年均风速2.7 m· s^-1。 灌木植被以油蒿、 白刺、 梭梭、 沙冬青、 柠条为主; 乔木植被以沙枣、 新疆杨、 旱柳、 樟子松、 胡杨、 欧美杨为主; 草本植物有芦苇、 蓼子朴、 猪毛菜等。 磴口县荒漠绿洲防护林既有大网格宽林带, 也有近些年建造的小网格窄林带, 其空间分布较自然林更为规则, 通常处在道路、 农田、 水系和房屋周边。 不同防护林带的长度、 宽度差异很大, 但均为条带状分布。 本研究选取磴口县荒漠绿洲区中约27 km²区域, 为图1中红色矩形区域, 用于最优分割参数的选取, 该区域同时包含多种防护林类型。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 研究区位置及验证样点空间分布图Fig.1 Location of the study area and the spatial distribution of sample points for verification

研究采用的遥感影像数据为GF-2影像, 包括4 m分辨率的多光谱影像和1 m分辨率的全色影像。 多光谱影像包含红、 绿、 蓝和近红外四个波段。 由于研究区面积较大, 影像获取难度高, 因此, 本研究镶嵌22景共三个时期的遥感影像作为研究区影像数据。 三个时期影像获取时间分别为2021年5月、 2020年7月和8月。 2020和2021年防护林空间分布变化不大, 且选择的月份, 防护林已过发芽期, 生长出叶片, 利于防护林的提取。 对于22景遥感影像, 首先进行辐射校正、 几何校正和影像融合操作, 经裁剪得到研究区影像^[8]。 然后, 结合实地调查及2019年7月— 8月份磴口县沙林中心提供的土地利用数据, 目视解译形成参考数据集, 用作分类器训练样本。 还对研究区实地调查367个防护林样地, 用于对识别结果的精度验证。

首先基于ROC-LV曲线与分割对比实验在27 km²小研究区获取最优分割参数, 随后应用于整体研究区影像分割, 根据分割结果进行特征空间优化, 然后优化分类器并训练, 最后得出分类结果与分类精度, 图2为本研究技术路线图。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 技术路线图Fig.2 Technology roadmap

1.3.1 影像分割

采用多尺度分割(multiresolution segmentation)方法进行多光谱影像分割, 目的是将影像分割成无数个矢量对象, 进而分别对这些矢量对象分类。 此过程分割结果由分割尺度(scale)、 形状(shape)、 紧凑度(compactness)三个参数确定。 分割尺度用来确定分割对象所允许的最大异质度, 值越大生成影像对象尺寸越大^{[9, 10, 11]}。

经过多次试验及多尺度分割参数的定义可知, 研究区大小不影响最优分割参数的确定, 故选取包含多种防护林种类的小区块作为确定最优分割参数的研究范围。 由于防护林形状不规则, 基于潜在分割误差、 分割强度、 欧式距离三个指标确定最优分割参数并不适用, 而ROC-LV(局部方差变化率)曲线方法广泛适用于不规则地物的分割, 故采用ROC-LV曲线方法进行影像多尺度分割。 首先, 通过统计影像同质性的局部方差LV(local variance)来表示分割结果内的标准差, 原则上LV会随分割尺度增大而增大。 然而, 只有LV不足以选择最优分割尺度(scale), 故引入LV的变化率ROC(rate of change), 绘制ROC-LV曲线图, 当ROC达到峰值时, 所对应的分割尺度极有可能为最优分割尺度^[12]。 ROC计算公式为

ROC=Li+1-LiLi×100%(1)

式(1)中, L_i+1为第i+1层对象的平均标准差, L_i为第i层中对象的平均标准差。

1.3.2 特征选择

采用的分类类别体系为防护林、 农田、 水体、 裸地和建设用地。 考虑防护林分布形状、 光谱特点, 共选取4大类共31个特征作为待筛选的特征空间, 用于影像提取防护林研究(表1)。 光谱特征反映了各对象间光谱信息的差异, 几何特征反应各对象间几何、 形状、 大小的差异。 纹理特征通过灰度共生矩阵(gray-level co-occurrence matrix, GLCM)提取, 描述具有某种空间位置关系两个像素灰度的联合分布。 通过对灰度共生矩阵(GLCM)的二次统计, 可以得到如对比度、 熵等的二次统计量, 来反映如防护林对象包括影子在内的纹理信息。

表1

Table 1

表1(Table 1)

表1 分类器待筛选的特征列表 Table 1 List of features to be selected for the classifier 光谱特征 几何特征 纹理特征 其他特征 均值 面积 熵 归一化植被指数 标准差 长度 角二阶矩 波段比 宽度 同质性 亮度 长宽比 异质性 密度 对比度 紧致度 平均值 圆度 标准差 边界指数 相关度 非对称性

表1 分类器待筛选的特征列表 Table 1 List of features to be selected for the classifier

选择的纹理特征为8个, 包括对比度(contrast, Con)、 熵(entropy, Ent)、 角二阶矩(ang.2^nd moment, ASM)、 同质性(homogeneity, Homo)、 均值(mean)、 异质性(dissimilarity, DISL)、 标准差(std dev)和相关度(correlation, Corr)。 这些纹理特征反映了对象图像清晰度和沟纹深浅、 灰度分布的复杂程度、 均匀程度和纹理粗细度、 灰度变化、 灰度值之间的稳定性与灰度关系的相似度。 其他特征选择了归一化植被指数(NDVI), 已有研究表明NDVI很好地反映了植被生长状态, 对区分植被与非植被有很大作用。

ASM=∑i=1n∑j=1nP2(i, j)(2)

Con=∑i=1n∑j=1n(i-j)2P(i, j)(3)

Ent=∑i=1n∑k=1nP(i, j)lnP(i, j)(4)

Homo=∑i=1n∑j=1n11+(i-j)2P(i, j)(5)

Mean=1n∑i=1n∑j=1nP(i, j)(6)

DISL=∑i=1n∑j=1n(i-j)P(i, j)(7)

Std=1n-1∑i=1n∑j=1n(xi, j-x̅)2P(i, j))(8)

Corr=∑i=1n∑j=1nijP(i, j)-μxμyσxσy(9)

NDVI=NIR-RNIR+R(10)

式中, σ _x和μ _x分别为P_x(i)的方差和均值; σ _y和μ _y分别为P_y(j)的方差和均值; P(i, j)为灰度共生矩阵; NIR和R分别为近红外波段反射率和红波段反射率。

多光谱影像的纹理特征、 光谱特征、 几何特征繁多, 形成了高维特征空间, 但特征并非越多越好, 大量冗余的特征不仅会增加计算机运行负担而且会造成“ Hughes现象” ^[13]。 因此, 需要对特征进行筛选。 基于随机森林算法获取特征重要性评分来筛选特征, 通过袋外误差率控制特征个数, 达到优化特征空间的目的。

随机森林算法是一种Bagging算法, 由Breiman提出^[14], 以决策树为基础, 通过Bootstrap自动抽样方法随机抽取样本子集, 未被选中的训练样本作为袋外数据(out of bag, OOB)用于验证模型。 袋外数据(OOB)约占总训练数据的1/3, 袋外数据误分个数占总样本总数的比率称为随机森林的OOB误分率^[15]。 随机森林计算特征重要性时方法主要有袋外误差率法和基尼系数法^[16]。 基尼系数用来表示一个随机选中的样本在其子集中被错分的可能性, 即计算样本集合的纯度和不确定性, 其定义如式(11)

Gini(p)=∑i=1npi(1-pi)=1-∑i=1npi2(11)

基尼系数法是在构建随机森林结构时计算得出特征重要性, 与每棵决策树的节点分裂有关, 而袋外误差率法使用OOB数据, 通过给特征加上随机噪音, 比较前后误差率得到特征重要性。 综合应用基尼指数(Gini)的下降量与OOB误分率的增加量来优化特征空间。

1.3.3 分类器模型与训练

分别采用支持向量机(SVM)^[17]、 CART决策树、 随机森林(RF)和K-近邻算法(KNN)四种分类器进行分类对比。 对防护林、 农田、 水体、 裸地和建设用地5个类别进行目视解译人工取样, 取样样本在研究区内随机分布, 类别与样本数量关系如表2所示, 所选取样本数量均大于特征数量, 样本数量充足。

表2

Table 2

表2(Table 2)

表2 分类类别与样本数 Table 2 Number of samples for each category 识别对象类别 样本数量 防护林 3 126 农田 573 水体 394 裸地 1 797 建设用地 958

表2 分类类别与样本数 Table 2 Number of samples for each category

随机森林是利用自助法重采样技术(Bootstrap), 从原始训练样本集中有放回地抽取训练样本, 并且决策树节点随机抽取特征, 构造多个互不相同且相互独立的决策树。 在判断新样本类别时, 根据所有决策树投票类别的众数, 进而决定最终类别^[18]。 随机森林的重要参数包括节点特征变量个数mtry和决策树的数量ntree。 Breiman建议对于分类问题, 二叉树节点使用的特征变量个数应该为总变量个数的平方根。 为提升分类效率, 应尽量减少决策树的数量。 本研究最优决策树数量(ntree)采用OOB误分率进行优化(图3), 在ntree约大于400后, 各项OOB误分率基本保持平稳, 总体误差在4%左右, 最小总体误差为ntree=369时取得。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 决策树数量(ntree)与OOB误分率(out of bag error)关系图Fig.3 The relationship between the number of decision tree (ntree) and out of bag (OOB) error

支持向量机(SVM)是一种机器学习算法, 在解决小样本、 非线性及高维模式识别中表现出运算速度快、 精度高、 泛化能力强的优点, 并能够推广应用到函数拟合等其他机器学习问题。 SVM分类器采用高斯核函数[式(12)]时有两个主要参数分别为惩罚系数(C)和Gamma系数(γ )。 本研究基于K-fold cross validation方法进行调参。 经过调整参数后, SVM算法最优Gamma=0.1, C=10, 训练数据的最优误差率约为3.2%。

k(xi, xj)=exp-‖xi-xj‖2δ2(12)

CART决策树是一种后剪枝决策树, 将训练样本划分为测试变量和目标变量, 由初始节点依次形成二叉决策树, 采用基尼系数(Gini)作为选择最佳分类特征和分割阈值的标准。 KNN最近邻分类算法以所有已知类别的样本作为参照, 计算未知样本与所有已知样本的距离, 选出距离最近的部分已知样本, 与类别占比较多的归为一类。 为了消除样本数量对分类器训练的影响, 对训练样本进行了标准化处理。

1.3.4 精度评价

为了评估各分类器对于防护林的识别准确度, 基于野外采集样点数据, 采用TTA mask方法计算混淆矩阵(confusion matrix)进行精度评价。 混淆矩阵(confusion matrix)法是通过对每个样本在地物实际的类别与其在分类结果影像中对应的类别进行比较并进行统计与计算的一种精度验证方法。 利用混淆矩阵可以计算出评价单个类别分类效果的指标包括用户精度(user’ s accuracy)[式(13)]、 生产者精度(producer’ s accuracy)[式(14)]与总体精度(overall accuracy)。

PAi=xiix+i(13)

UAi=xiixi+(14)

式中, x_ii为主对角线上的值, x_i+和x_+i为第i行和第i列的和。

遥感影像分割需要充分考虑分割对象的影像光谱特征。 在荒漠绿洲区域, 灌木、 乔木植被叶面普遍较小且分布稀疏, 部分植物体有白色绒毛, 有些单株生长, 在影像中光谱特征较弱, 对提取防护林造成很大干扰。 因此, 参考彭佳忆等^[19]有关荒漠植被提取的相关参数, 将形状参数设定为0.1, 紧致度设定为0.5。 初始Scale Parameter设为30, 步长设为2, 步频设为60。 由于本研究的分割用于提取荒漠绿洲防护林, 植被对近红外波段反映灵敏, 故在影像多尺度分割中将近红外波段权重设为2, 其余设为1。 先基于以上设置参数, 对GF-2号影像进行多尺度分割, 然后依据每个分割尺度对应的LV, 通过公式计算相应ROC值, 得到ROC-LV曲线图, 如图4所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 分割尺度与ROC-LV关系Fig.4 The relationship between Scale Parameter and ROC-LV

可以看出LV随分割尺度增大而增大。 因尺度越大, 分割出的对象面积越大, 相邻对象的差异越明显。 ROC总体呈下降趋势, 在下降过程中存在明显的峰值。 峰值的分割尺度分别为52, 68, 82, 104, 130和134, 利用其进行分割对比实验。 结合人工目视判断方法对ESP2得出的最优分割尺度可能值进行判别, 当尺度参数为82时, 防护林提取的效果最好。 因此, 研究区的最优分割参数组合设为Scale=82, Shape=0.1, Compactness=0.5。

在采用所有分类器分类之前, 本研究基于随机森林(RF)算法计算待筛选特征的特征重要性, 以提取重要特征提高分类效率与精度。 从图5中可以看出, 随着特征数量的增加, OOB误分率呈逐渐下降的趋势。 这主要是由于多种特征的加入增加了特征空间的维度, 提高了分类器性能。 当特征数量增加至17之后, OOB误分率趋于平缓, 此时分类器性能已经达到饱和。 随特征数量继续增加, 并没有出现明显的“ Hughes” 现象, 即冗余特征导致分类精度下降的现象, 这可能是由于特征数量不够多或者特征与特征之间相关性较强所致。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 特征数量与OOB误分率的关系Fig.5 The relationship between the number of features and OOB error

从31个待选特征中, 综合利用Gini系数及袋外数据误差评估特征重要性。 如表3所示。 从重要性排序可以看出波段比率特征(Ratio)与纹理特征综合重要性排名靠前, 表明这两类是分类器的重要特征, 对提取防护林有重要意义。 将重要性排名前17位的特征参与所有分类器计算, 这些特征中包含了4个波段的波段比及蓝、 近红外波段的均值, 这也说明了光谱特征是区分地类的重要标志。 并且, 边界指数、 宽度、 紧凑度、 密度、 对称度这些形状指数有效地区分了长条状林带、 类似长方形农田建设用地, 与不规则的水体及裸土。 在研究过程中, 通过特征空间优化, 显著减少了训练耗时, 并且通过优化决策树数量从1 000个降至369个, 随机森林训练器训练与分类总耗时从23 h缩减到了12 h左右, 提升了近一倍的分类效率。

表3

Table 3

表3(Table 3)

表3 Gini系数法与OOB误差法综合重要性得分排名前10的特征 Table 3 Top 10 features of overall importance ranking by Gini index method and OOB error method 特征 Gini系数法 重要性 OOB误差法 重要性 综合重要性 排名 Ratio BLUE 455.5 43.7 1 GLCM Correlation 223.3 44.6 2 GLCM Contrast 354.9 34.3 3 Ratio NIR 368.1 28.9 4 Standard deviation NIR 280.8 28.7 5 Border index 148.3 40.1 6 Width 149.9 34.3 7 Compactness 114.4 44.6 8 Mean BLUE 324.2 23.7 9 Density 109.0 34.5 10

表3 Gini系数法与OOB误差法综合重要性得分排名前10的特征 Table 3 Top 10 features of overall importance ranking by Gini index method and OOB error method

结合实测数据得到367个包含防护林及其他类别的空间分布数据, 通过构建TTA mask得到分类结果混淆矩阵, 经计算得到用户精度与生产者精度(表4)。 对于防护林总体识别精度, 四种分类器均得到了较高的分类精度, 均大于94%。 将CART、 随机森林、 支持向量机、 K-近邻法验证结果进行比较, 随机森林分类器得到的用户精度为四种分类器中最高, 达到93.5%, 说明随机森林分类器的分类存在更少的错分现象。 支持向量机法与K近邻法生产者精度均为97.14%, 说明其分类正确的防护林数量最多。 在两种基于决策树的算法中, 随机森林算法生产者精度与用户精度和CART决策树算法相近, 用户精度与SVM和KNN相差不多。 整体而言, SVM与KNN总体分类精度最高, 均为96.19%, 是更理想的分类器。

表4

Table 4

表4(Table 4)

表4 CART, RF, SVM, KNN分类精度计算结果 Table 4 Classification’ s accuracy results of CART, RF, SVM and KNN classifier 分类器 精度指标 生产者精度/% 用户精度/% 总体精度/% CART 93.6 92.9 94.8 RF 92.9 93.5 94.8 SVM 97.1 93.2 96.2 KNN 97.1 93.2 96.2

表4 CART, RF, SVM, KNN分类精度计算结果 Table 4 Classification’ s accuracy results of CART, RF, SVM and KNN classifier

选取4块典型小区块, 将RF, CART, SVM和KNN的分类结果与原始遥感影像图进行对比分析(图6)。 其中, 区块一为靠近城区的农田防护林网, 存在林带断带。 其分类结果显示, SVM(a), RF(b)和KNN(d)的防护林提取效果均远好于CART(c), 相较于SVM, KNN的错分现象更多, RF(b)存在一定的防护林带漏分, SVM结果中较好地显示了断带信息。 区块二为距城区较远的农田防护林网, 周边为裸露沙地, 部分林带树间距较大, 防护林稀疏。 在区块二内, 四种分类器均提取出了零散防护林。 SVM(a)与RF(b)的提取效果更好于CART(c)与KNN(d), 但SVM将部分细长型田间道路也划分为防护林带, 存在错分现象。 RF的林带提取结果中道路很少, 用户精度较高。 区块三为纳林湖景区的一部分。 在纳林湖中有大面积的芦苇(图中左下)。 四种分类器均能较好地区分水体与非水体, 而较难对芦苇与乔木防护林进行区分, 其原因可能是芦苇与防护林的光谱特征和纹理特征相近。 区块四为磴口县城区, 在道路两旁有行道树种植。 四种分类器均能提取城市道路防护林, 但也有将房屋阴影错分为林带的现象, RF(b)的分类结果优于其他分类器, 错分现象最少。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 四个区块的防护林分类结果 (a): SVM; (b): RF; (c): CART; (d): KNNFig.6 Shelterbelt’ s classification results in four regions (a): SVM; (b): RF; (c): CART; (d): KNN

为深入探究各分类器对芦苇与防护林错分原因, 分别对芦苇与防护林带提取特征均值并归一化, 进行对比分析。 在各分类器的分类结果中对芦苇和防护林分别采样300个样本, 计算样本的31个特征均值, 由图7可知, 四种分类器应用波段比特征、 相关度和熵等特征较难区分芦苇和防护林, 而这些特征较多是特征重要性综合排名较高的特征。 斑块长宽比、 灰度共生矩阵角二阶矩、 对比度和同质性是芦苇和防护林斑块相差较多的特征, 这些特征的特征重要性综合排名较低, 没有被选入特征空间。 对比四个分类器, KNN分类器提取所得的防护林特征值均较高, 而芦苇的特征值较其他分类器较低, 四种分类器均具有区分芦苇和防护林的潜力。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 芦苇与防护林带特征值对比 (a): SVM; (b): RF; (c): CART; (d): KNNFig.7 Feature value comparison of reeds and shelterbelts (a): SVM; (b): RF; (c): CART; (d): KNN

由分类结果的防护林分布图(图8)可知, 磴口县防护林主要分布在东部黄河流经地区及城镇周围, 以及中部地区及农田分布较密集的区域, 而西北部山区及南部乌兰布和沙漠几乎没有防护林存在, 显示出由荒漠到绿洲到城镇的防护林体系空间递变规律, 以及防护林与水系的关系, 其为维持荒漠绿洲稳定性提供了较大帮助。 四种分类器中K近邻分类器(d)得到的防护林面积最大, 为224.95 km², 随机森林算法(b)得到的防护林面积最小, 为143.89 km², SVM与KNN得到防护林带面积均在200 km²左右。 RF与CART两种基于决策树算法的分类器结果的边长面积比接近, 均小于SVM与KNN, 说明基于决策树算法的分类器更倾向于提取形状较圆润规则的斑块作为防护林。 根据2018年磴口县统计公报, 磴口县森林覆盖率为20.6%, 且灌木乔木与其他类型林地占比约为2:2:1, 乔木林地为防护林的主要林种, 去除西部阴山部分面积后, 其面积约为210 km²。 本研究得到2020年防护林(大多数为高大乔木, 包含小部分高植株灌木)面积约为200 km², 与实际情况相符。

图8

Fig.8

	Figure Option ViewDownloadNew Window
	图8 研究区内四种分类器的防护林分布图 (a): SVM; (b): RF; (c): CART; (d): KNNFig.8 Shelterbelts distributions in study area by different classifiers (a): SVM; (b): RF; (c): CART; (d): KNN

纵观整体磴口县分类结果, 综合分类精度、 小区块和宏观面积考虑, SVM分类器展示了较好的荒漠绿洲防护林提取能力, 其分类精度高, 微观小区块内准确度高, 总体防护林面积接近于实际面积, 为最优的防护林分类器。 但本研究基于上述方法, 对于防护林的提取效果仍可能存在不足, 例如, 对于边界规则的大面积防护林地, 大多数用于提供树苗或经济林, 其树高与农田中植株高的作物相似, 可能导致其误分为农田; 对于边界不规则防护林地, 其形状特征不明显, 内部结构杂乱, 并且树间距较大, 露出林下地表, 使得光谱特征存在干扰; 对于芦苇与防护林的混淆问题, 文中已做了分析, 若将芦苇作为单独的地类, 则分类效果可能会更好。 此外, 对于大范围防护林提取, 由于树种、 树龄以及防护林种植数量差异, 使得距离较远的两个地区防护林光谱特征、 形状特征、 纹理特征相差较大。 因此, 样本选取必需遵循多区域、 广覆盖、 多样本的规则, 随着未来传感器与航天平台的发展, 更高分辨率的多波段大尺度影像将为提取防护林提供更多可能。 在本研究中, 特征空间优化与参数优化有效提升了荒漠绿洲防护林提取精度与效率, 但在模型训练与应用时, 数据处理仍耗费大量时间(20 h左右)。 未来硬软件技术的发展, 将使识别提取防护林空间信息更加快速便捷。