图像分割是图像处理中的重点问题之一，它包括两个基本任务：确定图像类别数及分割同质区域。目前，灰度图像的分割算法相对成熟，主要分为阈值分割( 李积英 等，2014; 尹奎英 等，2011; 吴一全 等，2012)、聚类分割( 卢洁 等，2011; 田小林 等，2008; 于波 等，2013)和统计分割( 倪维平 等，2011; 杨学志 等，2014; 郑玮 等，2008)等。随着传感器技术的不断发展，人们可以获取越来越多的彩色遥感图像，由此，其分割技术的开发引起广泛的重视。对彩色遥感图像而言，受地物目标及其分布的多样性、复杂性和未知性等因素的限制，人为确定其类别数非常困难。因此，实现彩色遥感图像的可变类分割是相当必要的，同时也是十分困难的工作。

研究表明，针对回归变量选择、目标识别、图像分割等诸多可变类问题，基于统计的处理方法是最为有效的方法(Green，1995; Richardson和Green，1997; Kato，2008)。为了解决上述可变类问题，Green (1995)提出RJMCMC(Reversible Jump Markov Chain Monte Carlo)算法。该算法通过构建马尔可夫随机场链，使其状态能够在不同维度的参数空间跳变，最终实现目标概率密度函数的模拟。在Green(1995)工作基础上，Richardson和Green(1997)利用RJMCMC算法模拟贝叶斯混合模型，实现混合模型中组分数及各组分概率分布参数估计。但该研究只适用于小数据量混合模型。Kato(2008)将高斯混合模型与RJMCMC算法相结合，以完成可变类彩色图像分割。该方法虽能实现彩色图像分割，但图像分割细节模糊。

以上方法均为基于像素的图像分割，但基于像素的方法很难解决高分辨率遥感图像分割问题，分割精度较低，为此，研究者试图在图像建模过程中引入区域处理的思想，即，对图像域进行几何划分，并建立基于区域的图像模型。常用的几何划分方法有：规则划分(Askari 等，2013)、Voronoi划分(Lucarini 等，2009; 赵泉华 等，2013)和泊松划分(Schneider，2009，2010)等。由于规则划分在划分图像域时原理简单、操作方便，可快速地实现图像域的划分。

本文提出了一种融合规则划分和RJMCMC算法的彩色遥感图像可变类分割方法。首先，利用规则划分将图像域划分成若干个相同大小的规则子块，并建立基于区域的图像分割模型。然后，利用RJMCMC算法自动确定图像类别数，并实现粗分割。为了提高图像精度，提取分割图像同质区域的边界线，并以这些边界线为中心生成缓冲区，再利用邻域准则对缓冲区的每个像素重新分类，生成新的分割区域边界线。在RJMCMC算法中，设计的移动操作类型包括：分裂或合并实类、改变参数矢量、改变标号及生成或删除空类。为了验证提出的方法，分别对合成及彩色遥感图像进行可变类分割实验，定性和定量的评价结果表明提出方法的可行性及有效性。

彩色遥感图像 z ={ z _s ， s ∈ S }可以看作定义在图像域 S 上的离散随机场 Z ={ Z _s ， s ∈ S }的一个实现，其中， Z _s ={ Z _sd ， d=1， $ \cdots $ ， r}为光谱测度矢量， Z _sd 为 d 波段光谱分量， r 是图像波段的总数， s 为像素位置(在不至混淆情况下，亦作为像素索引)。

为了建模遥感图像，首先利用规则划分将 S 划分成 J 个含 m × m 个像素规则子块， S ={ P _j ， j=1， $ \cdots $ ， J}，其中， J= n / m ²， n 为图像总像素数。所有统计模型将建立在规则划分的图像域上。 图 1(a)为规则划分的图像域，每个规则子块含3×3个像素； 图 1(b)为 P _j 与它的邻域子块集合 NP _j 。

为每个规则子块 P _j 分配一个标号变量 X _j ∈{1， $ \cdots $ ， k}，表示其类属性，其中 k 为图像的总类别数，设其为随机变量，并假设服从某种先验分布。每个目标类的几何区域由一组规则子块拟合而成。显然，对所有规则子块的标号集合形成一个随机标号场 X ={ X _j ， j=1， $ \cdots $ ， J}，该随机场的每个实现对应于遥感图像的一种分割。

为了表达邻域子块标号的相关性，在规则划分的图像域上建立MRF模型，并采用改进的静态Potts模型定义其邻域关系( Besag，1986; Strauss，1975)。假设 X _j ( j=1， $\cdots$ ， J)相互独立，则 X 的概率密度函数可表示如下

$\begin{aligned} & p({\textbf{\emph{X}}}|k)=\prod\limits_{j=1}^{J}{p({{X}_{j}}|{{X}_{{{j}^{'}}}},{{P}_{j'}}\in {{N}}{{{P}}_{j}})}= \\ & \text{ }\prod\limits_{j=1}^{J}{\frac{\exp \left\{ -\eta \sum\limits_{{{P}_{j'}}\in {{N}}{{{P}}_{j}}}{t({{X}_{j}},{{X}_{j'}})} \right\}}{\sum\limits_{{{l}^{'}}=1}^{k}{\exp \left\{ -\eta \sum\limits_{{{P}_{j'}}\in {{N}}{{{P}}_{j}}}{t({{l}^{'}},{{X}_{j'}})} \right\}}}} \\ \end{aligned}$

(1)

式中， η 为邻域子块的空间作用参数； NP _j 为规则子块 P _j 的八邻域子块的集合(如 图 1(b))， t 为指示函数，若 X _j ≠ X _j' ， t ( X _j ， X _j' )=1；否则为0。

假设，同一个规则子块 P _j 中具有标号 X _j = l 的随机变量 Z _j ={ Z _s ， s ∈ P _j }服从独立的多值Gaussian分布，其概率密度函数具体如下

式中， Φ _l =( μ _l ， Σ _l )， μ _l =( u _{l 1}， $\cdots$ ， u _lr )为其分布的均值矢量、 Σ _l =[ Σ _lpq ] _{r x r} ( p 、 q ∈{1， $\cdots$ ， r})为协方差矩阵。

进一步假设，所有规则子块相互独立，则 Z ={ Z _j ， j=1，…， J}的概率密度函数为

$ \begin{aligned} p({\textbf{\emph{Z}}}|{\textbf{\emph{X}}},{{{\varPhi }}},k)=\prod\limits_{j=1}^{J}{p({{{\textbf{\emph{Z}}}}_{j}}|{{X}_{j}}=l,{{{{\varPhi}}}_{l}})} \end{aligned} $

(3)

假设类别数 k 为随机变量，服从均值为 λ 的泊松分布，则其概率密度函数为

$ \begin{aligned} p(k)=\frac{{{\lambda}^{k}}}{k!}\exp (-\lambda) \end{aligned} $

(4)

多值Gaussian分布的均值(协方差)参数 μ _l ( Σ _l )的各分量，服从均值 μ _μ ( μ _Σ )、标准差 σ _μ ( σ _Σ )的正态分布且相互独立( Yong 等，2006)，则其概率密度函数可分别写为

$\eqalign{ & p(\mu ) = \prod\limits_{l = 1}^k {\prod\limits_{i \in \{ 1, \cdots ,r\} } {p({\mu _{li}})} } \cr & = \prod\limits_{l = 1}^k {{\prod _{i \in \{ \{ 1, \cdots ,r\} }}{1 \over {\sqrt {2\pi } {\sigma _\mu }}}\exp \left\{ { - {{{{({\mu _{li}} - {\mu _\mu })}^2}} \over {2\sigma _\mu ^2}}} \right\}} \cr} $

(5)

$\eqalign{ & p(\sum ) = \prod\limits_{l = 1}^k {\prod\limits_{i \in \{ 1, \cdots ,r\} } {\prod\limits_{i' \in \{ 1, \cdots ,r\} } {p({\sum _{lii'}})} } } \cr & = \prod\limits_{l = 1}^k {\prod\limits_{i \in \{ 1, \cdots ,r\} } {\prod\limits_{i' \in \{ 1, \cdots ,r\} } {{1 \over {\sqrt {2\pi } {\sigma _\sum }}}\exp \left\{ { - {{{{({\sum _{lii'}} - {\mu _\sum })}^2}} \over {2\sigma _\sum ^2}}} \right\}} } } \cr} $

(6)

多值Gaussian分布的参数矢量 Φ =( μ ， Σ )，假设 μ 和 Σ 相互独立，且根据贝叶斯定理得到联合后验概率密度函数 p ( k， X ， Φ | Z )，即

$\begin{aligned} & p(k,{\textbf{\emph{X}}},{{{\varPhi }}}|{\textbf{\emph{Z}}})\propto p({\textbf{\emph{Z}}}|{{{\varPhi}}}{\textbf{\emph{X}}},{\textbf{\emph{X}}},k)p({\textbf{\emph{X}}}|k)p({{{\varPhi}}}|k)p(k) \\ &=p({\textbf{\emph{Z}}}|{{{\varPhi}}},{\textbf{\emph{X}}},k)p({\textbf{\emph{X}}}|k)p(k)p({{{\mu}}}|k)p({\textbf{\emph{∑}}}|k) \\ \end{aligned}$

(7)

假设总参数矢量 ${{\varphi}}$ =( Φ ， X ， k)且相互独立，利用RJMCMC算法对 ${{\varphi}}$ 进行相关采样。在相关采样过程中状态 ${{\varphi}}$ 到状态 ${{\varphi}}$ ^*，其状态的转换概率为 q ( ${{\varphi}}$ ， ${{\varphi}}$ ^*)接受率为( Green，1995)，

$a(\varphi ,{{\varphi }^{*}})=\min \left\{ 1,\frac{p({{\varphi }^{*}}|Z)q({{\varphi }^{*}},\varphi )}{p(\varphi |Z)q(\varphi ,{{\varphi }^{*}})} \right\}$

(8)

当参数空间维度不变时，式(8)可简化为( Green，1995)

$a({{\varphi }_{l}},\varphi _{l}^{*})=\min \left\{ 1,\frac{p(Z|\varphi _{l}^{*})p(\varphi _{l}^{*})}{p(Z|{{\varphi }_{l}})p({{\varphi }_{l}})} \right\}$

(9)

式中， ${{\varphi}}$ _l 为总参数矢量 ${{\varphi}}$ 的第 l 个元素， ${{\varphi}}_l^*$ 为 ${{\varphi}}$ _l 候选参数。

当参数空间维数改变时，其接受率表示为

$a(\varphi ,{{\varphi }^{*}})=\min \left\{ 1,\text{ }\frac{p({{\varphi }^{*}}|Z){{q}_{m}}({{\varphi }^{*}},\varphi )}{p(\varphi |Z){{q}_{m}}(\varphi ,{{\varphi }^{*}})} \right\}$

(10)

式中， m 移动操作的转换概率 q _m ( ${{\varphi}}$ ， ${{\varphi}}$ ^*)可在改变参数空间维数时保持细致平衡。若接受率 a 大于在[0，1]间生成的随机数，那么状态由 ${{\varphi}}$ 变到 ${{\varphi}}$ ^*；否则，状态保持不变。

在RJMCMC算法中，设计了如下的移动操作。对每次迭代需遍历所有的移动操作；在每次迭代中，每个移动操作的初始状态为上一个操作的结果状态；而前一次迭代参数结果是后一次迭代的初始参数。

(1) 分裂或合并实类

在图像域中有对应的规则子块隶属该类别，则该类别为实类， k _r 表示其实类别数。分裂操作就是将一组同标号的规则子块集合分解为两组标号不同的规则子块集合。在当前状态下，假设 S ={ S ₁， $\cdots$ ， S _l ， $\cdots$ ， S _kr }、 Z ={ Z ₁， $\cdots$ ， Z _l ， $\cdots$ ， Z _kr }和 X ={ X ₁， $\cdots$ ， X _j ， $\cdots$ ， X _J }，其中 S _l ={ P _j ， X _j = l}， Z _l ={ Z _j ， X _j = l}。首先，在实类中随机抽取任意实类 l，其概率为1/ k _r 。将 S _l 分解为标号分别为 l 和 k+1的规则子块集合，即 ${{{S}}_l} = {{S}}_l^* \cup {{S}}_{k + 1}^* $，其中 S _l ^*={ P _j ， X _j ^*= l}， S _k+1 ={ P _j ， X _j ^*= k+1}。将 Z _l 分解为 Z _l ^*={ Z _j ， X _j ^*= l}和 ${{Z}}_{k + 1}^* = \left\{ {{{{Z}}_j}, \ X_j^* = k + 1} \right\}$。其对应的新参数矢量分别为 Φ _l ^*=( μ _l ^*， Σ _l ^*)=( μ _l ， Σ _l )、 ${{\varPhi}} _{k + 1}^* = \left\{ {{{\mu}} _{k + 1}^*, \ {{\varSigma}} _{k + 1}^*} \right\}$，其中 μ ^* _k+1 、 ${{{\varSigma}} _{k + 1}^*}$服从各自的先验分布。分裂操作的接受率为

$ \begin{aligned} {{a}_{{{s}_{r}}}}({\varphi} ,{{\varphi}^{*}})=\min \left\{ 1,\text{ }{{R}_{{{s}_{r}}}} \right\} \end{aligned} $

(11)

式中，

$$

(12)

$\begin{align} & {{R}_{{{s}_{r}}}}=\frac{p(Z_{l}^{*}|\Phi _{l}^{*},X_{l}^{*},k+1)p(Z_{k+1}^{*}|{{\Phi }_{k+1}},{{X}_{k+1}},k+1)}{p({{Z}_{l}}|{{\Phi }_{l}},{{X}_{l}},k)}\times \\ & \frac{p(X_{l}^{*}|k+1)p(X_{k+1}^{*}|k+1)}{p({{X}_{l}}|k)}\times \frac{{{r}_{{{m}_{k+1}}}}}{{{r}_{{{f}_{k}}}}}\times \frac{p(k+1)}{p(k)}\times \\ & \frac{p({{\mathbf{\mu }}_{k+1}})p({{\sum }_{k+1}})}{q(u)}\times \left| \frac{\partial {{\varphi }_{k+1}}}{\partial ({{\varphi }_{k}},u)} \right| \\ \end{align}$

(12)

式中， r _{f k} = f _k /{ k _r (2 ^B -2)}，其中选择分裂操作的概率为 f _k ， B 为分裂前标号为 l 的所有规则子块的总个数； ${r_m}_{_{k + 1}}$ =2 m _{k+ 1}/{ k _r ( k _r +1)}，其中 m _k 为选择该操作的概率。

合并实类操作是分裂实类操作的对偶操作，即在现有实类中随机选择两个实类，将其合并成一个新的实类。因此，其接受率为

$ \begin{aligned} {{a}_{{{m}_{r}}}}({\varphi},{{\varphi}^{*}})=\min \left\{ 1,\text{ }1/{{R}_{{{s}_{r}}}} \right\} \end{aligned} $

(13)

(2) 改变参数矢量 Φ

等概率(1/ k)的从{1， $\cdots$ ， k}中任抽一标号 l，对应的规则子块集合 S _l ={ P _j ， X _j = l}；改变 μ _l ( Σ _l )通过改变其中的某些因素完成，候选参数 μ _l ^*( Σ _l ^*)从均值和方差分别为 μ _l ( Σ _l )、 ε _μ ( ε _Σ )的正态分布抽取，其中 ε _μ ( ε _Σ )为预定义参数；算法顺序改变 μ 和 Σ 。

首先判断 l 是空类还是实类，若是空类，则对应的接受率为

${{a}_{{{\varPhi }}}}({{{{\varPhi}}}_{l}},{{{\varPhi}}}_{l}^{*})=\min \left\{ 1,\text{ }\frac{p({{{\varPhi}}}_{l}^{*}|k)}{p({{{{\varPhi}}}_{l}}|k)} \right\}$

(14)

若 l 是实类，则其接受率为

${{a}_{\Phi }}({{\Phi }_{l}},\Phi _{l}^{*})=\min \left\{ 1,\frac{\prod\limits_{{{P}_{j}}\in {{S}_{l}}}{p({{Z}_{j}}|{{X}_{j}}=l,\Phi _{l}^{*})p(\Phi _{l}^{*}|k)}}{\prod\limits_{{{P}_{j}}\in {{S}_{l}}}{p({{Z}_{j}}|{{X}_{j}}=l,{{\Phi }_{l}})p({{\Phi }_{l}}|k)}} \right\}$

(15)

(3) 改变标号场

等概率(1/ J)从图像域 S ={ P _j ， j=1， $\cdots$ ， J}中以任抽一规则子块 P _j ，其对应标号 X _j = l ；然后从现有实类内等概率抽取候选标号 l ^*且 l ^*≠ l 。则改变 l 的接受率可计算为

${{a}_{X}}({{X}_{j}}=l,{{X}_{j}}={{l}^{*}})=\min \left\{ 1,\frac{p({{Z}_{j}}|{{X}_{j}}=l,{{\Phi }_{{{l}^{*}}}})p({{X}_{j}}={{l}^{*}}|k)}{p({{Z}_{j}}|{{X}_{j}}=l,{{\Phi }_{l}})p({{X}_{j}}=l|k)} \right\}$

(16)

(4) 生成或删除空类

如果一类别数在图像域中没有对应的规则子块集合，那么这个类别数为空类， k _e 表示其空类数。移动操作过程为生成一空类，设其标号为 k+1；由于该操作仅改变参数矢量 Φ 和标号场 X ，故令候选总参数矢量 ${{\varphi}}$ ^*=( Φ ^*， X ^*， k+1)且 Φ ^*、 X ^*、 k+1相互独立，其中 Φ ^*=( Φ ₁， $\cdots$ ， ${{\varPhi}}_l^*$ ， $\cdots$ ， ${{\varPhi}}_{k+1}^*$ )，假设 ${{\varPhi}}_l^*$ =( ${ {\mu}} _l^*$ ， ${{\Sigma}}_l^*$ )=( μ _l ， Σ _l )、 ${{\varPhi}}_{k+1}^*$ =( ${{\mu}}_{k+1}^*$ ， ${{\Sigma}}_{k+1}^*$ )， ${{\mu}}_{k+1}^*$ 、 ${{\Sigma}}_l^*$ 服从各自的先验分布。由式(8)可得到生成空类操作的接受率为

${{a}_{{{b}_{e}}}}({\varphi},{{\varphi}^{*}})=\min \{1,{{R}_{{{b}_{e}}}}\}$

(17)

式中，

$\begin{align} & {{R}_{{{b}_{e}}}}=\frac{p({{\varphi }^{*}}){{r}_{{{d}_{k+1}}}}({{\varphi }^{*}})}{p(\varphi ){{r}_{{{b}_{k}}}}(\varphi )q(u)}\left| \frac{\partial {{\Phi }^{*}}}{\partial (\Phi ,u)} \right|= \\ & \frac{p({{\Phi }^{*}}|k+1)p({{X}^{*}}|k+1)p(k+1)}{p(\Phi |k)p(X|k)p(k)}\times \frac{{{r}_{{{d}_{k+1}}}}({{\varphi }^{*}})}{{{r}_{{{b}_{k}}}}(\varphi )q(u)}\left| \frac{\partial {{\Phi }^{*}}}{\partial (\Phi ,u)} \right| \\ \end{align}$

(18)

式中， r _{d k+1} = d _{k +1}/( k _e +1)和 r _bk = b _k ，其中 d _{k +1}、 b _k 分别为选择其对应操作的概率。设 u =( ${{\mu}}_{k+1}^*$ ， ${{\Sigma}}_{k+1}^*$ )，与参数矢量 Φ =( μ ， Σ )=( μ ₁， $\cdots$ ， μ _k ， Σ ₁， $\cdots$ ， Σ _k )相互独立，所以， Φ ^*=( Φ ， u )，因此， p ( Φ ^*)= p ( Φ ) q ( u )，式(18)中Jacobian项为1( Green，1995)。

综上所述，式(18)可简化为

$\begin{align} & {{R}_{{{b}_{e}}}}=\frac{{{d}_{k+1}}p({{X}^{*}}|k+1)p(k+1)}{{{b}_{k}}({{k}_{0}}+1)p(X|k)p(k)}= \\ & \frac{{{d}_{k+1}}\underset{{{P}_{j}}\in S}{\mathop{\prod }}\,\sum\limits_{l'=1}^{k}{\exp \left\{ -\eta \sum\limits_{{j}'\in N{{P}_{j}}}{t(l',{{X}_{{{j}'}}})} \right\}}\lambda }{{{b}_{k}}({{k}_{e}}+1)\underset{{{P}_{j}}\in S}{\mathop{\prod }}\,\sum\limits_{l'=1}^{k+1}{\left\{ -\eta \sum\limits_{{j}'\in N{{P}_{j}}}{t(l',{{X}_{{{j}'}}})} \right\}\left( k+1 \right)}} \\ \end{align}$

(19)

删除空类操作是生成空类操作的对偶操作，即在现有的空类中随机选择一个空类，将其删除。因此其接受率为

$ \begin{aligned} {{a}_{{{d}_{e}}}}({\varphi},{{\varphi}^{*}})=\min \{1,\text{}1/{{R}_{{{b}_{e}}}}\} \end{aligned} $

(20)

图 2(b)为由Worldview-2上截取不同地物填充在 图 2(a)(模板)的、包含红、绿和蓝3个波段的合成遥感图像，其类别数为5，尺度为128×128像素。

为了合理划分图像域，采用规则划分，分别将图像域划分成若干个4×4和8×8的规则子块，并在此基础上，在Intel Core2 2.53Ghz/1G内存/Matlab7.0环境下，验证提出的方法，然后通过比较其结果，选择合适的规则子块尺寸。在实验中，进行了20000次迭代， 图 3为截取前300次迭代对应的图像实类别数 k _r 的变化。 图 3(a) (b)分别为在迭代过程中采用4×4和8×8子块规则划分对应的 k _r 变化。

实验表明：(1)在 k _r 稳定前，其值在[2， k _max] 之间跳变；通过 图 3(a) (b)可以看出，对应4×4和8×8子块的规则划分，分别在70和25次迭代后 k _r 达到稳态，即， k _r =5；(2)在跳变期间，没有任何一个类别数保持连续20次不变；因此，为稳健起见，假设在迭代过程中，如果一类别数连续50次不变，则认为 k _r 达到稳态，并设该数为图像实类别数；(3)8×8子块规则划分对应 k _r 的收敛速度要比4×4子块规则划分快，因为划分的规则子块尺寸越大，同一幅图像划分的规则子块数越少，其对应类别数收敛的也就越快。

图 4为分别采用4×4和8×8规则子块划分的粗分割结果，其中 图 4 (a)分别为 k _r 初收敛在5时的初分割结果； 图 4 (b)分别为在20000次迭代的粗分割结果，用它们的彩色均值表示其同质区域； 图 4 (c)分别为粗分割轮廓线与原图的叠加。

为了提高图像分割精度，设计了精细化操作。首先，以提取出的粗分割轮廓线为中心，生成宽为2个规则子块的缓冲区；然后采用邻域准则，即利用静态Potts模型定义像素间的邻域关系( Besag，1986)，对缓冲区内每个像素重新分类，形成新的边界。 图 5(a)分别为4×4和8×8子块规则划分的粗分割轮廓线(白色)与缓冲区轮廓线(黑色)的叠加； 图 5(b)分别为其细化结果，各分割区域用该区域像素彩色均值表示； 图 5(c)分别为对应的细分割轮廓线与原图叠加。从视觉效果来看，4×4子块规则划分提取的细分割轮廓线与实际轮廓线更加吻合。从 图 3、 图 4和图 5可以看出，提出的方法不仅能自动确定图像类别数，而且可以很好的实现图像分割。

为了对提出方法进行定量评价，将合成彩色遥感图像的模板( 图 2(a))作为标准分割数据，求其混淆矩阵，并据此分别计算用户精度、产品精度、总精度及Kappa值，见 表 1。从 表 1可以看出，4×4子块规则划分对应的产品精度、用户精度及总精度均大于等于99%，而Kappa值高达0.996；而8×8子块规则划分对应的第3个用户精度低达94.8%，而其Kappa值为0.981。综上所述，采用4×4子块规则划分图像域更为适合。

图 6为3幅分辨率为1.8 m的、包含Worldview-2红、绿和蓝3个波段的真彩色遥感图像，尺度128×128像素，人为判读其类别数分别为3，3，4。在Intel Core2 2.53Ghz/1G内存/Matlab7.0环境下，利用提出方法对 图 6进行分割实验， 图 7为图像的分割结果，其中 图 7(a)为其对应的粗分割结果，而 图 7(b)为精细化后的分割结果。通过 图 7可以看出，提出的方法不仅能自动确定图像类别数，还能实现图像分割。

为了验证提出方法的分割精度，分别对其分割结果进行了视觉评价，图 8为3幅真彩色遥感图像对应的粗轮廓线(白色)、细轮廓线(黑色)与原图叠加结果。通过图 8可以看出，3幅图像的粗轮廓线与实际边界的吻合不甚理想，而细分割轮廓线能很好的与实际边界吻合，因此，设计的精细化操作有效地提高了图像分割精度，且提出方法的分割精度很高。

图 9为4幅4波段IKONOS的整景图像，分辨率为10 m，尺度为256×256像素，人为判别其类别数为3、3、5和6，其中图 9(a)包含红、绿和蓝3个波段的真彩色图像，图 9(b)为近红外波段的灰度图像。在Intel Core2 2.53Ghz/1G内存/Matlab7.0环境下，利用提出方法对其进行可变类分割，图 10(a)为其对应的粗分割结果，而图 10(b)为精细化结果。为了验证提出方法的分割精度，分别对提取分割结果的粗轮廓线(白色)和细轮廓线(黑色)，并将其轮廓线叠加到原图中，见图 10(c)。通过图 10可以看出，提出方法对多波段图像的可行性及有效性。

为了对提出方法进行定量评价，以 图 6和 图 9中的遥感图像的手动标记模板 图 11为标为标准分割数据，分别求其对应最终分割结果的混淆矩阵，并依次计算其对应的产品精度、用户精度、总精度及kappa值(见 表 2)。由 表 2计算得，平均产品精度为92.9%，平均用户精度为91.1%，平均总精度为94.2%，平均Kappa值为0.911。通过这些参数值，验证了提出方法的可行性及有效性。

利用ENVI软件中的ISODATA算法对Worldview-2合成及彩色遥感图像和多光谱IKONOS图像进行可变类分割， 图 12为对应的分割结果。通过视觉评价可以看出，ISODATA算法在某些情况下虽然能正确确定图像类别数，但分割结果精度不高，甚至有的分割结果与实际情况不符。而提出方法不仅能正确确定图像类别数，且能很好地实现图像分割。因此，相对而言，本文方法更具有优越性。

利用 赵泉华等人(2013)提出的基于Voronoi划分和EM/MPM的分割算法分别对Worldview-2合成及彩色遥感图像和多光谱IKONOS图像进行不变类分割， 图 13为对应的分割结果。通过视觉评价可以看出，该方法虽然有效地提高区域分割精度，但由于某些Voronoi多边形跨越边界，导致区域边界分割精度降低。而提出方法不仅能自动确定类别数，还可以很好地实现分割。另外，从时间上角度看，文献算法实现分割的平均运行速度为45 min，而提出方法实现分割的平均运行速度为40 min，因此提出方法更有效率。

提出了一种基于规则划分和RJMCMC算法的可变类遥感图像分割方法，该方法不仅能自动确定图像类别数，还可以精确地实现同质区域分割。为了提高图像分割精度，利用规则划分将图像域划分成规则子块，并在此基础上建立统计模型，但是分割结果边界与实际边界的拟合不甚理想，为此，设计了精细化操作。利用提出方法、ISODATA算法和基于Voronoi划分和EM/MPM的分割算法分别对Worldview-2合成及彩色遥感图像和多光谱IKONOS图像进行分割实验，通过提出方法的分割结果，验证了该方法的可行性及有效性；通过对比实验，验证了提出算法的优越性。在今后的工作中，将采用其他几何划分方法划分图像域。