高空间分辨率遥感影像最优分割结果自动确定方法

程结海，黄中意，王建如，何 湜

河南理工大学测绘与国土信息工程学院，河南 焦作 454000

基于像元的影像分析方法一直是遥感界进行遥感影像分析的主流方法[1]。近些年，随着遥感影像空间分辨率不断提升，基于像元的影像分析方法受到越来越多的质疑[2-4]。原因之一在于高空间分辨率遥感影像中类内光谱差异增大[5-7]，影响遥感影像信息提取精度[8]。采用面向地理对象影像分析(GEOBIA)方法可以有效解决上述问题[1,9-10]。影像分割是GEOBIA中一个关键步骤，即：将影像分割成相对均质的区域，该均质区域被称作“分割对象”或“影像对象”[11]。分割对象取代像元成为影像后期分析的基本单元，其主要优势有：①减小类内光谱差异，且能融合形状、纹理和上下文信息等，使得GEOBIA具有处理复杂影像分析任务的潜力[12-14]；②克服传统分类方法易出现的椒盐现象[15]。在GEOBIA框架下，目前主流的影像分割算法有多尺度分割[11]、均值漂移分割[16]和分水岭分割[17]等。这些方法均通过设置分割参数来生成不同质量的影像分割结果。影像分割结果质量直接决定影像后期分类质量[18]。因此，确定影像最优分割结果至关重要[19]。

目前，确定影像最优分割结果方法主要分为监督评价法和非监督评价法[20]。监督评价法人为干预较多，费时费力，自动化程度普遍不高[21]。非监督评价法不需要选取参考对象，被认为是自动确定影像最优分割结果最具有潜力的方法[21-22]，该类方法采用以下标准：①同一分割对象内部是光谱均一的；②相邻分割对象间有明显光谱差异；③分割对象内部没有孔洞；④分割对象具有规则整齐的边界。由于高空间分辨率遥感影像普遍具有较强纹理信息，遥感影像对应的下垫面地物往往边界形状各异，使得标准③和标准④不具有适用性。对于遥感影像来说，只有标准①和标准②是符合实际需求的[21]。按照标准①和标准②，文献[10,13—14,22—26]采用方差表达分割对象内部光谱均质性，文献[10,13,23—25]采用Moran 's I空间自相关度量指标[27]表达分割对象与其相邻分割对象间光谱异质性。但方差和Moran 's I指标是在影像单一波段上完成求算，而高空间分辨率遥感影像往往存在多个波段数据。虽然文献[10,13,22]采用分波段求算，然后对所有波段上求算的方差和Moran 's I值分别取算术平均的方法来解决上述不足。但是，取平均值做法本身掩盖了波段间光谱差异，特别是对于下垫面复杂的城市区遥感影像来说，这种做法更不利于客观表达分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性。文献[5，18]提出了采用光谱角测量指标表达多波段遥感影像分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性的方法，但光谱角测量指标本身只对光谱曲线形状较敏感，对光谱距离特征不具有处理的能力，且该方法确定出的影像最优分割结果只是在某一分割参数下影像整体的最优分割结果，而影像整体的最优分割结果并不能确保影像范围内各类型地物都能同时达到各自的最优分割状态，经常会出现影像整体达到所谓的最优分割结果，而影像中某一些地物对象却是处于过分割或欠分割状态。这主要与不同类型地物的光谱特征、形状和尺寸大小等差异有关。

针对现有方法不能充分顾及遥感影像多波段光谱信息，以及确定出的影像最优分割结果只是在某一分割参数下影像整体最优分割结果，而忽视遥感影像中地理要素的多尺度特性等问题，本文提出了一种自动确定高空间分辨率遥感影像最优分割结果的非监督评价方法。

1 本文方法

图1 本文方法总体流程Fig.1 The flowchart of the proposed method

1.1 分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性的表达

分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性这两项指标用以表征影像分割结果质量，其计算如下。

(1) 定义分割对象内部两不同像元X和Y的光谱矢量

X=(x1,x1,…,xL)T

(1)

Y=(y1,y1,…,yL)T

(2)

式中，L是波段数；xl和yl表示像元X和Y在第l波段上灰度值(l=1，2，…，L)；上标T表示转置操作。

(2) 在光谱矢量基础上，定义像元X和Y的概率矢量

(3)

(4)

(3) 计算像元X和Y的信息熵

(5)

(6)

(4) 在信息熵基础上，计算像元X和Y之间相对熵

(7)

(8)

(5) 将像元X和Y之间相对熵相加，得到像元X和Y之间光谱信息离散度

SID(X,Y)=D(X‖Y)+D(Y‖X)

(9)

(6) 计算分割对象内部任意两个不同像元间光谱信息离散度，将得到的光谱信息离散度进行平均，即为该分割对象内部光谱均质性

(10)

式中，P是分割对象内部不同像元对的集合；I是P集合中像元对的数量。

(7) 利用分割对象内像元光谱矢量均值计算分割对象与其相邻的任一分割对象间光谱信息离散度，再根据分割对象与其相邻分割对象间公共边长度设置权重，得到的光谱信息离散度进行加权平均，来表达分割对象与其相邻分割对象间光谱异质性

(11)

1.2 影像分割结果质量的表达

分别基于θ和φ指标，以分割对象面积除以影像总面积作为权重，得到整个影像分割结果光谱均质性Θ和异质性Φ

(12)

(13)

式中，n为影像分割结果中分割对象的总个数；θi和φi分别为第i个分割对象内部光谱均质性和第i个分割对象与其相邻分割对象间光谱异质性；ωi为权重，即第i个分割对象面积占影像总面积的比例。

联合Θ和Φ，构建SQI指标，用于表达影像分割结果质量

SQI=e-Θ+eΦ

(14)

式中，SQI值越大，表示影像分割结果质量越好；反之，亦然。

1.3 影像最优分割结果的粗估计

极端欠分割结果对于找寻影像最优分割结果会造成干扰，没有实际意义。为此，构建指标剔除极端欠分割结果，以提高整个工作效率和精度

(15)

式中，sp表示分割参数；Δsp表示分割参数的增量；Θsp+Δsp表示影像在sp+Δsp分割参数下的影像分割结果的光谱均质性Θ值；Θsp表示影像在sp分割参数下的影像分割结果的光谱均质性Θ值；Φsp+Δsp表示影像在sp+Δsp分割参数下的影像分割结果的光谱异质性Φ值；Φsp表示影像在sp分割参数下的影像分割结果的光谱异质性Φ值。

1.4 影像最优分割结果的精确定

由于影像上地物尺寸和类型存在差异，R1st并不能确保各类地物对象同时达到各自的最优分割状态，往往会出现某些地物分割对象仍处于过分割或欠分割状态。为此，基于θ和φ指标，构建HI来标定处于过分割或欠分割状态的分割对象

HI=(θ×φ)/(θ+φ)

(16)

式中，HI越大，表示分割对象越倾向于欠分割状态；HI越小，表示分割对象越倾向于过分割状态。

对R1st中所有分割对象按HI值大小进行排序，按5%、10%、15%和20%比例从中选取HI值较大的欠分割对象进行再分割优化。将再分割优化后的分割对象与R1st中未参与再分割的分割对象进行联合，得到一系列优化后的影像分割结果，其中，SQI值最大的影像分割结果为该次优化后的影像最优分割结果，记为R2nd。

对R2nd中所有分割对象按HI值大小进行排序，按5%、10%、15%和20%比例从中选取HI值较小的过分割对象进行合并优化。将合并优化后的分割对象与R2nd中未参与合并的分割对象进行联合，得到一系列优化后的影像分割结果，其中，SQI值最大的影像分割结果为该次优化后的影像最优分割结果，即全面影像最优分割结果，记为R3rd。

根据商保公司经办城乡居民基本医保模式主要要素，总结分析商保机构经办城乡居民基本医保的试点情况。依据理论研究和实践结果，目前安徽省商保公司经办城乡居民基本医保试点初步形成的模式主要是委托经办模式，即：政府将业务委托商保公司，按照协议由后者获取适当的经办费用，提供经办服务。在委托经办模式中，按照经办费用来源已初步形成委托经办支付费用模式、委托经办不支付费用模式两种。每一大类中，按照经办服务内容、经办费用获取方式、基金风险承担实际运行情况又可分不同的模式。见表2

在上面的操作中，关键点是确定用于再分割和合并的分割参数值，其具体过程如下。

(1) 定分割参数范围。标记R1st对应的SQI值为SQIR1st，对应的分割参数为spR1st。由于SQIR1st为局部极大值，因此将SQIR1st×99.9%值对应的两个分割参数spL和spR构成的区间[spL，spR]确定为分割参数范围。其中，[spL，spR1st)范围内分割参数用于再分割操作，(spR1st，spR]范围内分割参数用于合并操作。若spL小于分割算法的最小分割参数，则spL直接取最小分割参数；若spR大于极端欠分割结果对应的分割参数spMax，则spR直接取spMax。

(2) 确定具体分割参数值。具体分割参数值的确定用式(17)—式(19)来完成

(17)

(18)

(19)

一个Isp值对应一个分割参数，用该分割参数对影像进行分割，可使影像上某一类地物达到最优分割状态。考虑到影像上有不同尺寸和类型地物，故将[spL，spR1st)和(spR1st，spR]分割参数范围内所有Isp值对应的分割参数全部确定出来，分别用于再分割和合并操作。

2 试验与分析

2.1 影像及分割

图2 研究区影像及影像欠分割结果Fig.2 The study area image and image under-segmentation results

2.2 试验结果与分析

将分割参数scale以50为间隔，在区间(0,500]内对影像进行分割、计算，得到A1、A2和A33个区域值随分割参数scale变化的第1个拐点spMax分别是150、250和200，对应的影像欠分割结果如图2所示。由图2可以看出，有大量的分割对象处于欠分割状态。若继续增大分割参数scale，这种欠分割现象将越来越明显。因此，将A1、A2和A33个区域的分割参数scale分别限定在(0,150]、(0,250]、(0,200]范围内，从而剔除极端欠分割结果。分别在3个区间范围内，以5为间隔设定分割参数scale，分别对3个区域影像进行分割，并计算各影像分割结果的SQI值(图3)。

图3中虚垂线与曲线相交的点分别为3个区域SQI最大值点，其对应的scale分割参数分别为80，90和180，对应的影像分割结果即为R1st(图4)。

图3 不同分割参数scale下影像分割结果的SQI值Fig.3 The SQI values of image segmentation results under different segmentation parameter scales

对R1st中所有分割对象按HI值进行排序，从中选取HI值较大的5%、10%、15%、20%数量比例的欠分割对象。图5为从R1st中选取HI值较大的15%数量比例的欠分割对象。由图5可以看出，欠分割对象大部分集中在林地与水体、裸地或道路与裸地混交区域。此类混交区域的土地利用类型光谱特征比较接近，容易被分割成同一个分割对象。

图4 初步影像最优分割结果R1stFig.4 The preliminary optimal-image segmentation results R1st

图5 选取的欠分割对象Fig.5 The extracted under-segmented objects

根据图3结果，R1st中欠分割对象的再分割参数范围[spL，spR1st)分别为[35,80)，[35,90)和[135,180)。求算再分割参数范围内所有Isp值对应的分割参数，见表1。

表1 各研究区用于再分割和合并的scale分割参数值

利用表1中再分割参数分别对R1st中欠分割对象进行再分割处理，并分别与R1st中未参与再分割的分割对象进行联合，得到一系列优化后的影像分割结果，分别为12个、16个和16个，其中，SQI值最大的影像分割结果即为一次优化后影像最优分割结果R2nd，如图6所示。

对R2nd中所有分割对象按HI值进行排序，从中选取HI值较小的5%、10%、15%、20%数量比例的过分割对象。图7为从R2nd中选取HI值较小的15%数量比例的过分割对象。由图7可以看出，过分割对象主要是道路区域和部分建筑物。原因在于由于道路本身尺寸比较大，且道路内部往往存在一些类似植被隔离带等地类，使得在实际影像分割中道路往往被分割成多个分割对象；在实际场景中，由于脊线和光照差异等原因，建筑物顶部经常被分割成多个分割对象。图7(c)中被挑选出的过分割对象数量比较少，这与该影像分割结果中分割对象总数量少有关。同时，观察发现图7(c)中所谓的过分割对象实际上并不处于过分割状态，而更倾向处于一种最优分割状态。原因在于对于形状相对规整和光谱均一的农用地区域来说，经过前面优化后的影像最优分割结果，其所有分割对象均倾向处于最优分割状态。在此情况下，按照HI值从小到大顺序规则来选取过分割对象，往往存在被挑选出的所谓过分割对象实际上并不处于过分割状态的现象。但这并不影响整个方法的有效性，原因在于这些被挑选出的过分割对象参与后期的合并优化，不会降低影像分割结果质量。这与本文所采取的合并优化方案有关，即对过分割对象间进行合并操作，并将合并后的结果与R2nd中未参与合并的分割对象进行联合，得到优化后的影像分割结果。类似图7(c)中被挑选出的过分割对象周围并没有其他过分割对象，使得这些过分割对象没有进行合并操作处理，进而保证了这些分割对象仍然处于各自的最优分割状态。

图6 一次优化后影像最优分割结果R2ndFig.6 The optimal image segmentation results after the first optimization R2nd

图7 选取的过分割对象Fig.7 The extracted over-segmented objects

根据图3结果，R2nd中过分割对象的合并参数范围(spR1st，spR]分别为(80,150]、(90,150]和(180,200]。求算合并参数范围内所有Isp值对应的分割参数，见表1。利用表1中合并参数对R2nd中过分割对象进行合并处理，并分别与R2nd中未参与合并的分割对象进行联合，得到一系列优化后的影像分割结果，分别为16个、16个和4个，其中，SQI值最大的影像分割结果即为二次优化后影像最优分割结果R3rd，如图8所示。

2.3 对比验证

采用文献[29]提出的ADI和PDI指标以及本文构建的SQI指标，对比本文方法与文献[18]方法、文献[10]方法和文献[13]方法确定出的影像最优分割结果质量。ADI和PDI指标值越小，表明影像分割结果质量越好。SQI指标值越大，表明影像分割结果质量越好。影像参考分割结果通过目视解译方法进行确定。对比结果见表2。

由表2可以看出，相比于文献[18]方法，本文方法所确定出的各级影像最优分割结果不论是ADI和PDI指标，还是SQI指标，均有更好的表现，特别是全面影像最优分割结果R3rd指标值基本上都优于文献[18]。虽然在A2研究区，文献[18]的PDI指标值优于本文方法的3个PDI指标值，但其ADI指标值和SQI指标值均比本文方法的对应指标值差。这些均表明，本文采用基于信息熵生成光谱信息离散度构建相应指标来表达分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性的做法要优于文献[18]采用光谱角测量指标的做法，本文方法更能充分顾及遥感影像多波段光谱信息。

图8 二次优化后影像最优分割结果R3rdFig.8 The optimal image segmentation results after the secondary optimization R3rd

表2 不同方法确定出的影像最优分割结果质量对比

与本文方法不同，文献[18]和文献[10]均没有采用“粗估计+精确定”的逐步多级优化策略，确定出的影像最优分割结果只是在某一分割参数下影像整体最优分割结果。对比本文方法中R1st、R2nd和R3rd发现，优化后指标值普遍要好于优化前指标值。进一步地，对比表2中R3rd与文献[18]和文献[10]的指标值可以看出：在A1和A2研究区，使用本文方法确定出的最终影像最优分割结果质量明显要好于另两种方法确定出的结果质量；在A3研究区，虽然R3rd的ADI指标值相比于另两种方法的ADI指标值不是最小的，但与最小ADI指标值比较接近，而且R3rd的PDI和SQI两指标值均比另两种方法相应指标值更优。试验结果验证了本文所采用的“粗估计+精确定”逐步多级优化策略的有效性。

与本文方法类似，文献[13]采用了逐步多级优化策略来提高影像最优分割结果质量。对比本文全面影像最优分割结果R3rd与文献[13]影像最优分割结果，由表2可看出：在A1研究区，无论是ADI和PDI指标，还是SQI指标，本文方法均更优；在A2和A3研究区，本文方法在PDI和SQI指标上表现更好。这进一步表明，本文充分顾及遥感影像多波段光谱信息进行影像最优分割结果确定的有效性和必要性。类似文献[13]，在处理遥感影像多波段光谱信息时，只简单地将各波段上求算的指标值取平均的方法势必会影响最终结果。同时，对比验证结果也证明了本文所采用的逐步多级优化策略更具有适用性。

3 结 论

本文方法基于信息熵生成光谱信息离散度，利用光谱信息离散度构建能表达分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性指标。基于分割对象内部光谱均质性和分割对象与其相邻分割对象间光谱异质性的指标，采取“粗估计+精确定”的策略，逐步得到一个多级优化后的影像最优分割结果。在整个影像最优分割结果确定过程中，不论是在粗估计阶段的初步影像最优分割结果确定，还是在精确定阶段的两级优化影像最优分割结果确定，本文方法都可以实现自动化操作，其有效性已经在不同下垫面影像区域得以证明。同时，与类似方法相比，本文方法确定出的影像最优分割结果质量更高，与参考分割结果更加贴近。在后期进一步研究中，可针对同一分割算法有多个分割参数情况，探讨本文方法的优化方案。