基于分形纹理的遥感影像土地覆盖的分类方法研究

测绘科学

ScienceofSurveyingandMapping

Vol󰀁33No󰀁2

Mar󰀁

基于分形纹理的遥感影像土地覆盖的分类方法研究

王󰀁娟

󰀂

,张󰀁军

󰀂

,吕兆峰

!

(󰀂辽宁工程技术大学,辽宁阜新󰀁123000; 中国测绘科学研究院,北京󰀁100039;!西安通信学院,西安󰀁710106)∀摘󰀁要#提出一种基于分形理论和改进模糊C均值聚类的遥感图像非监督分类方法,该方法尝试将图像的光谱

信息和纹理特征相结合。将图像进行主分量变换,根据第一主分量计算图像的布朗运动的各方向的分形维数,差分盒维数和∃空隙%等纹理特征作为分类依据。采用改进的模糊C均值聚类,并用混淆矩阵方法评定分类结果精度。通过对试验区的分类试验,说明该方法对改善土地覆盖分类精度行之有效。∀关键词#纹理;分形;布朗运动;差分盒维数;空隙;模糊C均值;混淆矩阵

∀中图分类号#TP751󰀁󰀁󰀁󰀁∀文献标识码#A󰀁󰀁󰀁󰀁∀文章编号#1009󰀂2307(2008)02󰀂0015󰀂04DOI:10󰀁3771/j󰀁issn󰀁1009󰀂2307󰀁2008󰀁02󰀁005

1󰀁引言

在自然资源调查中,遥感影像数据已成为重要的数据源和地物光谱分析识别的重要依据,是进行土地覆盖变化动态监测中的重要技术手段[1]。目前提取影像信息的方法是基于地物光谱特征信息,在遥感影像中,像元的亮度值代表该像元地物的平均辐射值,它随着地物的成分、纹理、状态表面特征及使用电磁波段的不同而变化。但在实际应用中常会遇到地物类别光谱反射率混合问题,难以区分光谱混淆的地物,影响地物分类精度,如林地与草地等。为了比较有效地解决遥感影像分类中的混淆地物信息问题,本文选择基于分形[2]的纹理特征作为分类方法。该方法具有尺度不变性的重要特点,即在一定范围内,分形特征是于分辨率、视角而稳定存在的物质表示量,利用这一特性作为分类的主要依据进行特征信息提取,提高分类精度。

近些年来这方面的研究从未间断过,如黄桂兰等提出了基于多尺度分形分析和基于纹理质地子特征的影像纹理分析方法[1995];张继贤等提出基于分形几何的多分辨率纹理分析[1995];曹燕等用三种计算影像分维值的算法和数学模型,结果表明:综合运用三种方法进行纹理分析,其提取精度达到75󰀁6%[1997];舒宁讨论了影像纹理分析的分维估计原理以及当前应用的主要方法,以分维值作彩色密度分割结果提取地物[1998]。但其中有的研究只是就方法本身进行可行性试验,还未达到实用化,如只是从影像中取出一些窗作分维计算和分类,而未对整幅影像每一处进行分维估计和分类。如何在遥感影像纹理中实施基于分形理论的分维估计方法来提高分类精度是本文所要讨论的重要问题。另一方面,绝大多数遥感影像纹理分析的文章中只是针对单一波段,并未考虑如何针对多波段影像纹理分析,这也是本文所要讨论的。

分形布朗模型能够反映影像纹理整体和局部的起伏程度,因此,本文引入了该模型的分数维作为纹理特征识别依据。然而自然界有很多视觉上差别很大的纹理其分数维近似相同情况,为克服这一缺陷,引入了差分盒维数及其∃空隙%特征值概念,可有效辨析具有相同分数维但差别较大的纹理,改善对地物区分精度具有明显作用。

2󰀁1󰀁离散分形布朗模型(DiscreteFractalBrownianModel

-DFBM)

遥感影像作为数字图像是呈二维离散的数据,根据这一特性我们可用离散分形布朗随机场模型来定义,即

lgE|BH(x+󰀁x)-BH(x)|2=2Hlg&󰀁x&+lgC(1)󰀁󰀁在(1)中通过取不同的󰀁x,运用最小二乘法拟合数据点对可求得直线的斜率H,图像的分形维数D=3-H。为提供足够信息加以描述和识别自然纹理,引入了DFBM模型的各方向维数作为分类特征描述,其各方向的分维数计算[3]如下:

1)水平分维D1󰀁在灰度图像上,当式(1)中󰀁x为水平方向,k为整数变量,图像窗口大小为M∋N

M-1N-k-1x=0y=0

F1k=

((

fx,y-fx,y+kMN-k2

(2)

󰀁󰀁则式(1)可改写为

logF1k=2󰀂1logk+c(3)

对式(3)做最小二乘线性回归计算,可得H1,水平分维D1=3-H1

2)垂直分维D2󰀁当式(2)中󰀁x为垂直方向,

M-K-1N-1x=0y=0

F2 =

((

fx,y-fx+k,yN(M-K)2

(4)

2󰀁遥感影像纹理分类的基本数学模型

目前利用分形的方法描述遥感影像纹理有很多不同的模型,如利用表面积-体积关系、分形交集公式、计盒法原理、分数布朗表面的统计性质都可进行分维估值。但用

󰀁󰀁可求得D2=3-H2

3)对角分维D3当式(2)中󰀁x为对角方向,

M-k-1N-k-1

F3k=

((

x=0

y=0

fx,y-fx+k,y+k2M-kN-k2

+

2

M-k-1N-k-1x=0

((

y=k

fx,y-fx+k,y-k2M-kN-k(5)

作者简介:王娟(1981󰀂),女,辽宁工

程技术大学硕士研究生,主要从事地理信息系统及遥感影像、数字图像处理方面等的研究。

E󰀂mai:lwangjuankai23@126󰀁com收稿日期:2006󰀂11󰀂09

󰀁󰀁可求得D3=3-H3

4)综合分维D4

M-k-1N-1x=0y=0M-1N-k-1

F4k=

((((

fx,y-fx+k,y2

N(M-k)+M(N-K)

fx,y-fx,y-k2

x=0

y=k

+󰀁󰀁󰀁

N(M-k)+M(N-K)

(6)

16

测绘科学󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁第33卷

󰀁󰀁可求得D4=3-H4

2󰀁2󰀁差分盒(DBC)分形维及其∃空隙%特征(lacunarity)

盒维数(boxdimension)是应用很广泛的分形维数概念之一,最早由前苏联数学家柯尔莫洛夫(Kolmogorov)提出,其方法的出发点是利用不同的观察尺度来近似分形图形,并由此得到图形的分形维数。为了更加有效地对遥感影像特征信息分割提取,对盒维数[4]进行改进,即在三维空间中,沿x,y,z方向堆放一堆并非固定大小的立方体的盒子,这堆盒子足够多,覆盖了整个图像曲面,所有与图像曲面相交的盒子数就是该曲面的盒维数Nr(F),计算方法如下:

Nr=

图2󰀁不同窗口分形图表1󰀁空隙特征值

8∋8窗口类型1(红色)类型2(墨绿色)类型3(鲜绿色)类型4(浅青色)类型5(深蓝色)类型6(青色)

空隙特征值

C(2)-0󰀁4090-0󰀁6692-0󰀁14230󰀁3045-0󰀁44180󰀁8857

C(4)0󰀁63140󰀁67810󰀁66710󰀁78690󰀁75250󰀁9828

C(8)0󰀁97450󰀁97420󰀁97870󰀁97310󰀁98860󰀁9986

(

i,j

maxIk-minIk/r

(7)

󰀁󰀁其中Ik(k=1,2,),n)表示r∋r单位面积内n个

像素灰度值的大小。对于不同的r值,由式(7)可求得一组点ri,Nr(i=1,2,3,),m),根据文献[4]利

i用最小二乘法进行线性回归即可得到图像曲面的DBC盒维数。针对遥感影像地物特征具有相同的分数维而差异很大的纹理情况,这里引入了∃空隙%[1]特征概念解决出现的问题。在前述DBC分形维计算方法的基础上,设P为遥感图像S所含像素点数,建立一种简单有效的计算空隙特征的公式:

P2-N2r

(8)CL=2

P+N2r

3󰀁基于纹理的分形模型的建立

3󰀁1󰀁影像波段的处理

为了使遥感影像与试验区时态同步,选取影像是一幅CBERS󰀂1卫星遥感影像图,地面分辨率为19󰀁5m。因CBERS󰀂1影像是多波段的,若每个波段都计算其分形维数并用于分类,会导致信息冗余,反而使分类精度下降。主分量变换[5]法就很好地解决了上述数据冗余问题,因其第一主分量包含了各个波段的绝大多数信息,具体结果见图1。

3󰀁3󰀁纹理特征值的提取

因对分类实验区没有先验知识,即属于非监督分类。依据地物的光谱特征的分布规律,在影像试验区中选取不同的地物类别(见图3)。

图3󰀁地物典型样本示例

每一地物类型分别采20个样本,运用第2节中所提到的数学模型,来提取各地类型的纹理特征值。并取其平均值,其具体数值见表2。

表2󰀁各地类纹理分形特征值及其空隙特征值

特征值各地类型类型1(红色)类型2(墨绿色)类型3(鲜绿色)类型4(浅青色)类型5(深蓝色)类型6(青色)

离散分形布朗模型(DFBM)水平分维

垂直分维

对角分维

综合分维

差分盒空隙(DBC)(lacunarity)维数特征值

-0󰀁4090-0󰀁6692-0󰀁14230󰀁3045-0󰀁44180󰀁8857

2󰀁40872󰀁50442󰀁51252󰀁44322󰀁33442󰀁36342󰀁32412󰀁36062󰀁302󰀁18372󰀁27962󰀁47662󰀁50802󰀁34092󰀁26682󰀁12042󰀁18412󰀁13452󰀁11762󰀁00052󰀁622󰀁48872󰀁552󰀁53522󰀁392󰀁21262󰀁30512󰀁14872󰀁12292󰀁0835

图1󰀁主分量变换图

3󰀁2󰀁滑动窗口的选取

图像中每一个像素点的分维数,都是基于纹理的。因此通常采用一定大小的滑动窗口,计算窗口内的分维数值,作为滑动窗口中心点的分维数。可见滑动窗口大小的选择十分重要,窗口选择得过小,不能充分反映出该局部范围的纹理特征;而选择得过大,会使其分形特征降低得到错误的结果,同时由于该算法是n2型的,即运算时间与窗口大小的平方成正比[6],选择过大的窗口会导致程序效率的急剧下降。要保证窗口的大小合适,就是要保证回归方程有较好的直线性。图2是以图1a图为例,按照不同的窗口大小,根据式(7)得到DBC盒维数分形图,从中可以看出选择合适窗口的重要性。

试验表明:窗口的大小选择6~10时,直线拟合性较好,从各种类型的纹理特征考虑直线拟合的优劣,对于DFBM维数和DBC盒维数来说,通过计算分析比较,当窗口选择8∋8的时候,直线拟合性最好。而对于由DBC维数衍生的空隙特征C(L)来说,当L值较小时的C(L)反映的纹理的细节和结构信息就越丰富。利用式(8)计算C(L)的值如表1,通过试验分析比较,选取L=2。

3󰀁4󰀁分形图像的生成

对遥感图像在基于分形的理论上对其进行分类,就须知道每个像素的分形维数。这里运用滑块邻域操作技术,计算图像上每个像素点的分形维数。考虑到窗口中心要遍历整个图像,当窗口在图像的周边移动时,就需在图像的周边填充像素。经过试验反复比较,用复制外边界[4]的方法来填充周边像素,并依据式(2)、(4)、(5)、(6)、(7)、(8)计算影像的各分形特征图。从图4中可以看到,各分形图从不同的角度反映了影像纹理的变化趋势。

4󰀁改进模糊聚类的影像识别

在分类过程中,由于遥感信息的不确定性和混和像元问题使部分像元很难进行非此即彼的划分。因此,需要考虑像元类别的隶属度问题进行软划分。模糊C均值聚类[7]就是进行非监督分类的一种很好的软化分方法。一方面,󰀁第2期󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁王󰀁娟等󰀁基于分形纹理的遥感影像土地覆盖的分类方法研究17

次实验的样区内的信息是由长江三峡库区(奉节-巴东段)-移民工程遥感动态监测图集提供的,其遥感数据源为ETM+,分辨率分别为15m。将其实际所属类别与分类结果进行比较,得出它们的混淆矩阵,统计结果为表3和表4。

表3󰀁基于分形方法的分类

精度分类类型林地草地旱地新开发用地

水面滩涂

生产者精度/%

88󰀁2484󰀁0081󰀁0291󰀁1682󰀁3479󰀁99

用户精度/%总体分类精度/%

85󰀁7187󰀁5081󰀁0088󰀁5782󰀁3380󰀁00

84󰀁19

图4󰀁遥感影像各分形特征图

该算法具有良好的局部收敛特性;另一方面,它适合在高维特征空间(本实验研究中分形模型具有6维特征空间)中进行像素的分类。

本文采用了改进的模糊C均值聚类法,更好地发挥其在遥感数据分类方面的优势,其基本思想为定义目标函数J

JmU,V=

k=1i=1

表4󰀁基于ERDAS系统非监督分类方法

精度分类类型林地草地旱地新开发用地

水面滩涂

生产者精度/%

59󰀁5680󰀁0076󰀁6775󰀁0082󰀁5161󰀁43

用户精度/%总体分类精度/%

56󰀁6376󰀁00󰀁2983󰀁7180󰀁4172󰀁33

72󰀁23

((

nc

uikm

dik2(9)

在dik的计算中,为适合遥感数据特点[8],本文采用Mahalanobis距离度量方式,其中m为模糊加权指数。由于篇幅关系这里不对该函数分解进行赘述。

试验区为三峡库区巫山县,见图1a,运用改进的分类器对其进行分类。通过事后对各类地物的光谱响应曲线进行分析,并参考长江三峡库区(奉节-巴东段)-移民工程遥感动态监测图集,可知图3中的各类型分别代表:类型1∗林地,类型2∗草地,类型3∗旱地,类型4∗新开发用地,类型5∗水体,类型6∗滩涂。其各地物类型用索引图像[5]来表示见图5a,分类结果为5b图,运用ERDAS系统中的非监督分类的结果得图5c。

从两表中可以看出,在表4中,由于仅仅依靠遥感图像

的光谱特征进行地物分类,地物类型便不能很好的区分。只有新开发用地和水面的精度较高,其他地物类型都比较低,总分类精度才72󰀁23%。而表3的分类结果显示,可较好的区分∃同谱异物%的像元问题,使得分类精度有很大提高。

6󰀁结束语

分形维数能有效地表现图像纹理的复杂度和粗糙度,能够揭示地物类型内在的自相似性。本文以分形理论为基础,提出一种遥感影像的地物分类方法。它将多波段遥感影像进行主分量变换,根据得到的第一主分量图计算分形布朗运动各方向的分形维数,DBC盒维数和∃空隙%等纹理特征,采用改进的模糊C均值聚类,解决了影像中的不确定性问题。与传统的基于光谱特征的分类方法相比,本文提出的方法具有一定的实用性和较好的分类效果。

该方法也存在许多有待改进的方面,如分形维数的稳定性问题,在实际遥感图像中,即使是同一类图像,所选样本图像不同,分形维数往往也有差异;另该方法会随着分类精度的提高而增加计算量。因此,在今后研究中可以考虑如何提高分形维数的稳定性并对计算进行优化,以进一步提高该方法的效率,使其更加完善。

参考文献

[1]

陈美球,赵小敏,吴次方󰀁城镇用地扩张的遥感监测研究进展[J].江西农业大学学报,1999,21(2):237󰀂241󰀁

MandelbrotBB󰀁TheFractalGeometryofNature[M].SanFrancisco:Freeman,1982󰀁[3]王思贤,江炳,寥孟扬󰀁基于小波变换及布朗分形的纹理分割[J].武汉大学学报,1998,44(1):118󰀂120󰀁SarkarN,ChaudhuriBB󰀁AnEfficientApproachtoEsti󰀂mateFractalDimensionofTexturalImages[J].PatternRecognition,1992,25(9):1035󰀂1041󰀁RafaelCGonzalez,RichardEWoods,StevenLEddins󰀁Dig󰀂italImageProcessingUsingMATLAB[Z].(下转第32页)

图5󰀁分类结果

在图5c中,说明水体、草地和新开发用地的分类效果

较好,而林地,旱地以及滩涂的误分较大,而且边缘刻画得很模糊,与实际相差较远。,在图5b中各地物类型均被正确识别出来,且地物边界清晰可见。

[2]

[4]

5󰀁精度评定

本文对遥感图像的分类结果进行精度评定的方法采用混淆矩阵法[10],它已被证明是一种有效的精度评定方法。方法如下:在遥感图像中的训练区内选取1000个采样点,此[5]

32

!当Kf<0.1时,计算S01

0.1∋Kf∋BTl 1

S01=2!∀0.01∋BT1B1+∀∋Kf+计算岭参数

测绘科学󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁󰀁第33卷

2

5󰀁结束语

(8)

由以上分析和算例可以看出:󰀂斜率法解算岭参数可

以和判断方程病态性一块进行,病态性的判定和解算形成一个整体; 要求病态曲线的负斜率不大于临界点(Kf=0.1)的负斜率是合理的,有效克服了系统病态;!斜率法计算简单、几何意义明确。

参考文献

[1]

陈希孺,王桂松󰀁近代回归分析/原理方法与应用[M].合肥:安徽教育出版社,1987󰀂9󰀁

[2]王桂松󰀁线性模型的理论及其应用[M].合肥:安徽

教育出版社,1994󰀁

[3]黄维彬󰀁近代平差理论及其应用[M].北京:

出版社,1992󰀁[4]Golub,etal󰀁GeneralizedCross󰀂validationasamethodfor

choosingagoodridgeparameter[J].Technometrics,1979,21(2):215󰀂223󰀁

[5]HansenPC󰀁Analysisofdiscreteill󰀂posedproblemsby

meansoftheL󰀂curve[J].SIAMReview,1992,34(4),561󰀂580󰀁

[6]HansenPCO.LearyDP󰀁TheuseoftheL󰀂curveinthe

regularizationofdiscreteill󰀂posedproblems[J].SIAMJ󰀁Sci󰀁Comput,1993,14(6):1487󰀂1503󰀁

[7]王振杰,欧吉坤󰀁用L曲线法确定岭估计中的岭参数

[J].武汉大学学报,2004,29(3):235󰀂238󰀁

[8]王振杰󰀁大地测量中不适定问题的正则化解法研究

[D].中国科学院测量与地球物理研究所,2003󰀁

[9]欧吉坤󰀁测量平差中不适定问题解的统一表达式与选

权拟合法[J].测绘学报,2004,33(4):283󰀂288󰀁[10]周江文,黄幼才,杨元喜,等󰀁抗差最小二乘法

[M].武汉:华中理工大学出版社,1997󰀁

[11]丁士俊,等󰀁平差不适定问题解性质与正则参数的确

定方法[J].测绘科学,2006,31(2)󰀁

∀=0

&y1&2

Skf==S01(9)

!∀󰀁󰀁当存在Kfi<0.1,i=1,2,)q;q∀=max∀i(10)

4󰀁算例

表1󰀁岭参数∀

病态

指标Kf

0󰀁010󰀁020󰀁030󰀁040󰀁050󰀁060󰀁070󰀁080󰀁090󰀁1

岭参数∀0󰀁5200󰀁6920󰀁7780󰀁8030󰀁7770󰀁7040󰀁50󰀁4320󰀁2360

不同病态指标下解算的岭参数∀见表1和图1,−.

所示,∀所在的点避开了曲线急剧变化的部分,有效地减小了&x&的值,证明岭参数∀的大小是合理的。

解算不同病态指标下岭估计&x∀&和LS估计&x&(图2),可以看出:󰀂岭估计改善了LS估计; 岭估计是一种压缩估计,病态性越严重,压缩比例越大。

表2󰀁不同岭参数和&x&2比较

方法斜率法L󰀂曲线法LS估计

岭参数∀0󰀁2922󰀁0940

&x&0󰀁1230󰀁06910󰀁366

图2󰀁不同Kf条件下岭估计

和LS估计的比较

Theslopemethodofridgeparameterdetermination

Abstract:Thekeytotheridgeestimationistoselectproperridgeparameter󰀁Theruleoftheparameters.normvaryingwiththeridgeparameterisanalyzedandstudiedindifferentill󰀂conditionedsystem󰀁Thentheslopemethodisputforwardtodeterminetheridgeparameter󰀁NumericalsimulationofslopemethodandL󰀂curve,LSestimation,demonstratesthattheslopemethodissimpleandexplic󰀂it,andtheridgeparameterdeterminedisreasonable󰀁

Keywords:ridgeparameter;parameterestimation;ill󰀂conditionedindex;slopemethodZHANGHui󰀂zhan,DUZhi󰀂xing(ShandongUniversityofScienceandTechnology,Qingdao2665110,China)

(上接第17页)

[6]斯林󰀁利用分形方法提取遥感图像空间结构信息的应用研究

[D].中国林业科学研究院,2000󰀁

[7]PothA,KlausD,VobM,etal󰀁OptimizainatMulti󰀂spectral

LandCoverClassificationwithFuzzyClusterningandtheKohonenFeatureMap[J].Int󰀁J󰀁RemoteSensing󰀁2001,22(8):1423󰀂1439󰀁

[8]

哈斯巴干,马建文,李启青,等󰀁模糊C-均值算法改进及其对卫星遥感数据聚类的对比[J].计算机工程,2004,30,(11)󰀁[9]廖克,等󰀁高分辩率卫星遥感影像在土地利用变化动态监测

中的应用[J].测绘科学,2006,31(6)󰀁

[10]StevenEFranklinandBradleyAWilson󰀁Athree󰀂stageclassifier

forremotesensingofmountainenvironments[J].PE&RS󰀁1992,58(4):449󰀂4󰀁

Studyonclassificationoflandcoverwithremotesensingimagebasedonfractaltexture

Abstract:ThispaperpresentsaunsupervisedclassificationmethodofremotesensingimagebasedonFractalTheoryandameliora󰀂tedFuzzyC󰀂MeansClustering(FCMC),whichtakesadvantagesofbothspectralandtexturalfeatures󰀁Firstly,principal󰀂componenttransformationisperformedwiththehelPofseveralclassificationelementssuchasvariousdirectionsfractaldimensionsofFractalBrownMotion,DifferenceBox󰀂Countingdimensionandlacunarityandsoon󰀁Then,withimprovedFCMCasclassifier,theprecisionofclassi󰀂ficationresultswithConfusion󰀂Matrixisassessed󰀁Finally,theexperimentofclassifyingthelandcoverofremotesensingimagehasbeendoneandtheresultissatisfactory,whichverifiestheeffectofthismethod󰀁

Keywords:texture;fracta;lBrownmotion;differencebox󰀂countingdimension;lacunarity;fuzzyC󰀂means;confusion󰀂matrix

󰀂 󰀂 !

WANGJuan,ZHANGJun,L󰀁Zhao󰀂feng(󰀂LiaoningTechnicalUniversity,Fuxin123000; ChineseAcademyofSurveyingandMapping,Beijing100039,China;!Xi.anCommunicationInstitute,Xi.an710106,China)