结合边界约束网络和分水岭分割算法的建筑物提取

罗壮，李明，张德朝

1.太原理工大学 数学学院,晋中 030600;2.太原理工大学 大数据学院,晋中 030600;3.清华大学 地球系统科学系,北京 100084

1 引 言

建筑信息广泛应用于灾害监测，数字城市建设和流行病监控等。高分辨率光学影像是建筑提取的重要数据源。高分辨率影像提供丰富的光谱特征，突出物体的结构、纹理等信息。同时高空间分辨率也带来了更多的图像噪声和“异物同谱”等问题（明冬萍等，2005）。此外因建筑结构的多样性和背景环境的复杂性，从高分辨率图像中准确检测建筑仍然是一个挑战。

中小尺寸的建筑物是城市的主体，在因人口聚集而土地稀缺的城市密集排列，常表现为大量大小相似的建筑物，位置非常接近。这类场景的遥感影像建筑物提取主要有两个问题：（1）精度不高，受建筑物结构的多样化和周边环境（特别是相邻建筑）的影响，容易误判，降低分类结果精度。（2）建筑边界不清晰，常出现多个建筑物“粘连”为一个。使用带掩模分支的目标识别框架（如Mask R-CNN（He 等，2017））和语义分割网络（如U-net（Ronneberger 等，2015））识别时各有不足。目标识别允许对象之间有重叠（Goldman等，2019），不适合建筑实例提取。不加约束和处理的单纯语义分割算法，常将多个建筑分割为一个区域对象。

建筑检测的算法和理论，可以分为传统面向对象的图像分析方法和深度学习方法。传统方法使用光谱、纹理（谭衢霖，2010；Wu等，2014）、形态轮廓（Lv 等，2014）等特征提取建筑，用多尺度方法（Tian 等，2016；Liu 等，2005）、先验知识（如阴影（Dosovitskiy 等，2015；Liow 和Pavlidis，1990））提升精度。相较于传统方法，深度学习的主要优势是学习特征表达的自适应性及较强的泛化能力。如DeepGlobe 挑战（Demir 等，2018）和INRIA 建筑检测比赛（Maggiori 等，2017） 的前5名方案均基于神经网络。当前建筑检测技术在精度上还有较大进步空间，特别是边界定位准确度需要提高。

针对密集建筑区域内，建筑提取精度不足和边界模糊两个问题，根据文献成果，加入边缘检测网络（Wu等，2018）或边缘检测分支（Marmanis等，2018），优化边界损失函数（Bokhovkin 和Burnaev，2019）可改进遥感影像的语义分割，受其启发，本文提出了一个完整的建筑检测方案。结合神经网络和自下而上分组的分水岭算法，对数据标签预处理，生成两类辅助标签-边界标签和相邻建筑分割线标签，在性能较优的建筑检测框架ICT-Net网络（Chatterjee 和Poullis，2019） 基础上改进，使用更适合建筑特性和3类标签的新损失函数，提升精度，以基于分水岭算法和梯度提升回归树（Ke等，2017）的后处理，从“粘连”的多个建筑掩模中，提取单体建筑。实验证明本文检测方案容易训练，既显著提高建筑检测精度，又使提取的建筑边界更准确。

2 方 法

本文的建筑检测总流程如图1 所示，分3 个阶段，数据预处理、监督分类和后处理。数据预处理就是转换原始建筑样本标签，生成建筑、边界和相邻建筑分割线3类标签。监督分类使用转换后的训练样本训练本文网络，预测3类别概率。后处理阶段，对网络预测结果应用基于分水岭算法和梯度提升回归树的后处理，从密集建筑掩模中提取建筑实例。

图1 建筑提取流程图Fig.1 The workflow of building detection

2.1 数据集与预处理

本文使用的数据源是INRIA 数据集中的芝加哥（Maggiori等，2017），共36张0.3 m空间分辨率的彩色正射校正航空影像和对应的建筑掩模标签，每景影像为5000像素×5000像素，面积约81 km2。

数据预处理就是对数据标签的预处理，形成3 类标签：原建筑标签、边界标签、相邻建筑分割线标签，前两类存在于每一个建筑，第3类分割线仅存在于紧凑相邻的建筑物之间。

图2为数据与预处理示例，每个建筑标签区域（建筑像素为1，非建筑像素为0）对应一个单体建筑，数据预处理的流程如下（其中形态学腐蚀尺寸、膨胀尺寸和邻域范围是根据图像特点设置的经验值）：

图2 3类标签Fig.2 The label of three categories

（1）对建筑标签应用形态学腐蚀操作（腐蚀尺寸设置为3×3），腐蚀后的建筑掩模与原建筑掩模的差，作为边界标签；

（2）对建筑标签应用形态学膨胀操作（膨胀尺寸设置为9×9），并减去原建筑标签区域，形成分割线标签的候选点集合；

（3）对原建筑标签做区域标号，遍历候选点集合，若该像素邻域范围（尺寸设置为7×7）内存在至少2 个不同建筑物，则该像素作为分割线标签。

2.2 网络结构与损失函数

本文考虑建筑物的特点，改进原ICT-Net网络。ICT-Net（Chatterjee 和Poullis，2019）属于完全卷积网络，基于U-net网络（Ronneberger等，2015），结合Dense blocks 与Squeeze-and-Excitation 的组合（这里简称FRDB 模块，图3 显示ICT-Net 网络结构，细节参考Chatterjee和Poullis（2019）。

图3 ICT-Net网络结构（无箭头虚线表示由FRDB，串联，下（上）采样组成的模块的重复）Fig.3 The structure of ICT-Net（where dotted line without arrow is repetition of the module consisting of FRDB，Concatenation and Down（Up）sampling）

首先修改网络最后一层为三通道，分别对应建筑、边界、分割线3 个类，应用Sigmoid 激活函数获得每个类的概率，并裁剪预测结果的4条边各24 个像素宽的边缘，以避免边界效应（Iglovikov等，2017），被裁剪的边缘不参与网络优化，减少后向传播的计算量。

其次，本文使用了更适合建筑提取和更有利于边界准确的损失函数，焦点损失（Lin 等，2017）、骰子损失和Jaccard 损失（Iglovikov 等，2017）。各损失计算公式如下：

式中，真实像素标签ygt∈{0，1}，预测概率ypd∈[0，1]。由计算公式，焦点损失从像素层次优化网络，骰子损失和Jaccard 损失作为IOU 微分形式的变体，从图像层次优化网络，根据Iglovikov 等（2017），使用3个损失加权比单个损失有优势，因此本文对每个类应用3个损失并加权，更好地优化网络。

语义分割领域常用的交叉熵损失，骰子损失等对线状对象（如边界）的错误匹配不敏感：当预测结果与真实标签偏移5—10个像素，上述损失数值无显著变化。为正确解释边界信息，指导网络优化，本文对边界类和分割线类添加边界损失（Bokhovkin 和Burnaev，2019），记作LBoundary，计算如下：

式中，掩模标签Ygt由元素ygt∈{0，1}构成，预测掩模Ypd由元素ypd∈[0，1]构成， *，*F表示矩阵点积求和运算，pool(*，pool_size)表示最大池化，pool_size为池化尺寸。建筑类、边界类、分割线类的损失分别是

式中，β1，β2，β3，β4为各损失项的权重。权衡建筑类和两个辅助类的重要性，网络总损失为

式中，α1，α2，α3表示3个类的权重。

2.3 后处理

后处理是将网络输出的建筑、分割线和边界3 类概率值转化为建筑区域标记的过程。首先合成3 类概率为一个值，使预测结果合成为一个灰度图，合成公式为

式中，Pdurban，Pdline，Pdborder∈[0，1]分别表示建筑概率、分割线概率和边界概率；(1 -Pdline)表示分割线约束，建筑不应与分割线重叠；(1 -kPdborder)表示边界约束，处于建筑边界的像素属于建筑的可能性较小，参数k控制边界约束的权重。

后处理也分为训练和预测两个部分。训练阶段的流程如下：

（1）对Gp图，以较小的阈值l0二值化，形成建筑标记L0（像素值大于l0标记为建筑，其余标为非建筑）。对建筑标记L0和Gp图采用基于标记的分水岭分割算法，获得预测图P0。

（2）在P0建筑标记范围内，分别以从小到大的阈值l1，l2，l3对Gp图二值化，形成初步建筑标记L1，L2，L3，并分别结合Gp图分水岭分割，形成相应的分割结果P1，P2，P3，都是在P0建筑标记范围内的建筑细分。

（3）对分割结果P0，P1，P2，P3中每一个连通区域（视作一个建筑），计算与真实标签的IOU=若IOU 小于0.5 的，IOU 值改为0；其中TP 表示正确检测（True Positive）；FP 表示误检（False Positive）；FN表示漏检（False Negative）。

（4）以上一步的建筑区域的IOU 为目标值，对每个建筑区域，分别经验性选取面积、面积周长比、长短轴、邻域信息等41 个常用区域及邻域相关的描述特征，训练梯度提升回归树。特征计算依赖于skimge 库（Van Der Walt 等，2014）的regionprops函数，该函数为每一个连通区域计算囊括上述特征的数十个属性值。

预测阶段，后处理流程如下：

（1）对Gp图，经训练阶段的第1—2 步处理，获得分割结果P0，P1，P2，P3，对结果的每一个连通区域，用训练阶段第4步的同样特征和训练好的回归树预测IOU。

（2） 对于P0的每一个连通区域，选取P0，P1，P2，P3中的IOU 平均值最大的分割结果，替换原连通区域，形成最终分割结果。

（3）对最终分割结果的每一个区域，若对应的IOU 大于接受阈值t1，或该区域概率均值大于概率阈值t2，则该区域标记为建筑。

3 实验结果与分析

本论文实验环境是Linux系统下以Tensorflow为后端的keras框架（Chollet，2015）。根据多次实验结果，设置建筑、边界、分割线的损失权重α1，α2，α3=1，0.4，0.4，焦点损失、骰子损失、Jaccard损失、边界损失的权重β1，β2，β3，β4=1/3，1/3，1/3，1。网络使用Nadam（Dozat，2016）优化器，初始学习率设置为5×10-4，在第10、50 个周期，学习率除以5。网络总共训练80 个周期，每个周期迭代250 次，共2万次迭代。网络的输入尺寸为196×196，批尺寸为16。后处理中，融合处理参数k= 1，建筑标记阈值l0，l1，l2，l3= 0.55，0.65，0.75，0.85，阈值t0，t1，t2=0.45，0.5，0.85。实验数据已在2.1 节介绍。本文采用四折交叉验证比较实验效果，即每一次交叉验证使用27 张尺寸5000×5000 的图像为训练集，剩余9张为测试集，分别训练和评估，并且所有实验使用相同的随机数。

3.1 指标

根据经验（Maggiori 等，2017），本文使用召回率Recall，准确率Precision，F1- score 和IOU 等指标衡量网络性能，计算公式如下：

式中，TP表示正确检测（True Positive）；FP表示误检（False Positive）；FN表示漏检（False Negative）。召回率反映真实样本遗漏比率，准确率反映预测结果误检比率。F1- score 是精确度和召回率的综合测量值。IOU是预测结果与真实标签的重叠比率。

本文使用新指标评估后处理结果。对建筑标签中每一个建筑，计算其与预测结果中每一个建筑的IOU，并将最大值作为该建筑实例的IOU，若IOU＞0.5，视为该建筑实例被检测；计算其与预测结果中每一个建筑的重叠面积，占该建筑面积的比例，并将最大值记为该建筑实例的Cover 值。由此延伸实例召回率Recallc（Hosang 等，2016）、覆盖比例Pc及IOUunit指标：

式中，Rc为实例召回率，Pc为覆盖比例，gti表示真实标签的第i个建筑掩模，pdj表示预测结果的第j个建筑掩模，N为建筑总数，Ⅱ(*)为示性函数。根据IOUi与Coveri的定义，若真实标签的建筑A与预测结果的建筑B的IOU＞0.5，则建筑A至少50%像素与建筑B重叠，即覆盖比例是实例召回率的上界。

3.2 数据预处理与边界损失

本节探究数据预处理、边界损失与网络性能的关系，比较建筑类网络，无边界损失的网络，本文网络和基准网络的建筑预测图。由Chatterjee和Poullis（2019），基准网络使用原始的ICT-Net，网络最后一层仅预测建筑类别，使用交叉熵作为损失函数；建筑类网络，与基准网络结构相同，使用Lurban作为损失函数；未应用边界损失的网络，与本文网络结构相同，修改边界损失权重β4为0。

表1显示比较结果，与基准、建筑类和无边界损失的网络相比，本文网络在IOU、召回率和F1-score 表现较优，其中本文网络相比基准网络提高IOU 约1%，对于IOUunit和实例召回率，本文网络高于前三者至少2.8%和5.7%，表明应用了数据预处理和边界损失的本文网络能更准确地检测建筑实例。

表1 不同方案下网络性能的指标比较Table 1 The results of different networks

本节选取两个子区域进一步比较，如图4，相较于其他3种网络，本文网络提取的建筑边界形态更接近真实标签。

图4 不同方案的预测结果比较Fig.4 The results of different networks

3.3 后处理

本文选取一典型示例，阐明后处理的具体步骤。

如图5，将图5（a）遥感影像输入网络，预测3 类概率，其中对图5（c）网络预测的建筑类概率，以0.5 二值化形成图5（d）建筑预测图。同时合成3类概率形成图5（e）Gp图，以标记阈值l0形成图5（f）预测图P0。以标记阈值l1形成图5（g）建筑标记L1。采用基于标记的分水岭分割算法以Gp图和L1分割P0，形成建筑实例，每个实例以回归树预测IOU。将每一个IOU 大于接受阈值的实例标记为建筑，形成图5（h）后处理预测图，其中不同颜色表示不同建筑实例。

图5 后处理示例Fig.5 An example of postprocess

为验证后处理性能，比较本文网络、合成处理和后处理的结果：建筑预测图、预测图P0和后处理预测图，并与基准网络的建筑预测图比较。对比结果如表2，后处理与本文网络相比，提高实例召回率10.5%，提高IOUunit4.5%，与基准网络相比，提高实例召回率22.9%，提高IOUunit10.7%，表明后处理在本文网络的基础上进一步提高建筑检测精度。

表2 本文方法不同阶段的结果比较Table 2 The results at different stages

本文以两个典型例子展示本文方法各阶段的结果，进一步说明后处理的作用，如图6，红色框展示结果改善的区域。图6（e）中蓝色建筑显示后处理从“粘连”的掩模中区分建筑实例，其检测结果与参考标签基本吻合。绿色建筑处于阴影或树荫下，未能被本文网络和后处理检测。

图6 本文方法不同阶段的结果比较Fig.6 The results at different stages

4 讨 论

4.1 网络结构

本节讨论网络结构的性能差异，选用U-net和FC-DenseNet56 （Jégou 等，2017），如2.2 节对ICT-Net 修改，改变网络最后一层结构，并应用3类标签和边界损失。

表3， 显示ICT-Net相较于U-net 和FCDenseNet56，在准确率、F1- score 和IOU 上较优，图7呈现训练过程中的损失和指标，间接证明了上述结论。

表3 网络结构比较Table 3 Comparison of different networks

图7 训练过程的损失图和IOU指标图Fig.7 The loss and IOU in training phase

4.2 参数选取

由2.3 节所述，合成公式参数k和建筑标记阈值l0，l1，l2，l3对后处理性能有较大影响，本节讨论参数的选取。

针对建筑标记阈值threshold，固定合成公式参数k= 1，选取建筑标记阈值threshold ∈{0.05，0.1，…，0.95 共19 个值，生成预测图P0，计算实例召回率和覆盖比例。如图8（a）所示，随着建筑标记阈值的提高，覆盖比例逐渐减小，实例召回率逐渐增大。当建筑标记阈值较小时，实例召回率较低，后处理无法充分利用网络输出的概率信息，且覆盖比例与召回率的差距过大，易带来过多噪音和预测负担；当建筑标记阈值过大时，后处理倾向于对预测图过度分割，覆盖比例衰减较大，覆盖比例作为实例召回率的上界限制了后处理的提升。权衡实例召回率和覆盖比例，本文选取四个建筑标记阈值l0，l1，l2，l3=0.55，0.65，0.75，0.85，以l0= 0.55 对应的建筑标记形成预测图P0，以l1，l2，l3= 0.65，0.75，0.85 形成的建筑标记对预测图P0进一步分割。

针对参数k，固定建筑标记阈值threshold =0.55，选取k∈{0，0.1，…，1}共11个值，生成预测图P0，计算实例召回率和覆盖比例。如图8（b），实例召回率随k逐渐增大，而覆盖比例下降较小，本文选取k=1。

图8 参数选取Fig.8 Parameters selection

5 结 论

本文针对密集建筑区域的建筑检测，提出了一种结合神经网络自下而上分组的分水岭算法。该方法通过数据预处理，利用原建筑标签生成两类辅助掩模标签；在性能较优的ICT-Net 基础上，修改网络结构和损失函数，针对两类辅助标签添加边界损失；然后在网络预测结果的基础上，通过一个结合分水岭算法和梯度提升回归树的后处理完成建筑实例的准确检测。实验表明，在中小建筑密集的城市区域，本文检测方案能有效解决密集建筑掩模“粘连”的问题，准确定位建筑边界，检测建筑实例，检测效果明显优于原始的ICT-Net网络。本文方法的优势如下：

（1）数据预处理生成的两类辅助标签，注重边界特征，可丰富边界信息。

（2）经修改的网络结构、添加的边界损失函数充分考虑建筑物特点和边界优化的需要，相较于基准网络，本文网络提取的建筑边界形态更接近真实标签。

（3）后处理借助分水岭分割算法可利用网络预测的概率信息进一步优化建筑边界，采用梯度提升回归树可筛选最优的分割结果，因此后处理在本文网络的基础上可使建筑实例的检测精度更高。

本文方法目前仅引入了建筑边界信息，未考虑建筑的内部结构，今后将在这一方面深入探索。