基于半监督流形学习的多光谱遥感数据土地利用分类方法
技术领域
[0001] 本发明涉及一种基于多光谱遥感数据土地利用分类方法,特别涉及一种半监督流形学习的适用于多光谱数据土地利用分类的方法。
背景技术
[0002] 人们从20世纪60年代就开始用多光谱技术获取地球表面信息。多光谱成像是利用多光谱摄影系统或多光谱扫描系统对电磁波谱不同谱段做同步摄影遥感,分别获得植被及其他地物在不同谱段上的影像的遥感技术。多光谱扫描仪(Multi-spectral Scanner,MSS)是一种光谱测量传感器,它通过记录地物对不同光谱波段的响应来获取地物信息。多光谱成像是通过扫描的方式完成的。对于每一个分辨率大小的地面区域,MSS一次可获取多个光谱响应值,每个值对应一个光谱波段。国际遥感界的共识是光谱分辨率在λ/10数量级范围的称为多光谱(Multi-spectral),这样的遥感器在可见光和近红外光谱区只有几个波段,如美国Landsat MSS,TM,法国的SPOT等。多光谱遥感图像的特点:①工作波段宽,从近紫外、可见光到热红外波段,波长范围达0.35~20微米;②各波段的数据容易配准。这两个特点非其他遥感器所能具有,因而多光谱扫描仪是气象卫星和“陆地卫星”的主要遥感器。随着传感器技术的进步,新型多光谱传感器的光谱分辨率达到了纳米级,空间分辨率达到了米级。在相关领域具有广泛的应用和发展空间。在遥感图像的实际使用中,数据集较以往有了显著的变化,其主要特点可以归纳为:高数据量、高维数、高数据增长率以及高分辨率。
[0003] 随着多光谱传感器的空间分辨率不断地提高,人们可以从遥感图像中获得更多的有用的数据和信息。多光谱图像数据也越来越多的被用于土地利用和覆被监测领域中来,这就不可避免的遇到一个关键环节---多光谱图像分类,用计算机对多光谱遥感图像进行土地利用分类时,多光谱图像数据一般数据量较大,把不同传感器得到的遥感图像进行融合使用时数据量更大。多光谱图像数据的土地利用分类,通常是建立在不同地物在各波段反射的电磁波谱差别的基础上的。若以各波段接收到的电磁波强度为坐标,则n个波段可形成n维波谱空间。各波段上同一像元对应于N维空间的一个点,而同类地物可形成一个点集。多光谱遥感数据土地利用分类的基本原理在于把波谱空间中的点集区分开来,多光谱图像分类中,由于数据维数较高、有类别标签的地面样本获取代价比较大,所以有类别标签的训练样本数量通常很有限,传统的模式识别方法(如贝叶斯方法,神经网络等)大多是以经典统计学为基础,在假设训练样本数目足够多的前提下进行研究的,这些方法只有当训练样本数目趋于无穷大时其性能才能达到理论上的最优。目前多光谱遥感数据土地利用分类主要集中两步,首先对多光谱遥感数据进行分类前的处理---特征提取,这一步主要目的是既可以消除原图像中对分类或其他处理不利的外部影响因素,又可以消除原图像各波段间的数据相关性从而避免许多无效的数据处理。其次再选择合适的分类器进行多光谱遥感数据土地利用分类。
[0004] 针对高光谱图像分类问题,研究者提出如下两方面技术,一是把机器学习引入到多光谱图像分类中,力图在多光谱数据中通过机器学习发现数据内部的变化规律;二是寻求或改进更有效的分类算法;按照地物分类时是否需要训练样本,这些方法可大致分为两类:监督分类和非监督分类。由于多光谱图像有类别标签的训练样本获取的代价比较大,样本数目非常有限,所以使用监督分类效果就有些局限性,虽然非监督分类不使用任何关于数据的先验信息,但是分类效果通常较差。因此,把半监督流形学习方法引入进多光谱遥感数据土地利用分类中,即通过对训练样本的学习,掌握数据内部隐含的规律,实现对未来样本的准确判决。
[0005] 基于数据的机器学习是现代智能技术中的重要方面,研究从观测数据(训练样本)出发寻找尚不能通过原理分析得到的规律,利用这些规律对未来数据或无法观测的数据进行预测。现实世界中存在大量目前无法准确认识但却可以进行观察的事物,因此,把机器学习引入到多光谱遥感数据土地利用分类中来是十分有必要的,数据经过投影,通过半监督流形学习法来发现数据内部的变化规律,提高识别精度,这样既比非监督方法精度要好,又可以弥补使用监督方法时样本获取代价过高的缺点。基于数据的机器学习在现代科学、技术、经济、社会等各个领域都发挥着十分重要作用。
发明内容
[0006] 有鉴于此,为了解决上述问题,本发明提出一种把机器学习引入到遥感图像分类中,融合半监督学习和流形学习的方法,对投影后的多光谱数据的进行训练来找出数据内部的流行变化规律;同时将用户先验知识提供的标注信息融入到半监督流形学习中,并用于多光谱遥感图像土地利用分类。
[0007] 本发明的目的是提出一种基于半监督流形学习的多光谱图像土地利用分类方法;
[0008] 本发明的目的是这样实现的:
[0009] 本发明提供的基于半监督流形学习的多光谱图像土地利用分类方法,包括以下步骤:
[0010] (1)读入多光谱遥感图像数据;
[0011] (2)读入多光谱遥感图像中选取多个样本数据点,对每一个样本数据点根据其波段生成一个向量,然后将所得向量组成矩阵,作为样本数据集;
[0012] (3)根据先验知识从样本数据集中选取部分样本数据进行已知地物类别的标注,生成样本类别标签;
[0013] (4)通过样本类别标签来构建度量数据点的相似性与相异性的相似图和相异图;
[0014] (5)利用热核构建相似图的权重矩阵;
[0015] (6)根据样本数据集和权重矩阵来计算局部相似结构矩阵和全局相异结构矩阵;
[0016] (7)通过全局相异结构矩阵和优化目标函数来计算多光谱遥感图像的投影矩阵;
[0017] (8)通过投影矩阵对整幅多光谱遥感数据进行投影;
[0018] (9)对投影后的整幅多光谱遥感数据用K-近邻分类算法进行分类;
[0019] (10)输出多光谱遥感数据地物分类图。
[0020] 进一步,所述步骤(2)中的训练样本集通过以下方式实现:
[0021] 根据不同波段对地物的光谱反射特性,把多光谱遥感图像数据变为三维数据,X=T
{x1,x2,…,xM×N},其中M×N是图像的空间尺寸,B表示波段数,矩阵的每一列表示一个波段的数据值,每一行表示多光谱遥感图像数据中的一个点。
[0022] 进一步,所述步骤(3)中的样本类别标签通过以下方式实现:根据先验知识从训练样本集中选取部分样本数据进行已知地物类别的标注,生成样本类别标签:
[0023] X={(x1,l1),(x2,l2),…,(xC,lC),xC+1,xC+2,...,xM×N}T
[0024] 式中li为对数据点xi标注的地物类别标签,前C个点具有类别信息,其余的M×N-C个为无类别标签样本点。
[0025] 进一步,所述步骤(4)中的相似图和相异图包括以下步骤:
[0026] 1)将样本数据点分成相似图GS和相异图GD,在相似图GS和相异图GD的节点之间采用k最临近方法k-NN方法定义两个图中的边;
[0027] 2)每一点xj对xi,若xj∈knnS(xi),即每一点对(xi,xj),两者为近邻且相似度较高,则用一条边连接图GS中xi和xj两点;
[0028] 3)若xj∈knnD(xi),即相似度较低,则用一条边连接图GD中xi和xj两点;
[0029] 式中 为点xi的k个近邻点组成的子集,对每一点xi,近邻数
据集knn(xi)可以分为两部分:knnD(xi)与knnS(xi),其中:knnD(xi)表示近邻点来自于不同类样本点,knnS(xi)则为knn(xi)中其余部分,需要进一步说明的是,knnS(xi)中部分样本可能不包含类别信息,但是这些样本距离xi足够的近,knn(xi)同knnD(xi)和knnS(xi)之间的关系可用下式来表示:
[0030]
[0031] knnS(xi)=knn(xi)-knnD(xi)
[0032] 显然, knnS(xi)∪knnD(xi)=knn(xi);
[0033] 根据knnD(xi)和knnS(xi),就可以构建相似图GS和相异图GD。
[0034] 进一步,所述步骤(5)中权重矩阵通过以下方式实现:
[0035]
[0036] 权重矩阵为Wij=Ws,ij*Aij
[0037] 其中,Aij定义如下:
[0038]
[0039] 其中d2(xi,xj)是从Xi到Xj的距离,t为一常数。
[0040] 进一步,所述步骤(6)包括以下步骤:
[0041] (1)计算相似局部结构矩阵Srlw
[0042] 根据输入的多光谱样本数据集X和权重矩阵Wij,通过
[0043]
[0044]
[0045]
[0046]
[0047] 式中DS是一个对角阵,且Ds,ii=∑jWij,LS=DS-WS为Laplacian矩阵;
[0048] 然后,防止Slw为奇异矩阵采取如下方式,
[0049] Srlw=Slw+SE
[0050] 式中SE为β×Ε即一个较小的常数β乘以单位矩阵;
[0051] (2)计算全局相异结构矩阵Srlb
[0052] 根据输入的多光谱样本数据集X和权重矩阵Wij,通过
[0053] Srlb=Slb+ST
[0054] 其中ST,Slb定义如下:
T
[0055] ST=E{(X-μ)(X-μ)}
[0056]
[0057] 式中μi为每类的均值;μ为总样本的均值;n为每类含有的样本数;C是总类别数。
[0058] 进一步,所述步骤(7)通过以下方式实现:
[0059] 多光谱遥感图像数据投影后最理想的目标:使同类样本数据间散度尽可能小、不同类数据间散度尽可能大,基于上述考虑,投影矩阵A通过目标优化函数:
[0060] Maximize J(A)=tr(ATSrlbAT)
[0061] subject to tr(ATSrlwAT)=1
[0062] 通过采用Lagrange函数ATSrlbAT-λ(ATSrlwAT-1),上述最优化问题可以很容易地转换为广义的特征值求解问题,基于上述考虑,投影矩阵A通过目标优化函数:
[0063]
[0064] 求得,即A由以下特征方程的d个最大特征值λ1>λ2>…>λd对应特征向量v1,v1,…,vd构成
[0065] Srlbv=λSrlwv。
[0066] 进一步,所述步骤(8)通过以下方式来实现:
[0067] Y=ATZ;
[0068] 式中Y为整幅多光谱样本数据集Z通过投影矩阵A投影到嵌入空间的数据集。
[0069] 本发明的优点在于:本发明利用半监督流形学习方法,且有区别的对待标注数据与无标注数据,首先对样本数据进行训练,发现其内部数据变化规律,找出投影矩阵,把整幅图像进行投影,使不同的类尽可能的分开相同类尽可能的聚拢,再用最近邻方法去实现多光谱遥感数据土地利用分类;具有以下优点:
[0070] (1)通过加入不知类别信息的样本数据,来得到投影矩阵,这样在没有降低分类正确度的情况下,大大节约了选取训练样本的成本;
[0071] (2)利用半监督流形学习方法对多光谱遥感数据进行学习,发现隐藏在多光谱遥感数据中的低维流形结构;
[0072] (3)通过有区别的对待用户标注数据和无标注数据,人的视觉解译信息融入半监督流形学习方法,从而得到一个结合了用户语义理解的多光谱遥感数据的投影矩阵,在投影后的空间中实现同类数据点之间依然保持其近邻关系,不同类数据点之间的距离尽可能最大化,最大程度的增加了不同类型地物目标之间的可分性,从而实现多光谱遥感数据土地利用分类。
[0073] 本发明的其它优点、目标和特征在某种程度上将在随后的说明书中进行阐述,并且在某种程度上,基于对下文的考察研究对本领域技术人员而言将是显而易见的,或者可以从本发明的实践中得到教导。本发明的目标和其它优点可以通过下面的说明书,权利要求书,以及附图中所特别指出的结构来实现和获得。
附图说明
[0074] 为了使本发明的目的、技术方案和优点更加清楚,下面将结合附图对本发明作进一步的详细描述,其中:
[0075] 图1为基于半监督流形学习的多光谱遥感数据土地利用分类方法流程图;
[0076] 图2算法分类精度图样本数据投影前的数据分布情况;
[0077] 图3训练样本集投影前二维分布图样本数据经过投影矩阵投影后的分布情况;
[0078] 图4训练样本集投影后二维分布图算法实现分类后输出效果图;
[0079] 图5本专利方法最终输出河流湖泊(水域)效果图;
[0080] 图6本专利方法最终输出森林用地效果图;
[0081] 图7本专利方法最终输出建筑用地效果图;
[0082] 图8本专利方法最终输出水田耕地效果图;
[0083] 图9本专利方法最终输出绿地效果图;
[0084] 图10本专利方法最终输出边界效果图;
[0085] 图11本专利方法最终输出整体效果图。
具体实施方式
[0086] 以下将结合附图,对本发明的优选实施例进行详细的描述;应当理解,优选实施例仅为了说明本发明,而不是为了限制本发明的保护范围。
[0087] 实施例1
[0088] 图1为基于半监督流形学习的多光谱遥感数据土地利用分类方法流程图;图2算法分类精度图样本数据投影前的数据分布情况;本发明提供的基于半监督流形学习的多光谱图像土地利用分类方法,包括以下步骤:
[0089] (1)读入多光谱遥感数据;
[0090] (2)读入选取的多光谱遥感图像中每一个样本数据点,根据其波段生成一个向量,从而将多光谱遥感图像中选取的样本用一个矩阵来表示,作为训练样本集;
[0091] (3)根据先验知识从训练样本集中选取部分样本数据进行已知地物类别的标注,生成样本类别标签;
[0092] (4)在部分数据类别信息已知的前提下,通过样本数据集构建相似图和相异图来度量数据点的相似性与相异性;
[0093] (5)利用热核构建来构建权重矩阵;
[0094] (6)计算局部相似结构矩阵和全局相异结构矩阵;
[0095] (7)通过目标优化函数计算投影矩阵;
[0096] (8)通过投影矩阵对整幅多光谱遥感数据进行投影;
[0097] (9)对投影后的整幅多光谱遥感数据用K-近邻分类算法进行分类;
[0098] (10)输出多光谱遥感数据地物分类图。
[0099] 作为上述实施例的进一步改进,所述步骤(2)中的训练样本集通过以下方式实现:
[0100] 根据不同波段对地物的光谱反射特性,把多光谱遥感图像数据变为三维数据X=T
{x1,x2,…,xM×N},其中M×N是图像的空间尺寸,B表示波段数,矩阵的每一列表示一个波段的数据值,每一行表示多光谱遥感图像数据中的一个点,如图2所示,基于半监督流形学习的多光谱遥感数据土地利用分类方法分类后算法的精度,取十次平均值为94.545%。
[0101] 作为上述实施例的进一步改进,所述步骤(3)中的样本类别标签通过以下方式实现:根据先验知识从训练样本集中选取部分样本数据进行已知地物类别的标注,生成样本类别标签:
[0102] X={(x1,l1),(x2,l2),…,(xC,lC),xC+1,xC+2,...,xM×N}T
[0103] 式中li为对数据点xi标注的地物类别标签,前C个点具有类别信息,其余的M×N-C个为无类别标签样本点。
[0104] 作为上述实施例的进一步改进,所述步骤(4)中的相似图和相异图包括以下步骤:
[0105] 1)将样本数据点分成相似图GS和相异图GD,在相似图GS和相异图GD的节点之间采用k最临近方法k-NN方法定义两个图中的边;
[0106] 2)每一点xj对xi,若xj∈knnS(xi),即每一点对(xi,xj),两者为近邻且相似度较高,则用一条边连接图GS中xi和xj两点;
[0107] 3)若xj∈knnD(xi),即相似度较低,则用一条边连接图GD中xi和xj两点;
[0108] 式中 为点xi的k个近邻点组成的子集,对每一点xi,近邻数
据集knn(xi)可以分为两部分:knnD(xi)与knnS(xi),其中:knnD(xi)表示近邻点来自于不同类样本点,knnS(xi)则为knn(xi)中其余部分,需要进一步说明的是,knnS(xi)中部分样本可能不包含类别信息,但是这些样本距离xi足够的近,knn(xi)同knnD(xi)和knnS(xi)之间的关系可用下式来表示:
[0109]
[0110] knnS(xi)=knn(xi)-knnD(xi)
[0111] 显然, knnS(xi)∪knnD(xi)=knn(xi);
[0112] 根据knnD(xi)和knnS(xi),就可以构建相似图GS和相异图GD。
[0113] 作为上述实施例的进一步改进,所述步骤(5)中权重矩阵通过以下方式实现:
[0114]
[0115] 权重矩阵为Wij=Ws,ij*Aij
[0116] 其中,Aij定义如下:
[0117]
[0118] 其中d2(xi,xj)是从Xi到Xj的距离,t为一常数。
[0119] 作为上述实施例的进一步改进,所述步骤(6)包括以下步骤:
[0120] (1)计算相似局部结构矩阵Srlw
[0121] 根据输入的多光谱样本数据集X和权重矩阵Wij,通过
[0122]
[0123]
[0124]
[0125]
[0126] 式中DS是一个对角阵,且Ds,ii=∑jWij,LS=DS-WS为Laplacian矩阵;
[0127] 然后,防止Slw为奇异矩阵采取如下方式,
[0128] Srlw=Slw+SE
[0129] 式中SE为β×Ε即一个较小的常数β乘以单位矩阵;
[0130] (2)计算全局相异结构矩阵Srlb
[0131] 根据输入的多光谱样本数据集X和权重矩阵Wij,通过
[0132] Srlb=Slb+ST
[0133] 其中ST,Slb定义如下:
T
[0134] ST=E{(X-μ)(X-μ)}
[0135]
[0136] 式中μi为每类的均值;μ为总样本的均值;n为每类含有的样本数;C是总类别数。
[0137] 作为上述实施例的进一步改进,所述步骤(7)通过以下方式实现:
[0138] 多光谱遥感图像数据投影后最理想的目标:使同类样本数据间散度尽可能小、不同类数据间散度尽可能大,基于上述考虑,投影矩阵A通过目标优化函数:
[0139] Maximize J(A)=tr(ATSrlbAT)
[0140] subject to tr(ATSrlwAT)=1
[0141] 通过采用Lagrange函数ATSrlbAT-λ(ATSrlwAT-1),上述最优化问题可以很容易地转换为广义的特征值求解问题,基于上述考虑,投影矩阵A通过目标优化函数:
[0142]
[0143] 求得,即A由以下特征方程的d个最大特征值λ1>λ2>…>λd对应特征向量v1,v1,…,vd构成
[0144] Srlbv=λSrlwv。
[0145] 作为上述实施例的进一步改进,所述步骤(8)通过以下方式来实现:
[0146] Y=ATZ;
[0147] 式中Y为整幅多光谱样本数据集Z通过投影矩阵A投影到嵌入空间的数据集。
[0148] 实施例2
[0149] 重庆市大渡口区多光谱遥感图像数据包含3个波段,分辨率为14.25m,主要包含建筑物、林地、草地、水域、耕地水田五类地物。本发明的实现流程如图1所示,具体实施方案按以下步骤进行:
[0150] (1)多光谱遥感数据读入:读入大渡口区多光谱遥感图像数据,波段分别为band1,band2至band3;
[0151] (2)选取多光谱遥感图像训练样本数据:每类选取50个已知类别信息的样本点,转化成300行3列的矩阵,另外在随机选取600个不知类别的样本点,表示为X=T
{x1,x2,...,x300,x301,...x900} ;
[0152] (3)用户利用先验知识标注部分样本:对每类地物标注50个数据点,即在样本训练数据X中前面300个点具有类别信息,余下600个样本点无类别标注;
[0153] (4)利用标注后的样本数据构建相似图GS和相异图GD:本发明首先找到每一个数据点的8个近邻点: 然后将knn(xi)分为两部分:knnD(xi)与
knnS(xi)。再在数据点和相似图GS和相异图GD的节点之间定义一一对应关系,然后根据k-NN方法定义两个图中的边,具体为:考虑每一点xj对xi,若xj∈knnS(xi),则用一条边连接图GS中xi和xj两点;若xj∈knnD(xi),则用一条边连接图GD中xi和xj两点;
[0154] (5)计算权重矩阵:
[0155] 多光谱数据的权重矩阵计算方法如下:
[0156]
[0157] 权重为Wij=Ws,ij×Aij
[0158] 其中,Aij表示xi和xj基于局部比例试探法的相似程度,定义如下:
[0159]
2
[0160] 这里的d(xi,xj)是从Xi到Xj的距离,t为一常数。在本发明中,t=1或25。
[0161] (6)计算局部相似结构矩阵和全局相异结构矩阵
[0162] 利用X和Wij计算局部相似结构矩阵Srlw和全局相异结构矩阵Srlb,具体如下:
[0163]
[0164]
[0165]
[0166]
[0167] 然后,防止Slw为奇异矩阵采取如下方式,Srlw=Slw+SE
[0168] 式中SE为β×Ε即一个较小的常数β乘以单位矩阵;
[0169] Srlb=Slb+ST 其中ST,Slb定义如下:
[0170] ST=E{(X-μ)(X-μ)T};
[0171] 式中μi为每类的均值;μ为总样本的均值;n为每类含有的样本数;C是总类别数。
[0172] (7)利用目标优化函数计算投影矩阵A
[0173] 投影矩阵A通过目标优化函数:
[0174] Maximize J(A)=tr(ATSrlbAT)
[0175] subjct to tr(ATSrlwAT)=1
[0176] 通过采用Lagrange函数ATSrlbA-λ(ATSrlwXTA-1),上述最优化问题可以很容易地转换为广义的特征值求解问题。基于上述考虑,投影矩阵A通过目标优化函数:
[0177]
[0178] 求得。利用Lagrange函数,上述最优化问题可以很容易地转换为广义的特征值求解问题
[0179] Srlbv=λSrlwv
[0180] 上述特征方程的d个最大特征值λ1>λ2>…>λd对应特征向量v1,v1,…,vd,构成投影矩阵A=(v1,v1,…,vd)。本发明中,根据以下方法确定d
[0181]
[0182] (8)利用投影矩阵A对整幅多光谱遥感图像数据Z进行投影:Y=ATZ;
[0183] (9)对投影后的整幅多光谱遥感数据用K-近邻分类算法进行分类;如图3和图
4所示,其中图3为训练样本集投影前二维分布图样本数据经过投影矩阵投影后的分布情况;图4为训练样本集投影后二维分布图算法实现分类后输出效果图;其中将河流湖泊(水域)、森林用地、建筑用地、水田耕地和绿地分类后如图所示。
[0184] (10)输出多光谱遥感数据地物分类图,至此,本发明完成了多光谱遥感数据土地利用分类的全过程。效果图见附图5至图11所示,图5红色代表河流湖泊(水域)效果图;
图6绿色代表森林用地效果图;图7蓝色代表建筑用地效果图;图8紫色代表水田耕地效果图;图9白色代表绿地效果图;图10黑色代表边界效果图;图11整体效果图。
[0185] 以上所述仅为本发明的优选实施例,并不用于限制本发明,显然,本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样,倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内,则本发明也意图包含这些改动和变型在内。