CN102254187A - 一种基于元胞自动机的遥感影像非监督分类方法 - Google Patents

Abstract

本发明公开了一种基于元胞自动机的遥感影像非监督分类方法，本发明方法在传统元胞自动机的基础上，采用地图代数的距离变换和模糊系统计算邻居元胞对中心元胞的隶属度，并基于影像非监督分类的基本原理，对遥感影像进行非监督分类。与传统遥感影像非监督分类方法相比，本发明方法在分类精度上有一定程度的改善；与基于模糊c均值的遥感影像分类方法相比，本发明方法对分辨率较大的遥感影像进行分类具有效率上的优势。

Description

一种基于元胞自动机的遥感影像非监督分类方法

技术领域

本发明属于遥感影像分类方法，尤其涉及一种基于云胞自动机的遥感影像非监督分类方法。

背景技术

非监督分类算法正是为解决难以获取先验样本这一问题而产生的，主要有K-均值(K-Means)、迭代自组织数据分析技术方法(ISODATA)和模糊c均值三种方法。K-均值是以各类中元素到其类中心距离平方和最小为聚类准则，并根据这一准则将待分元素分到各类的一种较典型的逐点修改迭代的动态聚类算法[1-4]。迭代自组织数据分析技术方法是K-均值的改进版本，它在迭代过程中能动态修改类别数^[5，6]，而K-均值法的类别数从始至终是固定的。这两种算法都将待分元素绝对的分到某一类中，而模糊c均值^[7-10]采用模糊集代替了这种非“黑”即“白”的绝对划分，多采用样本点与聚类中心间的距离并用隶属度的平方加权表示，把类内误差平方和目标函数变成带加权的误差平方和目标函数，可在一定程度上改善K-均值和迭代自组织数据分析技术方法的精度度^[11，12]。

地物之间是相互影响的，而基于上述三种方法的遥感影像分类方法仅根据待分元素与各类中心的距离来判断其归属问题，并没有考虑这种地物之间的相互作用。元胞自动机是一时间、空间和状态都离散的动力系统，它根据元胞及其邻居的当前状态确定该元胞下时刻的状态^[13]，能够考虑相邻地物之间的影响。自冯诺依曼于1947-1955年的早期研究和后来的在图像处理和模式识别等等方面的研究工作^[14]，元胞自动机可用在研究包括通信、信息传递、计算、构造、生长、复制、竞争与进化等一般现象^[15，16]，广泛应用于社会、经济、军事和科学研究的各个领域，元胞自动机在图像分类方面的改进，多是利用其进行边缘检测或对分类后的结果图进行去噪声处理、增强等^{[17，18，19]}，鲜少有人利用二维元胞自动机根据非监督分类的原理直接进行多波段遥感影像的分类。

而且根据地理学第一定律，地物之间存在着越相邻的地物相似程度越大。因而距离中心元胞距离不同的邻居，应该对其有不同的影响。但现今常用的二维元胞邻居(冯诺依曼型、摩尔型和扩展的摩尔型)虽然定义简单，容易实现，却会随着扩散步数的增加存在各向异性的缺陷，而且邻居元胞对中心元胞的影响也存在着“影响”和“不影响”、非0即1这种以绝对的衡量标准处理问题的不足。

本文涉及的参考文献主要如下：

[1]SHOKRI Z.SELM and M.A.ISMAIL. K-Means-Type Algorithms：AGeneralized Convergence Theorem and Characterization of Local Optimality.IEEETRANSACTIONS ON PATTFRN ANALYSIS AND MACHINE INTELLIGENCE，，1984，6(1)：81-87.

[2]HUANG ZX.Extensions to the k-Means Algorithm for ClusteringLarge DataSets with Categorical Values.Data Mining and Knowledge Discovery.1998(2)：283-304.

[3]Anderberg，M.R.Cluster  Analysis for Applications.USA：AcademicPress.1973.

[4]D.Hall and G. Ball.“Isodata：A novel method of data analysis andpatternclassification，”Stanford Res.Inst.，Stanford，CA，Tech.Rep.，1965

[5]J.MacQueen.“Some  methods for  classification and analysis ofmultivariateobservations，”in Proc.5th Berkeley Symp.Mathematics，StatisticsandProbability，1967，pp.281-296.

[6]Duda，R.O.，Hart，P.E.Pattern Classification and SceneAnalysis.WileyInterscience，New York.1973.

[7]Dunn，J.C.Some recent investigations of a new fuzzypartition algorithm andits application to pattern classificationproblems.J.Cybernetics.1974，4：1-15.

[8]Bezdek，J.C.Pattern Reeognition with Fuzzy ObjectiveFunction Algorithms.New York：Plenum Press.1981.

[9]Keh-Shih Chuang.Hong-Long  Tzeng.Sharon Chen.Jay Wu.Tzong-JerChen.Fuzzy c-means clustering with spatial information for imagesegmentation.Computerized Medical Imaging and Graphics.2006，30：9-15.

[10]Xizhao  Wang.Yadong  Wang.LijuanWang.Improving  fuzzy  c-meansclustering basedon feature-weight learning.Pattern  Reeognition Letters.2004，25：1123-1132.

[11]YanfeiZhong，Liangpei Zhang，Bo Huang，and Pingxiang  Li，AnUnsupervised Artificial Immune Classifier forMulti/Hyperspectral Remote SensingImagery.IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING.2006，44(2)：420-431.

[12]Jiao LC，GongMG，Wang S，et al.Natural and Remote Sensing ImageSegmentation Using Memetic Computing.IEEE COMPUTATIONAL INTELLIGENCEMAGAZINE.2010，5(2)：78-91.

[13]周成虎，孙战利，谢一春.地理元胞自动机研究.科学出版社，2001.

[14]Kendall Preston，M.J.B.Duff.Modern cellular automata：theory andapplications.Springer Press，1984.

[15]Smith A.Cellular Automata Theory.Technical Report 2，Stanford ElectronicLab.Stanford University，1969.

[16]Perrier，J.Y.，Sipper，A.G. and Zahnd，J.，Toward a Viable，Self-reproducingUniversal Computer Physica D，1996，97：335-352.

[17]王海军，张文婷，贺三维.利用元胞自动机和模糊c均值进行图像分割.武汉大学学报：信息科学版，2010，35(11)：1288-1291.

[18]陈炎雄，郝志峰，温尚清，杨晓伟.基于二维细胞自动机和中值运算的图像滤波方法.计算机辅助工程，2006，3：54-58.

[19]李伟，杨绍清，刘松涛.元胞自动机在图像高密度椒盐噪声去噪中的应.激光与红外，2010，5：563-567.

[20]胡鹏，游涟，杨传勇，吴艳兰.地图代数.武汉大学出版社，2006年。

发明内容

针对现有技术的不足，本发明在传统元胞自动机邻居度量方法的基础上、采用地图代数的距离变换[20]和模糊隶属度衡量周围邻居对中心元胞的影响程度，并在此基础上提出了一种基于元胞自动机的遥感影像非监督分类方法。

为解决上述技术问题，本发明采用如下技术方案：

一种基于元胞自动机的遥感影像非监督分类方法，包括以下步骤：

S1、在栅格空间下计算元胞自动机邻居元胞对中心元胞的隶属度，该步骤包括如下子步骤：

S1-1、确定待分类遥感影像中栅格空间的大小，并利用地图代数的距离变换计算栅格空间中各邻居元胞对中心元胞的距离，其中，栅格空间的一个栅格对应遥感影像的一个像素，即为一个元胞；中心元胞为栅格空间中的每一像素；

S1-2、根据地物之间的相互作用，通过引进模糊系统，计算栅格空间中邻居元胞对中心元胞的隶属度；

S2、对待分类遥感影像进行非监督分类，该步骤进一步包括以下子步骤：

S2-1、读取待分类遥感影像的灰度值到矩阵img(M×N)中，根据待分类遥感影像的特征确定其聚类数目n，并随机确定各聚类中心的初始值

且 $\min \left(\mathrm{img}(i,j)\right)\≤{\mathrm{Center}}_k^0\≤\max \left(\mathrm{img}(i,j)\right);$ 其中：img(M×N)即为待分类遥感影像的元胞空间，M、N分别为待分类遥感影像的长、宽像素大小，n为自然数且n∈[1，+∞]，k取1、2、...n，i取1、2、...M，j取1、2、...N，min(img(i，j))和max(img(i，j))分别表示待分类遥感影像所有像素灰度值中的最小灰度和最大灰度；

S2-2、根据公式

计算t时刻待分类遥感影像中各中心元胞(i，j)距各类聚类中心的距离

如果中心元胞(i，j)在t时刻到第m类聚类中心的距离

小于其到其他聚类中心的距离，则该元胞(i，j)属于第m类聚类中心，即t时刻元胞(i，j)的状态值S^t(i，j)＝m；根据公式

计算t时刻待分类遥感影像中各中心元胞(i，j)的邻居元胞到各类聚类中心的距离

中心元胞(i，j)的邻居元胞在t时刻到第m类聚类中心的距离

小于到其他聚类中心的距离，则该邻居元胞属于第m类；其中：img(i，j)为各中心元胞(i，j)的灰度值，NB(i，j)为各中心元胞(i，j)的邻居元胞的灰度值，为t时刻第k类聚类中心的灰度值，i取1、2、...M，j取1、2、...N，k取1、2、...n，m为自然数且m∈[1，n]，||…||表示距离运算；

S2-3、根据步骤S1中得到的邻居元胞对中心元胞的隶属度，并根据转换规则确定各中心元胞在t+1时刻的状态，所述的转换规则如下：

$S^{t+1}(i,j)=\left\{\begin{array}{cc}k& \mathrm{Max}\_U_{{\mathrm{NB}}_{(i,j)}^k}>\max (\mathrm{Max}\_U_{\mathrm{NB}})\×{\mathrm{\ϵ}}_1^{\′}\\ S^t(i,j)& \mathrm{Max}\_U_{{\mathrm{NB}}_{(i,j)}^k}\≤\max (\mathrm{Max}\_U_{\mathrm{NB}})\×{\mathrm{\ϵ}}_1^{\′}\end{array}\right.;$

其中：S^t+1(i，j)和S^t(i，j)分别为t+1时刻和t时刻中心元胞(i，j)的状态；

为中心元胞(i，j)的圆形邻居元胞中属于第k类聚类中心的最大隶属度，i取1、2、...M，j取1、2、...N，k为自然数且k∈[1，n]；Max_U_NB为所有元胞的圆形邻居属于某一类聚类中心的最大隶属度和的集合；ε′₁为根据遥感影像质量和分类目所设定的阈值，其值介于0到1之间；

S2-4、根据公式 ${\mathrm{Center}}_k^{t+1}=\frac{{\mathrm{\Σ}}_{\mathrm{img}(i,j)\∈{\mathrm{img}}_k^{t+1}}\mathrm{img}(i,j)}{N_k^{t+1}}$ 求取t+1时刻待分类遥感影像的新的聚类中心

其中：

为t+1时刻第k类聚类中心的灰度值，k取1、2、...n；img(i，j)为中心元胞(i，j)的灰度值，i取1、2、...M，j取1、2、...N：

为t+1时刻属于第k类聚类中心的中心元胞灰度集，k取1、2、...n；

为在t+1时刻属于第k类聚类中心的元胞数。

S2-5、判断t+1时刻和t时刻待分类遥感影像中同一类聚类中心的接近程度，如果公式 $\underset{k=1}{\overset{n}{\mathrm{\Σ}}}{({\mathrm{Center}}_k^{t+1}-{\mathrm{Center}}_k^t)}^2\≤{\mathrm{\ϵ}}_2$ 成立，k取1、2、...n，或者当迭代次数超过预先设定的最大迭代次数时，则迭代结束，执行步骤S3，否则，重新执行步骤S2-1；其中：ε₂为设定的收敛阈值；

S3、根据步骤S2所得的待分类遥感影像的聚类中心结果，分别将各类聚类中心值取整后赋给不同类别的所有元胞，即得到所有元胞的状态矩阵，并输出状态矩阵，即得到已分类的结果影像。

上述步骤S1-2中栅格空间中邻居元胞对中心元胞的隶属度通过隶属度函数得到，其中，隶属度函数为高斯型隶属度函数。

上述步骤S2-2中，中心元胞(i，j)距各类聚类中心的距离和中心元胞(i，j)的邻居元胞距各类聚类中心的距离

均为欧式距离。

与现有技术相比，本发明具有以下优点和有益效果：

1、本发明根据非监督分类的基本原理，将用地图代数的距离变换和模糊系统改进的元胞自动机用于遥感影像分类，与传统遥感影像非监督分类方法相比，在分类精度上有一定程度的改善；

2、本发明充分利用了矩阵运算，与基于模糊c均值的遥感影像分类方法相比，本发明方法对分辨率较大的遥感影像进行分类具有效率上的优势。

附图说明

图1为冯诺依曼邻居型示意图；

图2为摩尔邻居型示意图；

图3为基于地图代数距离变换的3×3栅格空间示意图；

图4为基于地图代数的距离变换的栅格空间(d＝4和d＝9)；

图5为基于距离变换的2×2栅格空间中邻居元胞对中心元胞的高斯隶属度；

图6为武汉市部分区域三波段遥感影像图；

图7为采用本发明方法对图6所示遥感影像图进行分类的结果；

图8为采用基于摩尔邻居的元胞自动机遥感影像分类方法对图6所示遥感影像图进行分类的结果；

图9为采用基于K-均值的遥感影像分类方法对图6所示遥感影像图进行分类的结果；

图10为采用基于模糊c均值的遥感影像分类方法对图6所示遥感影像图进行分类的结果；

图11为应用1中所述四种分类方法的执行效率和迭代次数的对比分析图；

图12为广东省揭西县城区部分区域三波段遥感影像图；

图13为采用本发明方法对图12所示遥感影像图进行分类的结果；

图14为采用基于摩尔邻居的元胞自动机遥感影像分类方法对图12所示遥感影像图进行分类的结果；

图15为采用基于基于K-均值的遥感影像分类方法对图12所示遥感影像图进行分类的结果；

图16为采用基于模糊c均值的遥感影像分类方法对图12所示遥感影像图进行分类的结果。

具体实施方式

基于K-均值(K-Means)和模糊c均值方法的传统的遥感影像非监督分类方法，仅以单个地物为分类对象，没有考虑地物之间的相互作用；到目前为止，利用元胞自动机改进的遥感影像非监督分类方法，一般是利用一维或者二维元胞自动机对分类后的结果影像进行去噪声处理，以达到更高分类精度的目的，但是，这种改进方法并没有将非监督分类的基本原理与元胞自动机模型直接结合来构造非监督分类器，而是在通过传统分类方法对图像分类的基础上再利用元胞自动机对分类结果进行的改进，并且，现有的元胞邻居构型以等权的方式衡量周围邻居对中心元胞的影响的，忽视了“越相邻的地物相似程度越大”这一客观事实，也是造成分类精度不够高的原因之一。

本发明在传统元胞自动机邻居度量方法的基础上、采用地图代数的距离变换和模糊隶属度衡量周围邻居对中心元胞的影响程度，并在此基础上提出了一种基于元胞自动机的遥感影像非监督分类方法，该方法包括以下步骤：

S1、在栅格空间下计算元胞自动机邻居元胞对中心元胞的隶属度，传统的计算邻居元胞对中心元胞隶属度的方法，是在确定邻居构型的基础上，将各个邻居对中心元胞的距离看成是一样的，即以一种等权的方式计算周围元胞对中心元胞的贡献；而本发明方法在距离变换后的栅格空间中，计算不同邻居元胞对中心元胞的影响程度，具体步骤如下：

S1-1、确定待分类遥感影像中栅格空间的大小，并利用地图代数的距离变换计算栅格空间中各邻居元胞对中心元胞的距离，其中，栅格空间的一个栅格对应遥感影像的一个像素，即为一个元胞；中心元胞为栅格空间中的每一像素；

在本发明中，邻居元胞和中心元胞都是指遥感影像的像素，在用元胞自动机处理遥感图像时，每一个图像像元都会成为中心元胞，从而对全图像进行处理，而中心元胞周围符合特定邻居构型的元胞在此时就是邻居元胞。

传统的二维元胞自动机邻居构型主要有冯诺依曼、摩尔和扩展的摩尔邻居型，它们定义简单，容易操作，但存在着随着扩散步数的增加而导致元胞演变的空间格局趋于不平衡状态和各向异性的缺陷，如冯诺依曼邻居型(如图1所示)在正东南西北四个方向上扩展较快，而摩尔邻居型(如图2所示)在东北、东南、西南和西北四个方向上扩散速度快于其他方向。

本发明对元胞自动机邻居构型进行改善，主要是利用地图代数的距离变换计算栅格空间中邻居元胞对中心元胞的距离，从而根据已经确定的邻域半径d得到一个n×n的栅格空间，再根据需要保留该栅格空间中距中心元胞元距离相同的邻居作为改进后的元胞自动机的邻居，其中，n是一个栅格空间在一条边长方向上的格子数，当每个栅格空间的边长为1个单位标准栅格时，d＝n；当栅格空间的边长为r个单位标准栅格时，d＝n*r。如确定相对于当前像元的偏移量即邻域半径为3，可以得到如图3所示的基于距离变换的邻居元胞对中心元胞在3×3的栅格空间中的距离关系，而图中虚线范围内的像素即为改进的2×2邻居元胞构型，粗实线范围内的为改进的3×3邻居元胞构型。利用这种方法得到的邻居构型在定义时避免了标准邻居定义的方向不等权特点，可以让元胞自动机快速实现类似于几何圆形的邻居扩散方式，如图4所示。

S1-2、根据地物之间的相互作用，通过引进模糊系统，计算栅格空间中邻居元胞对中心元胞的隶属度；

在S1-1的基础上，根据地物之间的相互作用，通过引进模糊系统，来快速计算栅格空间中邻居元胞对中心元胞的隶属度，从而得到一种改进的二维元胞自动机构型。

对于“越相邻的地物相似程度越大”这种客观现实，本发明引进模糊系统，使得空间范围可以一定程度属于某些中心点，空间属于中心的程度由“0”与“1”之间的一个数值来描述。具体刻画不同距离的邻居元胞对中心元胞的影响程度，即隶属度，是由隶属度函数来实现的。

隶属度函数可以是任意形状的曲线，取什么形状的曲线取决于分类的地物特征(地物之间实际的相互影响情况)、目的和允许的时间消耗等实际情况。常用的隶属度函数包括高斯型隶属度函数、钟形隶属度函数、神经元的非线性作用(sigmoid)函数型隶属度函数、梯形隶属度函数和三角形隶属度函数。图5是基于地图代数的距离变换的高斯型隶属度函数在2×2栅格空间中各邻居对中心元胞的隶属度值，图5所示的高斯型隶属度函数为

式中：c、σ为高斯型隶属度函数的两个特征参数，此处取c＝0，σ＝1；x为邻居元胞到中心元胞的基于地图代数的距离变换的栅格空间距离。

S2、对待分类遥感影像进行非监督分类，该步骤进一步包括以下子步骤：

S2-1、读取待分类遥感影像的灰度值到矩阵img(M×N)中，根据待分类遥感影像的特征确定其聚类数目n，并随机确定各聚类中心的初始值

且 $\min \left(\mathrm{img}(i,j)\right)\≤{\mathrm{Center}}_k^0\≤\max \left(\mathrm{img}(i,j)\right),$ 其中，img(M×N)即为待分类遥感影像的元胞空间；M、N分别为待分类遥感影像的长、宽像素大小；n为自然数且n∈[1，+∞]；k取1、2、...n；i取1、2、...M；j取1、2、...N；min(img(i，j))和max(img(i，j))分别表示待分类遥感影像所有像素灰度值中的最小灰度和最大灰度；

S2-2、根据公式

计算t时刻待分类遥感影像中各中心元胞(i，j)距各聚类中心的距离，记为根据公式

计算t时刻待分类遥感影像中各中心元胞(i，j)的邻居元胞距各聚类中心的距离，记为

如果中心元胞(i，j)在t时刻到第m类聚类中心的距离

小于其到其他聚类中心的距离，则该元胞(i，j)属于第m类聚类中心，即t时刻元胞(i，j)的状态值S^t(i，j)＝m；中心元胞的邻居的状态的确定与中心元胞的确定方法相同，也即是如果该邻居到第m类聚类中心的距离

小于到其他聚类中心的距离，则该邻居元胞属于第m类。其中，img(i，j)为各中心元胞(i，j)的灰度值；NB(i，j)为各中心元胞(i，j)的邻居元胞的灰度值；

为t时刻第k类聚类中心的灰度值；i取1、2、...M；j取1、2、...N；k取1、2、...n；m为自然数且m∈[1，n]；||…||表示距离运算，在后面的两个应用中均使用欧式距离；

S2-3、根据步骤S1中得到的邻居元胞对中心元胞的隶属度，并根据转换规则确定待分类遥感影像中各中心元胞在下一时刻的状态，所述的转换规则如下：

$S^{t+1}(i,j)=\left\{\begin{array}{cc}k& \mathrm{Max}\_U_{{\mathrm{NB}}_{(i,j)}^k}>\max (\mathrm{Max}\_U_{\mathrm{NB}})\×{\mathrm{\ϵ}}_1^{\′}\\ S^t(i,j)& \mathrm{Max}\_U_{{\mathrm{NB}}_{(i,j)}^k}\≤\max (\mathrm{Max}\_U_{\mathrm{NB}})\×{\mathrm{\ϵ}}_1^{\′}\end{array}\right.;$

其中，S^t+1(i，j)和S^t(i，j)分别为t+1时刻和t时刻中心元胞(i，j)的状态；

为中心元胞(i，j)的圆形邻居元胞中属于第k类聚类中心的最大隶属度，i取1、2、...M，j取1、2、...N，k为自然数且k∈[1，n]；Max_U_NB为所有元胞的圆形邻居属于某一类聚类中心的最大隶属度和的集合；ε′₁为根据遥感影像质量和分类目的所设定的阈值，其值介于0到1之间，当ε′₁为1时表明下一时刻中心元胞的状态完全取决于上一时刻该中心元胞的状态，而ε′₁为0时表明该中心元胞在下一时刻的状态完全由其邻居元胞所决定；

S2-4、根据公式 ${\mathrm{Center}}_k^{t+1}=\frac{{\mathrm{\Σ}}_{\mathrm{img}(i,j)\∈{\mathrm{img}}_k^{t+1}}\mathrm{img}(i,j)}{N_k^{t+1}}$ 求取t+1时刻待分类遥感影像的新的聚类中心

其中，

为t+1时刻第k类聚类中心的灰度值，k取1、2、...n；img(i，j)为中心元胞(i，j)的灰度值，i取1、2、...M，j取1、2、...N；

为t+1时刻属于第k类聚类中心的中心元胞灰度集，k取1、2、...n；

为在t+1时刻属于第k类聚类中心的元胞数；

S2-5、当情况1或情况2存在时，迭代结束，执行步骤S3，否则，重新执行步骤S2-1；

情况1为：判断t+1时刻和t时刻待分类遥感影像中同一类聚类中心的接近程度，如果公式 $\underset{k=1}{\overset{n}{\mathrm{\Σ}}}{({\mathrm{Center}}_k^{t+1}-{\mathrm{Center}}_k^t)}^2\≤{\mathrm{\ϵ}}_2$ 成立，k取1、2、...n，其中，ε2为根据影像大小和分类速度条件等要求设定的收敛阈值，会影响聚类的精度和速度，若ε₂数值设定过大，会导致聚类结果不稳定，ε₂数值设定过小，可能会导致过度计算，浪费时间；

情况2为：迭代次数大于等于预先设定的最大迭代次数时，其中，最大迭代次数由实验者根据实验结果的精度和执行时间要求而定，最大迭代次数是用来限制迭代次数的，有的时候如果一直达不到

成立这个条件，但是如果继续返回步骤S2-1会造成大量时间的浪费但又不能带来精度明显改善时，可以用最大迭代次数来使程序跳出循环，这个值的作用根据精度和时间要求而定)该值不能太小，否则将使程序在没有接近最优解时就执行完毕。

S3、根据步骤S2所得的待分类遥感影像的聚类中心结果，分别将各类聚类中心值取整后赋给不同类别的所有元胞，即得到所有元胞的状态矩阵，并输出状态矩阵，即得到已分类的结果影像。

下面将进一步说明本发明的应用及其有益效果。

应用1

如图6所示为湖北省武汉市城区某部分2005年分辨率为28.5米的多光谱美国陆地卫星(Landsat)的增强型转用绘图仪(ETM+)影像图，像素为403×299，波段为432三波段。

通过观察可知，图6所示影像主要分为四类：植被、水、建筑物和裸地，其中，植被类包括树木、草地以及农作物等；水类包含河流、湖泊以及没有栽种水生农作物的水田等；裸地类指道路广场等没有建筑物和植被覆盖的裸露地面。首先，在图6所示的原始影像上随机均匀采样，并通过目视解译和实地调查获取这些样本的真实类归属，实测数据结果的类归属统计值见表1所示。

表1各类样本的随机抽取情况表

分别采用本发明方法、基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法对图6所示的原始遥感影像进行分类。为了使结果具有可比性，上述四种方法均选取相同的迭代次数和收敛阈值，且本发明方法和基于摩尔邻居的元胞自动机遥感影像分类方法还选取了相同的转换阈值。基于摩尔邻居的元胞自动机非监督分类方法中采用的邻居结构为标准结构，其邻居元胞对中心元胞的隶属度均为1。上述四类方法的主要参数最大迭代次数、转换阈值ε₁、收敛阈值ε₂详见表2。

表2四种分类方法的参数选取情况

图7-10分别为采用本发明方法、基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法对图6所示的原始遥感影像进行分类的结果。从视觉角度直观对比分析图7-10，上述四种方法均能很好的区分出水类，但基于K-均值和模糊c均值的遥感影像分类方法对建筑物类、裸地类、植被类交错分布复杂的区域的分类结果很不理想，建筑物类用地容易被误分类为植被类、水类和裸地类。而采用本发明方法和基于摩尔邻居的元胞自动机遥感影像分类方法所得到的分类结果类似，这两种分类方法也存在建筑物类、裸地类、植被类被误分的情况，但在建筑物类的分类精度上明显优于基于K-均值和模糊c均值的遥感影像分类方法。

为了更有效的比较不同分类方法的分类结果，按照表1的抽样方式比较由上述四种分类方法获得的分类结果(如图7-10)和实测数据结果，并利用总精度和kappa系数定量评价了四种分类方法的分类准确度，如表3、4所示。表3所示为采样上述四种分类方法所得到的遥感影像分类结果的分类精度的误差矩阵，是目前反映分类结果好坏的一种表示方法，各列中粗体字的值越大则表明分类结果越好。表4中可以看出，本发明方法的分类效果明显优于其他三种方法：本发明方法的总精度最高，为83.57％，分别比基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法高3.49％、4.52％和7.80％；本发明方法的kappa系数也是最高，为0.7738，分别在基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法基础上改进了0.0455、0.0606和0.1012。从表3和表4中可以得出结论，本发明方法的分类效果明显优于其他三种方法。

表3四种分类方法的结果比较

表4四种分类方法的精度比较

传统的分类方法是以单个像素为对象进行分类，未考虑到相邻地物之间的相互影响；而基于摩尔邻居的元胞自动机遥感影像分类方法是基于标准的摩尔邻居，该邻居定义方式在东北、东南、西北、西南四个方向扩散较快而且将周围距中心元胞距离不同的八个邻居对中心元胞的影响等同，也即是采用相同的隶属度来衡量周围邻居对中心元胞的影响程度；本发明方法则同时考虑了地物之间的相互影响和周围元胞对中心元胞的随距离变化而产生的不同影响程度，可根据地物间的实际的相互作用情况确定邻居元胞对中心元胞的合适隶属度函数，从而保证分类结果更加精确。

在matlab中，当遥感影像图的分辨率较大时若采用遍历的方法实现元胞自动机的效率一般都较低，但矩阵运算效率很高。基于此，本应用中利用矩阵来实现本发明方法和基于摩尔邻居的元胞自动机遥感影像分类方法，运算速度比较理想。图11显示了在CPU为英特尔Core 2Duo E75002.93GHz、内存为4G，Window XP SP2配置环境下，不同分辨率下，分别采用上述四种方法处理图6所示原始遥感影像的运行时间和迭代次数。从运行效率来看，由于基于K-均值的遥感影像分类方法每次迭代过程简单且迭代次数较少，因而其执行时间最短；本发明方法在每一步迭代过程中由于邻居数量比基于摩尔邻居的元胞自动机遥感影像分类方法多、且引入的模糊系统造成计算复杂，故执行时间较长，但从整体来讲本发明方法的效率也是可以接受的，因为利用基于矩阵的运算方法计算一般大小的三波段遥感影像(例如分辨率为40dpi(360×285))仅需要半分钟的时间。当影像尺寸较小时，采用基于模糊c均值的遥感影像分类方法比基于元胞自动机的分类方法(本发明中，基于元胞自动机的分类方法指本发明方法和基于摩尔邻居的元胞自动机非监督分类方法)快；随着影像尺寸的增大，基于模糊c均值的遥感影像分类方法耗时的增速明显快于基于元胞自动机的分类方法，即，随着遥感影像分辨率的增加，基于元胞自动机的分类方法的时间的增速较基于模糊c均值的遥感影像分类方法慢，因而可以得出结论：基于元胞自动机的分类方法对分辨率较大的遥感影像图分类在时间方面具有一定的优势。但是，由于本发明方法的邻居元胞较多，当影像图较大时，对内存要求比较高。

应用2

图12所示为广东省揭西县某部分2008年的多光谱美国陆地卫星(Landsat)的专题绘图仪(TM)遥感影像图，像素为378×299，分辨率为10米。通过观察可知，图12所示影像主要分为三类：植被、建筑物和裸地。首先，在图12所示的原始影像上随即均匀采样，取样结果见表5。本应有中各主要参数的选择与应用1所选参数相同，但转换阈值有所不同，具体值见表6。

表5各类样本的随机抽取情况表

表6四种分类方法的参数选取情况

图13-16分别为采用本发明方法、基于摩尔邻居的元胞自动机遥感影像分类方法(CAS)、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法对图12所示的原始遥感影像进行分类的结果。从视觉角度直观对比分析图13-16，上述四种方法均能较好的区分出裸地类，与应用1相比，上述四种方法均容易将建筑物类用地误分类为其他类型用地，但本发明方法的分类结果还是由于其他三种分类方法。

利用总精度和kappa系数定量评价了采用四种分类方法对图12遥感影像进行分类的准确度，如表7所示。表7的数据可以充分说明了本发明方法提高了遥感影像的分类精度：本发明方法总精度和kappa值均最高，其总精度分别较基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法提高了0.47％、1.75％和3.7％，其kappa值分别在基于摩尔邻居的元胞自动机遥感影像分类方法、基于K-均值的遥感影像分类方法、基于模糊c均值的遥感影像分类方法的基础上改善了0.0064、0.0333和0.0574。

表7四种分类方法的精度比较

Claims (3)

1.一种基于元胞自动机的遥感影像非监督分类方法，其特征在于，包括以下步骤：

S1、在栅格空间下计算元胞自动机邻居元胞对中心元胞的隶属度，该步骤包括如下子步骤：

S1-1、确定待分类遥感影像中栅格空间的大小，并利用地图代数的距离变换计算栅格空间中各邻居元胞对中心元胞的距离，其中，栅格空间的一个栅格对应遥感影像的一个像素，即为一个元胞；中心元胞为栅格空间中的每一像素；

S1-2、根据地物之间的相互作用，通过引进模糊系统，计算栅格空间中邻居元胞对中心元胞的隶属度；

S2、对待分类遥感影像进行非监督分类，该步骤进一步包括以下子步骤：

S2-1、读取待分类遥感影像的灰度值到矩阵                                                

中，根据待分类遥感影像的特征确定其聚类数目n，并随机确定各聚类中心的初始值

，且

；其中：即为待分类遥感影像的元胞空间，M、N分别为待分类遥感影像的长、宽像素大小，n为自然数且

，k取1、2、…n，i取1、2、…M，j取1、2、…N，和

分别表示待分类遥感影像所有像素灰度值中的最小灰度和最大灰度；

S2-2、根据公式

计算t时刻待分类遥感影像中各中心元胞

距各类聚类中心的距离

，如果中心元胞

在t时刻到第m类聚类中心的距离

小于其到其他聚类中心的距离，则该元胞

属于第m类聚类中心，即t时刻元胞

的状态值

；根据公式

计算t时刻待分类遥感影像中各中心元胞

的邻居元胞到各类聚类中心的距离

，中心元胞

的邻居元胞在t时刻到第m类聚类中心的距离

小于到其他聚类中心的距离，则该邻居元胞属于第m类；其中：

为各中心元胞

的灰度值，

为各中心元胞的邻居元胞的灰度值，

为t时刻第k类聚类中心的灰度值，i取1、2、…M，j取1、2、…N，k取1、2、…n，m为自然数且

，||…||表示距离运算；

S2-3、根据步骤S1中得到的邻居元胞对中心元胞的隶属度，并根据转换规则确定各中心元胞在t+1时刻的状态，所述的转换规则如下：

；

其中：

和分别为t+1时刻和t时刻中心元胞

的状态；

为中心元胞

的圆形邻居元胞中属于第k类聚类中心的最大隶属度，i取1、2、…M，j取1、2、…N，k为自然数且

；

为所有元胞的圆形邻居属于某一类聚类中心的最大隶属度和的集合；

为根据遥感影像质量和分类目所设定的阈值，其值介于0到1之间； 

S2-4、根据公式求取t+1时刻待分类遥感影像的新的聚类中心

；其中：

为t+1时刻第k类聚类中心的灰度值，k取1、2、…n；

为中心元胞

的灰度值，i取1、2、…M，j取1、2、…N；

为t+1时刻属于第k类聚类中心的中心元胞灰度集，k取1、2、…n；

为在t+1时刻属于第k类聚类中心的元胞数；

S2-5、判断t+1时刻和t时刻待分类遥感影像中同一类聚类中心的接近程度，如果公式

成立，k取1、2、…n，或者当迭代次数超过预先设定的最大迭代次数时，则迭代结束，执行步骤S3，否则，重新执行步骤S2-1；其中：

为设定的收敛阈值；

S3、根据步骤S2所得的待分类遥感影像的聚类中心结果，分别将各类聚类中心值取整后赋给不同类别的所有元胞，即得到所有元胞的状态矩阵，并输出状态矩阵，即得到已分类的结果影像。

2.根据权利要求1所述的基于元胞自动机的遥感影像非监督分类方法，其特征在于：

所述的步骤S1-2中栅格空间中邻居元胞对中心元胞的隶属度通过隶属度函数得到，其中，隶属度函数为高斯型隶属度函数。

3.根据权利要求1或2所述的基于元胞自动机的遥感影像非监督分类方法，其特征在于：

所述的步骤S2-2中，中心元胞距各类聚类中心的距离

和中心元胞的邻居元胞距各类聚类中心的距离

均为欧式距离。

Images