龚 龑，舒 宁，王 琰，陶 醉，尹淑玲

1.武汉大学 遥感信息工程学院，湖北 武汉430079；2.中国科学院 遥感应用研究所，北京100101；3.武昌理工学院 信息工程学院，湖北 武汉430223

1 引 言

从遥感影像解译的角度来看，像斑不是无序的待识别模式集合，而是在空间呈现特定分布规律的解译基元集合。文献［1］的地理学第一定律认为，任何事物在空间相关，此种影响随距离的加大而递减。从获取服务于人类应用的有效解译信息的目的出发，符合人类认知规则的地物类别概念有着显著的空间关联特征。当前的研究实践表明，单纯依靠像斑自身的固有遥感数据特征，不足以对像斑的类别信息进行合理判断［2］。在面向对象的影像分析中，在光谱特征分析的基础上，综合考虑像斑的空间关系特征，对于获取符合认知语义规则的各类专题解译结果具有重要意义，这一点已引起国内外研究人员的重视［3－5］。

基于马尔可夫随机场（Markov random field，MRF）的分析方法从数学角度提供一种空间关系分析途径，使得空间关系特点可以用数值进行描述。文献［6］证明吉布斯分布率与之等价，进而为空间关系分析提供一种概率分布表达形式。最初的MRF模型是针对二维空间中标量状态进行分析的，在医学灰度影像分析中已经得到诸多实践［7］。近年来针对遥感影像的MRF模型研究也逐渐深入［8－11］。国内外近十几年的研究表明，MRF模型对二维影像空间中的栅格邻域关系有着良好的描述能力，这种特点对于分析遥感影像中像元空间分布规律较为适合。

但是，在面向对象的遥感影像分析中，像斑作为特征的承载基元，其空间特征的非规则性与经典的MRF模型有较大差异。在经典MRF的定义中，点位呈规则分布，相邻点位之间具备等间隔特征。在连续条件下，表现为坐标空间的平直特征，离散条件下则表现为栅格数据的行列规则性。建立在上述基础上的MRF模型可以充分借助概率分布函数等工具予以描述。与此不同，地物目标的复杂性决定了像斑形状大小各异，像斑及其分布通常不具备规则特征［12－14］。虽然在数学领域，文献［15］曾证明MRF可以定义在图上，以区域邻接图（region adjacency graph，RAG）为代表的像斑 MRF定义也初步展开［16－17］，但是相关的应用研究亟待深入。同时，对于MRF应用中重要的吉布斯参数问题［18］，在像斑数据中，可以描述某种像斑空间分布规律的数据并不如像素条件下充足，这势必造成样本数据相对待估计吉布斯参数数目的不足，加之像斑子团的非规则性，参数估计变得困难甚至不可行。

本文的非规则MRF模型与后向传播（back propagation，BP）神经元网络相结合，采取像斑类别共生规律统计的方法实现无参数形式势函数构建，避免吉布斯参数显式估计过程，在保留MRF优良特性的同时，实现针对遥感影像像斑的空间关系特征应用。

2 非规则MRF模型

设遥感影像Q定义在有限二维空间区域L上，经过影像分割［19］，在区域L上形成N个像斑Segi（1≤i≤N），即并且Segi∩Segj＝∅（i≠j）。设每个像斑有M个特征，构成一个M维随机矢量，任意像斑Segi的随机矢量记为Xi＝［Xi1Xi2…XiM］，则这N个随机矢量构成区域L上的矢量特征随机场X＝［X1X2…XN］，x＝［x1x2…xN］为X的一个实现。同时，像斑Segi的解译标号记为Si，所有N个像斑标号构成L上的标号随机场S＝［S1S2…SN］，记s＝［s1s2…sN］是随机场S的一个实现。像斑与临近像斑属性的相关性，决定了其解译标号存在马尔可夫特性。本文立足于对随机场X和S的分析，力图在状态空间ΩS中寻找出最符合空间关系规律的随机场S的实现形式。

设ηi是L上与像斑Segi有相互作用像斑组成的L的子集，即Segi的邻域，这种子集构成的子集类η＝｛ηi｜1≤i≤N，Segi∉ηi，ηi⊂L｝，称为L的邻域系统。本文以拓扑相邻及其传递定义像斑Segi的邻域ηi。如图1所示，像斑SegA的一阶邻域ηA（1）是与SegA拓扑相邻的像斑（含岛）集合，即ηA（1）＝｛SegB1，SegB2，SegB3，SegB4，SegB5｝。像斑SegA的二阶邻域记为ηA（2），表示在ηA（1）的基础上增加ηA（1）中像斑的拓扑相邻像斑，ηA（2）＝ηA（1）∪｛SegC1，SegC2｝，3阶及以上的邻域系统以此类推。

图1 非规则MRF模型邻域系统示例Fig.1 Example of neighborhood system in IR－MRF

基于上述定义，将遥感影像像斑空间分布规律视为全局分布特性，像斑及其邻域的解译标号视为局部特性。上述随机场，可以用如下吉布斯分布［9］予以描述

式中，Z为配分函数，针对随机场所有可能的实现进行归一化，可作常数处理；U（s）为吉布斯能量函数，像斑解译标号集s与空间分布规律的不一致性程度决定了能量函数U（s）的大小。从兼顾像斑自身特征和像斑空间关系的角度出发，本文的非规则马尔可夫随机场（irregular－Markov random field，IR－MRF）模型的能量函数形式如下

该能量函数等于单像斑能量函数U0（s）和邻域能量函数Uη（s）加权之和，p为权值。

3 吉布斯能量函数

3.1 单像斑势能

单像斑势能基于像斑自身特征，衡量其被赋予某一类别标号的合理性。本文以下式表示单像斑能量函数

即在给定像斑特征x和先验知识Θ的条件下，标号随机场以s形式实现所产生的能量。f0（s；x，Θ）表示衡量s、x和类别特征知识Θ相容性的测度函数，当s、x和Θ的相容性越高时，U0（s）取值越小，否则取值越大。本文采用增强全局寻优能力的改进BP神经元网络［20］来实现单像斑势能计算。BP神经元网络的输入层为单像斑的特征向量x，输出层节点个数为解译标记s的种类数K，每输出数值为像斑Segi相对与各种解译标记k的隶属概率Pik。引入上述BP神经元网络后，f0（s；x，Θ）为网络所输出的类别隶属概率。同时，BP神经元网络输出的类别隶属概率，将作为像斑势能计算中的动态定权参考。对于像斑Segi，在所有K个输出的类别隶属概率中，最大值记为Pmaxi，则式（2）中单像斑势能权重按下式计算

上式的含义为，若仅利用自身特征，像斑Segi相对于某一类别已具有明显的可分性优势，则单像斑势能权重较大，反之则较小。

3.2 无参数邻域势函数

邻域能量函数从像斑邻域状态角度，评估像斑标号空间分布的合理性。本文定义下式表示邻域能量函数

即在给定邻域空间分布知识Θ的条件下，标号随机场以s形式实现所产生的能量。s与空间分布知识Θ越相容，Uη（s）取值越小，邻域能量函数由邻域子团势函数VI（c）（s；Θ）累加而成。式（5）中，c表示子团，即中心像斑与其邻域系统内的像斑所构成的像斑组合。其中，对n阶邻域系统ηA（n），中心像斑和任意一个邻域像斑的组合称为n阶双像斑子团，中心像斑和任意两个邻域像斑的组合称为n阶三像斑子团。本文主要针对一阶双像斑子团进行讨论，以图1为例，中心像斑SegA所对应的邻域子团共有15个，包括5个双像斑子团C（2）和10个三像斑子团C（3），双像斑子团C（2）如图2所示。

图2 一阶双像斑子团示意图Fig.2 Pair－segments cliques of first－order neighborhood

对子团c中像斑，V（c）I（s；Θ）从标号分布规律角度，度量标号实现形式s和先验知识Θ的不一致程度。通常将吉布斯分布应用于像素分析的主要方法，是把不同的像素子团c和不同的吉布斯参数相对应起来，通过这些参数来反映像素在空间分布上的特征。针对像斑空间关系分析的特点，本文不再以确定性的吉布斯参数作为像斑空间分布规律的描述方式，而是直接建立子团状态统计结果与势能V（c）I（s；Θ）之间的联系，从不同类别像斑共生概率角度设计如下无参数吉布斯邻域势函数

式中，Cntc表示中心像斑所具有的子团个数；PΘ（c）（s）表示针对像斑标号随机场S的一种实现s，当前中心像斑及其子团c相应有一种实现形式s0和Rc（s），在具备空间分布规律Θ的像斑数据中，以s0为中心像斑的子团c以Rc（s）形式出现的概率即为PΘ（c）（s）。各类子团相应的PΘ（c）（s）基于对已标记像斑数据的统计得到。

3.3 像斑标号分析流程

基于上述非规则无参数随机场模型，进行遥感影像像斑空间分析主要包括如下环节：①通过影像分割获取像斑，结合具体应用背景确定解译标记s的集合，在此基础上，针对每个通道，分别统计每个像斑内所有像元的光谱均值等信息，建立特征随机场X；②基于样本进行BP神经元网络单像斑分类训练，获取单像斑相对于解译标记的隶属概率，以此为基础计算单像斑子团势能U0（s），并按照最大隶属概率所对应的类别对像斑赋予解译标号；③基于上述已标号的像斑数据，统计双像斑子团的空间共生概率P（c）Θ（s）；④依据像斑空间特征统计结果和网络输出层数值，依式（4）确定单像斑和多像斑的势能权重，将当前待分析像斑作为中心像斑，分别赋予不同的标号，按照式（6）计算像斑共生势能，然后按式（5）计算能量函数。最后，依据 式（1）、式（2）计算吉布斯概率并依据最大后验概率，对像斑标号进行更新；⑤如果像斑标号更新比例小于阈值或迭代次数超过阈值，则结束，否则转至第③步。

4 试验与分析

4.1 多光谱影像分析

4.1.1 试验数据

为分析单像斑特征和邻域特征在像斑分类中的影响，本文利用武汉地区2002年TM影像进行试验。首先利用基于多波段影像进行区域合并法影像分割，获取像斑数据，如图3所示。然后对每个波段进行像斑内灰度均值计算，形成每个像斑的7维光谱特征向量。针对区域实际情况，定义了江水、河水、湖水、农田、裸地、新建筑区、旧建筑区、树木8种类别。人工选取两组样本像斑数据分别作为训练样本和检测样本。

图3 TM影像和分割后的像斑数据（局部）Fig.3 TM image and segments data（local）

4.1.2 像斑邻域特征分析

针对上述像斑，进行基于BP神经元网络的初始分类和邻域特征分析。像斑邻域特征共生势能的计算立足于子团类别共生概率P（c）Θ（s），本文通过对全体像斑标号数据进行统计获取该值。表1是对像斑初始分类标号的第一轮统计结果，该表在行方向具备概率的归一化特点。

分析上表可以发现，除了同类像斑相邻概率普遍较大外，在武汉地区，江水的邻域中出现农田、旧建筑区概率较大，分别为21.8%和16.2%，这反映该地区的城市格局特点，即旧城区沿江而建的历史原因，以及江边滩涂地区农业利用的现状。湖水邻域较多出现农田（21.5%）与树木（15.6%），一定程度上由远城区围湖造田的空间格局决定，并且与主城区湖泊普遍存在沿湖景观绿化树木有关。农田邻域出现裸地（12.7%）以及裸地邻域出现农田（40.3%）的情况比较普遍，这一方面与农作物生长的周期间歇特点有关，另一方面还与2000—2002年建设期间，近郊田间存在较多的待开发裸露地块有关。树木的邻域像斑中，旧建筑区出现概率（25.54%）远大于新建筑区的出现概率（3.56%），说明该时期旧城区绿化状况相对较好。新建筑区邻域为旧建筑区（42.3%）和农田（41.6%）的概率较大，反映新建筑区的两个主要产生位置，即老城区改造和城乡结合部的城市扩张。以上分析说明子团类别共生概率较真实地反映了像斑空间邻域特征。

表1 基于双像斑子团初始标号的共生概率Tab.1 Co－occurrences of initial pair－segments clique labels

4.1.3 整体试验效果分析

基于上述统计参数计算共生势能，结合BP神经元网络获取的单点势能，可获取不同标号的后验概率，依据最大后验概率准则对像斑标号进行更新。本文采用两种势能组合形式进行分类，图4是分类效果图。

其中图4（a）仅考虑单像斑势能U0，图4（b）考虑单像斑能量U0和邻域共生能量Uη（s）。图5显示了局部细节区域的分类效果对比。

经由实地调查数据及人工判读检验，对比图5（a）、图5（b）可以发现，江水与湖水的混淆情况在U0＋Uη（s）方法中得到改善，这反映了仅依靠单像斑自身光谱特征尚不足以将湖水和江水有效区分，但是在考虑湖水与江水的空间关系后，在湖水中被误判为江水的像斑明显减少。图5（c）、图5（d）反映旧建筑区中的地块误分为农田的情况在考虑像班邻域特征能量后得到显著改善。

为定量评价上述两类方法的分类精度，本文基于检测样本进行总体分类精度评估，针对像斑特点设计如下像斑分类精度指标

式中，K为类别数；Nk表示第K类的测试像斑数，当测试像斑i分类正确时，Pki为1，否则为0；CntOfPtki表示第K类的第i个测试图斑包含的像元数；Sk为第K类全部测试图斑总像元数。基于单像斑势能分析的分类结果精度0.886，加入邻域势能分析后，分类的精度为0.915，整体精度提高2.9%。

图4 像斑分类结果Fig.4 Segments classification results

图5 分类效果细节比较Fig.5 Classification results detail comparison

4.1.4 类别空间依赖性分析

为分析不同地物类别对空间关系的依赖程度，分别对引入空间势能后各类精度ρk进行了分析

式中，等号右端中各变量含义与式（7）相同。

表2列出了针对不同地物，分别采用U0和U0＋Uη（s）分析的精度以及精度变化量Δρk。

从表2中可以看出，对于在光谱特征上相似但在语义概念上存在差异的湖水和河水，在考虑邻域关系势能后分类精度分别提高17.2%和1.6%，而同为水体的江水由于含沙量引起的光谱特征差异，在单像斑分析时已能达到94.7%的较好精度，所以引入空间关系后精度改善有限，为0.6%。农田和新旧筑区的分类精度分别提高4.0%、3.2%和3.3%，说明这几类地物也存在明显的空间相关性。

值得注意的是，对于树木和裸地，在考虑邻域势能后精度反而下降了4.1%和2.5%，这说明对这两类地物而言，空间关系依赖性并非解译的优势信息，反而对结果产生负面影响。虽然本文的模型中已经采用了基于式（4）的单像斑势能动态定权方法，但对于这类地物其邻域势能权重最好置零。

4.2 高分辨率影像分析

4.2.1 试验数据

为检验IR－MRF模型在高分辨率影像中的应用效果，本文以房屋信息提取为应用进行了如下试验。试验数据为天津某地区2005年World－View多波段影像，采用0.5m全色波段影像进行影像分割，获取像斑数据，如图6所示。再基于2m多光谱波段数据进行像斑内灰度均值计算，形成每个像斑的蓝、绿、红和近红外4维光谱特征向量。针对数据特点和房屋提取的应用背景，定义了楼房顶、墙面、平房顶、阴影、草坪、水泥地、水体、路面、树木、车辆10种类别标号。

图6 WorldView影像和分割后的像斑数据（局部）Fig.6 WorldView image and segments data（local）

表2 各类像斑分类精度对比Tab 2 Classify accuracy comparison of different class

4.2.2 像斑标记与房屋提取

分别以仅考虑单像斑势能和兼顾邻域势能两种方式，对所有像斑进行分类。在此基础上，针对房屋提取，选取与之相关的楼房顶、墙面、平房顶和阴影等4种标号类型，即得到候选房屋目标提取结果，经人工干预后将属于这4类标号的像斑标记为房屋，输出结果。表3反映了依据人工判读结果进行检验的房屋像斑提取精度分析，引入邻域势能后漏检率和误检率分别降低了5.6%和16.6%，与单像斑分析相比，尤其在误检率上明显改善。图7（a）和（b）分别为本文IR－MRF模型U0＋Uη（s）分析的像斑标记结果和房屋目标结果。试验证明在房屋提取中通过IR－MRF模型考虑像斑空间关系将有效提高分析精度。

表3 房屋像斑提取精度分析Tab.3 Building segments extraction accuracy

图7 IR－MRF模型像斑标记结果和房屋提取结果Fig.7 Segments labeled results and building extraction results

5 结 论

空间关系特征是遥感影像分析时的重要信息，本文的非规则MRF模型以全新的方式定义了像斑状态下随机场，采用新的无参数吉布斯邻域势能函数，在BP神经元网络的辅助下，实现了像斑空间关系特征的定量化描述，并且使MRF的应用摆脱了复杂的吉布斯参数估计的制约，是一种适用于遥感影像像斑空间关系分析的实用模型。

试验分析表明，这种无参数的非规则随机场，能够有效地表达不同地类之间存在的空间关系，与仅考虑单像斑的分析方法相比，这种兼顾单像斑和邻域像斑的分析模型，非常适合房屋、建筑等存在显著空间关系依赖特点的地类或目标提取，以及光谱相似但人类语义概念存在差异的目标之间的区分。对于绝大多数地类目标来讲，这种空间关系的应用能够提高整体分类精度的。本文的研究也表明，对于少数空间关系特点存在随机性的地类（如树木、裸地等），将空间关系作为自动约束有可能会降低分类精度，对于这类目标，在分析时应当依据先验知识增加单像斑势能权重。

在本文的非规则无参数MRF模型中，势函数定义具备开放性特点，除了像斑类别相邻规律以外，针对双像斑及三像斑子团，可进一步设计有效的邻接、方位及距离等空间关系势能函数，丰富空间关系分析层次。

［1］ TOBLER W.A Computer Movie Simulating Urban Growth in the Detroit Region ［J］.Economic Geography，1970，46：234－240.

［2］ JOHN A.Analysis of Remotely Sensed Data：The Formative Decades and the Future［J］.IEEE Transaction on Geoscience and Remote Sensing，2005，43（3）：422－432.

［3］ BLASCHKE T.Object Based Image Analysis for Remote Sensing ［J］.ISPRS Journal of Photogrammetry and Remote Sensing，2010，65（1）：2－16.

［4］ TAN Qulin.Urban Building Extraction from VHR Multispectral Images Using Object－based Classification［J］.Acta Geodaetica et Cartographica Sinica，2010，39（6）：618－623.（谭衢霖.高分辨率多光谱影像城区建筑物提取研究［J］.测绘学报，2010，39（6）：618－623.）

［5］ JIANG Chongya，LI Manchun，LIU Yongxue.Full－automatic Method for Coastal Water Information Extraction from Remote Sensing Image［J］.Acta Geodaetica et Cartographica Sinica，2011，40（3）：332－337.（江冲亚，李满春，刘永学.海岸带水体遥感信息全自动提取方法［J］.测绘学报，2011，40（3）：332－337.）

［6］ GEMAN S，GEMAN D.Stochastic Relaxation，Gibbs Distribution，and Bayesian Restoration of Images［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1984（6）：721－741.

［7］ SCHERRER B，FORBES F，GARBAY C，et al.Distributed Local MRF Models for Tissue and Structure Brain Segmentation［J］.IEEE Transactions on Medical Imaging，2009，28（8）：1278－1295.

［8］ YANG Honglei，PENG Junhuan.Remote Sensing Classification Based on Markov Random Field and Fuzzy c－means Clustering［J］.Acta Geodaetica et Cartographica Sinica，2012，41（2）：213－218.（杨红磊，彭军还.基于马尔可夫随机场与模糊C均值的遥感影像分类［J］.测绘学报，2012，41（2）：213－218.）

［9］ LIU Guoying，QIN Qianqing，MEI Tiancan，et al.Supervised Image Segmentation Based on Tree－Structured MRF Model in Wavelet Domain［J］.IEEE Geoscience and Remote Sensing Letters，2009，6（4）：850－854.

［10］ SERPICO B，MOSER G.Weight Parameter Optimization by the Ho－Kashyap Algorithm in MRF Models for Supervised Image Classification［J］.IEEE Transactions on Geoscience and Remote Sensing，2006，44（12）：3695－3705.

［11］ KOMODAKIS N，PATAGIOS N，TZIRITAS G.MRF Energy Minimization and Beyond via Dual Decomposition［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2011，33（3）：531－552.

［12］ LIU Yu，GUO Qinghua，KELLY M.A Framework of Region－based Spatial Relations for Non－overlapping Features and Its Application in Object Based Image Analysis［J］.ISPRS Journal of Photogrammetry and Remote Sensing，2008（63）：461－475.

［13］ GALLEGUILLOS C，RABINOVICH A，BELONGIE S.Object Categorization Using Co－occurrence，Location and Appearance［C］∥Proceedings of IEEE Conference on Computer Vision and Pattern Recognition.California：［s.n.］，2008：1－8.

［14］ XU Zhu，MENG Yanzi，LI Zhilin，et al.Recognition of Structures of Typical Road Junctions Based on Directed Attributed Relational Graph［J］.Acta Geodaetica et Cartographica Sinica，2011，40（1）：125－131.（徐柱，蒙艳姿，李志林，等.基于有向属性关系图的典型道路交叉口结构识别方法［J］.测绘学报，2011，40（1）：125－131.）

［15］ ZHENG Zhaobao.Image Analysis of Markov Random Field Method［M］.Wuhan：Wuhan University Press，2000.（郑肇葆.图像分析的马尔可夫随机场方法［M］.武汉：武汉大学出版社，2000.）

［16］ NIU X.A Semi－automatic Framework for Highway Extraction and Vehicle Detection Based on a Geometric Deformable Model［J］.Journal of Photogrammetry and Remote Sensing，2006，58（3－4）：166－186.

［17］ FLORENCE T，MICHEL R.Markov Random Field on Region Adjacency Graph for the Fusion of SAR and Optical Data in Radargrammetric Applications ［J］.IEEE Transactionson Geoscience and Remote Sensing，2005，43（8）：1920－1928.

［18］ ALY A F，REFAAT M，AYMAN E.A Unified Framework for MAP Estimation in Remote Sensing Image Segmentation［J］.IEEE Transaction on Geoscience and Remote Sensing，2005，43（7）：1617－1634.

［19］ WANG Aiming，SHEN Lansun.Research on Image Segmentation［J］.Measurement and control technology，2000，19（5）：1－6.（王爱民，沈兰荪.图像分割研究综述［J］.测控技术，2000，19（5）：1－6.）

［20］ YIN Shuling，SHU Ning，LIU Xinhua.Classification of Remote Sensing Image Based on Adaptive Genetic Algorithm and Improved BP Algorithm［J］.Geomatics and Information Science of Wuhan University，2007，32（3）：201－204.（尹淑玲，舒宁，刘新华.基于自适应遗传算法和改进BP算法的遥感影像分类［J］.武汉大学学报：信息科学版，2007，32（3）：201－204.）