本发明公开了一种系统性的高光谱草原群落划分方法，涉及遥感图像应用领域，包括：步骤1、结合多种数据选取样本标签数据；步骤2、利用线性预测方法选取适当的波段子集并构建光谱特征；步骤3、将选择后的波段利用扩展形态学方法提取空间特征；步骤4、将光谱特征与空间特征进行堆栈，利用主成分分析在低维空间中融合，并利用随机森林(RF)得到预分类结果；步骤5、利用标签相似概率过滤器对分类结果进行分类后处理，进一步利用RF进行分类，得到最终草原分类结果。本发明充分利用高光谱图像的空‑谱特征，解决了特征提取中的“维度灾难”问题，并通过分类后处理的方法缓解了噪声对分类结果的影响，从而提高草原群落识别精度。
1.一种系统性的高光谱草原群落划分方法，其特征在于，包括以下步骤：步骤S1、样本标签数据的选择在一定程度上决定了草原群落划分的准确性，因此为了确保样本选取的合理性，以内蒙古锡林郭勒盟为研究区，结合实测样本点数据、中华人民共和国植被图(1：1000000)数据，以下简称：植被图数据、国家级自然保护区功能区划数据，以下简称：保护区数据和ZY1‑02D高光谱数据共同构建样本标签数据，具体过程：步骤S1.1、在草原植被生长茂盛的7、8月份进行野外实地调查，并记录调查点的经纬度坐标和草原优势植被类型，继而选取云量小于20％且成像时间为7、8月份的4幅内蒙古锡林郭勒盟ZY1‑02D高光谱影像，在实测点周围做200m*200m的缓冲区，并将缓冲区内的像素作为样本标签数据；
步骤S1.2、植被图数据是2008年由地质出版社出版，包含我国11个植被类型组，55个植被型的900多个群系和亚群系以及约2000多个群落优势种的地理空间分布，选取覆盖锡林郭勒盟的植被图数据；保护区数据是国家对有代表性的自然生态系统、珍稀濒危野生动植物物种的天然集中分布区，包含核心区、缓冲区和实验区，核心区是生态系统保护的核心区域，禁止任何单位和个人进入，也不允许从事科学研究活动，核心区外围划定一定面积的缓冲区，只准从事科学研究观测活动，缓冲区外围划为实验区，从事科学试验、教学实习、参观考察、旅游以及繁殖珍稀、濒危野生动植物活动，选取锡林郭勒盟保护区数据；将植被图数据、保护区数据和影像数据在ArcGIS中进行叠加，核心区和缓冲区范围内的群落认为在多年内没有发生改变，但实验区内的群落需要结合当前影像进行样本标签选择；
步骤S1.3、在植被图数据、保护区数据和影像数据三者叠加的前提条件下，对于特殊的区域，包括建成区，即市行政区范围内经过征收土地和实际建设发展起来的非农业生产建设地段，以市区为主并参考植被图数据，认为建成区适宜生长的群落为市区所对应的群落；
步骤S1.4、通过步骤S1.1‑1.3，共选取了10种草原群落类型，分别为大针茅草原(13655)、芨芨草盐生草甸(7775)、农作物(1384)、丛生禾草(3822)、大果榆疏林(2249)、糙隐子草(2137)、锦鸡儿灌丛(1266)、贝加尔针茅(1168)、杂类草草甸草原(1354)和丛生隐子草(1146)，其中，括号内为选取的高光谱影像像素数目，共计35956；
步骤S2、利用线性预测LP方法从高光谱数据中选取适当的波段子集，组成高光谱数据的光谱特征；
步骤S2.1、从原始高光谱数据中随机选择一组波段B1和B2初始化算法，并生成波段子集Φ＝{B1B2}；
步骤S2.2、估计和线性预测的方法能够评估波段之间的相似性，因此能够从所有的波段中选出与当前Φ最不同的波段B3，并且更新波段子集Φ＝Φ∪{B3}，a0+a1B1+a2B2＝B′ (1)
其中，B′是B3的估计和线性预测结果，a0，a1，a2是使线性误差e＝||B3‑B′||最小的参数；
步骤S2.3、重复步骤S2.2，直到Φ达到预先设定的波段子集数，并得到光谱特征y ＝Φ；
步骤S3、将步骤S2中光谱特征Φ利用扩展形态学方法EMP提取空间特征，降低噪声对高光谱图像分类的影响，包括开变换和闭变换；
步骤S3.1、对高光谱像素x执行先腐蚀后膨胀开变换，操作定义：
其中，MPγ表示经过开变换的影像；γ (x)表示对影像x进行的开操作；l为开运算的次数；
步骤S3.2、对高光谱像素x执行先膨胀后腐蚀闭变换，操作定义：
其中，MPΦ表示经过闭变换的影像；Φ (x)表示对影像x进行的闭操作；p为闭运算的次数；公式(2)和(3)中k为结构元素大小，当k＝0时，表示未对原始图像进行操作；
步骤S3.3、将原始高光谱像素与开、闭操作后的结果进行合并，则得到EMP结果，也即空spa间特征y
y ＝EMP＝{MPγ，x，MPΦ } (4)利用扩展形态学EMP特征提取方法能够获取丰富的空间信息，并且能够对结构性物体或物质进行区分；
步骤S4、将高光谱图像的光谱特征和空间特征进行堆栈，以获取高光谱图像空‑谱特征，为了避免高维特性导致的“维度灾难”问题，利用主成分分析法PCA在低维空间中进行融合；
步骤S4.1、将步骤S2和S3中的光谱特征和空间特征进行堆栈并归一化得到融合特征步骤S4.2、对数据 利用PCA从n维空间投影到m维空间：其中，P＝[P1，...，Pm]，m＜n，表示根据特征值排列顺序选择前m个特征值即m维空间对应的特征向量组成的变换矩阵；
步骤S4.3、将得到的空‑谱融合特征G和样本标签数据输入到RF分类器中进行训练，利用训练得到的模型对影像进行测试得到草原群落预分类结果；
步骤S5、利用标签相似概率过滤器LSPF方法对分类结果做进一步分类后处理，利用与周围像素的空间相关性估计所有类别中所有像素的标签概率，然后将概率数据和样本标签数据输入到分类器中，以细化预先获得的分类图，减少椒盐噪声对分类结果的影响；
步骤S5.1、根据像素是否属于类别c，将初始分类结果映射转换为一系列的二值标签图，属于该类则设置为1，否则为0；影像场景中有10类地物类型，则预分类结果图分为10层二值标签图；
步骤S5.2、依据现实世界中二维高斯分布的简单性和有效性，在局部空间窗口中利用指数函数建模LSP：其中，g(x，y)为同一空间窗口中的中心像素与其相邻像素之间的LSP，e为自然常数，σ为标准偏差；
步骤S5.3、进一步的在每个空间窗口内，对于每个类c的LSPF表示为：其中，fc(i，j)为中心像素(i，j)类别c的累积概率，⊙为哈达玛积，Gij＝[g(x，y)]为在(i，j)_th空间窗口中的LSP矩阵， 为预分类结果图中的布尔矩阵；
步骤S5 .4、进一步的遍历所有类，得到所有像素的累积概率数据为
表示为一个1×C的向量，继而将所有像素的累积概率和样本标签数
据输入到RF分类器中，利用训练得到的模型对Fc(i，j)进行测试得到草原群落高光谱数据的分类图。
一种系统性的高光谱草原群落划分方法
技术领域
[0001] 本发明涉及遥感图像应用领域，具体的说是一种系统性的高光谱草原群落划分方法。
背景技术
[0002] 草原退化的主要特征是植被群落结构的改变，基于高光谱的草原群落识别与划分是利用遥感进行大面积高精度草原退化监测与治理的基础与前提。混杂度高且植被低矮的草原，多光谱遥感识别与划分难度较大，而高光谱遥感图像中的每个像元由几十个甚至几百个窄波段组成，通过分类操作可用于精细的草原群落划分。
[0003] 高光谱图像包含丰富的空间信息和光谱信息，有效地特征提取是实现高光谱图像精细分类的关键步骤，否则不仅会浪费大量的计算资源，而且具有较低的分类精度。传统的方法，一般是只利用光谱特征或者空间特征，但这会忽略对应的空间特征或者光谱特征对分类的贡献，此外将空‑谱特征进行堆栈融合的方法会遇到高光谱数据“维度灾难”的问题。近年来，以卷积神经网络为代表的深度学习方法以其端到端的训练方式被广泛的应用到高光谱图像分类中，但因其缺乏可解释性，以及需要耗费大量的人力、物力进行调参，限制了其在实际生产中的应用。因此本发明从样本标签数据采集到特征提取、融合，再到分类结果后处理，提供了一个系统性的操作流程，可在特定场景中进行应用实践。
发明内容
[0004] 本发明针对高光谱草原群落划分中传统方法特征利用不足以及深度学习方法训练时间长、调参困难等问题，提出一种系统性的高光谱草原群落划分方法，用来提高特征利用率和分类精度。
[0005] 这种系统性的高光谱草原群落划分方法，包括以下步骤：
[0006] 步骤S1、利用实测样本点数据、中华人民共和国植被图(1：1000000)数据(以下简称：植被图数据)和国家级自然保护区功能区划数据(以下简称：保护区数据)，选取研究区范围内部分像素作为样本标签数据；
[0007] 步骤S2、输入高光谱影像，利用基于线性预测(LP)的方法选取适当的波段子集并提取光谱特征；
[0008] 步骤S3、将选择后的波段利用扩展形态学方法(EMP)得到高光谱图像的空间特征；
[0009] 步骤S4、将提取的光谱特征与空间特征进行堆栈，利用主成分分析法(PCA)在低维空间中进行融合，结合样本标签数据并利用随机森林(RF)得到预分类结果；
[0010] 步骤S5、利用标签相似概率过滤器(LSPF)方法对分类结果进行分类后处理，进一步利用RF进行分类，减少椒盐噪声对分类结果的影响。
[0011] 本发明的特点还在于：
[0012] 执行步骤S1时，主要是通过实测样本点数据、植被图数据和保护区数据，并结合影像数据选取样本标签数据，具体过程：
[0013] 步骤S1.1、实测样本点数据是人工野外实地考察得到的数据，包括经纬度坐标和1m*1m格网内的草原优势种信息，然后选取与实测样本点数据采集时间相近的ZY1‑02D高光谱影像，继而依据实测的样本点数据做200m*200m的缓冲区，并将缓冲区内的像素作为样本标签数据；
[0014] 步骤S1.2、植被图数据是2008年由地质出版社出版，包含我国11个植被类型组，55个植被型的900多个群系和亚群系以及约2000多个群落优势种的地理空间分布；保护区数据是国家对有代表性的自然生态系统、珍稀濒危野生动植物物种的天然集中分布区，包含核心区、缓冲区和实验区，核心区是生态系统保护的核心区域，禁止任何单位和个人进入，也不允许从事科学研究活动，核心区外围可以划定一定面积的缓冲区，只准从事科学研究观测活动，缓冲区外围划为实验区，可以从事科学试验、教学实习、参观考察、旅游以及繁殖珍稀、濒危野生动植物等活动；将植被图数据、保护区数据和影像数据在ArcGIS中进行叠加，核心区和缓冲区范围内的群落认为在多年内没有发生改变，但实验区内的群落需要结合当前影像进行样本标签选择；
[0015] 步骤S1.3、在植被图数据、保护区数据和影像数据三者叠加的前提条件下，对于特殊的区域，如建成区(市行政区范围内经过征收土地和实际建设发展起来的非农业生产建设地段)，以市区为主并参考植被图数据，认为建成区适宜生长的群落为市区所对应的群落。
[0016] 执行步骤S2时，主要包括对高光谱影像进行波段选择并提取光谱特征，具体过程：
[0017] 步骤S2.1、通过选择高光谱图像中的一组波段B1和B2初始化算法，然后生成波段子集Φ＝{B1B2}；
[0018] 步骤S2.2、通过估计和线性预测的方法从所有的波段中选出与当前Φ最不同的B3，并且更新波段子集Φ＝ΦU{B3}；
[0019] a0+a1B1+a2B2＝B′ (1)
[0020] 其中，B′是B3的估计和线性预测结果，a0，a1，a2是使线性误差e＝||B3‑B′||最小的参数；
[0021] 步骤S2.3、重复步骤S2.2，直到Φ达到预先设定的波段子集数，并得到光谱特征spey ＝Φ。
[0022] 执行步骤S3时，主要包括对波段选择后的高光谱影像进行开变换和闭变换，具体过程为：
[0023] 步骤S3.1、对高光谱像素x执行先腐蚀后膨胀开变换，操作定义：
[0024]
[0025] 其中，MPγ表示经过开变换的影像；γk(x)表示对影像x进行的开操作；l为开运算的次数；
[0026] 步骤S3.2、对高光谱像素x执行先膨胀后腐蚀闭变换，操作定义：
[0027]
[0028] 其中，MPΦ表示经过闭变换的影像；Φk(x)表示对影像x进行的闭操作；p为闭运算的次数；公式(2)和(3)中k为结构元素大小，当k＝0时，表示未对原始图像进行操作；
[0029] 步骤S3.3、将原始高光谱像素与开、闭操作后的结果进行合并，则得到EMP结果，也spa即空间特征y
[0030] yspa＝EMP＝{MPγ，x，MPΦ} (4)
[0031] 执行步骤S4时，主要包括对高光谱影像的光谱特征和空间特征进行堆栈，并利用PCA在低维空间中融合，以利用RF进行分类，具体过程：
[0032] 步骤S4.1、将步骤S2和S3中的光谱特征和空间特征进行堆栈并归一化得到融合特征
[0033] 步骤S4.2、对数据 利用PCA从n维空间投影到m维空间：
[0034]
[0035] 其中，P＝[P1，...，Pm](m＜n)表示根据特征值排列顺序选择前m个特征值(即m维空间)对应的特征向量组成的变换矩阵。
[0036] 步骤S4.3、将得到的空‑谱融合特征G和样本标签数据输入到RF分类器中进行训练，利用训练得到的模型对影像进行测试得到草原群落预分类结果。
[0037] 执行步骤S5时，对得到的预分类结果利用LSPF(Label Similarity Probability Filter)进行分类后处理，继而再利用RF分类器得到草原群落分类结果，以解决光谱相似导致的椒盐噪声。
[0038] 步骤S5.1、根据结果像素是否属于类别c，将预分类结果映射转换为一系列的二值标签图，属于类别c则设置为1，否则为0，例如，如果影像场景中有5类地物类型，则预分类结果图分为5层二值标签图。
[0039] 步骤S5.2、在每个空间窗口内，像素(i，j)为类别c的累积概率被计算为fc(i，j)，然后遍历所有类，得到所有像素的累积概率数据为 可以表示为一个1×C的向量，继而将所有像素的累积概率和样本标签数据输入到RF分类器中，利用训练得到的模型对Fc(i，j)进行测试得到草原群落高光谱数据的分类结果图。
[0040] 本发明针对高光谱遥感图像提取的光谱特征和空间特征，利用PCA的方法在低维空间中融合，以此挖掘高光谱图像的空‑谱特征，并缓解高光谱图像“维度灾难”问题，随后利用RF对融合特征进行分类得到预分类结果，继而利用LSPF进行分类后处理，以解决分类结果中椒盐噪声现象，最后利用RF得到最终分类结果。本发明能够有效地弥补光谱特征或者空间特征带来的分类性能差、精度低等问题，为草原群落识别提供了一套完整的分类流程。
附图说明
[0041] 附图1为本发明的实施例流程图。
具体实施方式
[0042] 为使本发明的技术方案、解决的技术问题更加清楚明白，以下结合具体实施例图1，对本发明做进一步的描述。
[0043] 本发明提供了一种系统性的高光谱草原群落划分方法，草原作为生态系统的重要组成部分，在保护水土、防风固沙及保护生物多样性方面发挥着重要的作用，但近年来由于气候变化和人类活动的影响，中国已经出现大面积的草原退化，因此草原退化的监测显的尤为重要。利用高光谱遥感进行草原退化监测得到了广泛的应用，但由于其特征提取不足以及“维度灾难”的问题，限制了进一步的应用，同时为大面积的草原退化监测带来了挑战。
[0044] 步骤S1、样本标签数据的选择在一定程度上决定了草原群落划分的准确性，因此为了确保样本选取的合理性，本发明以内蒙古锡林郭勒盟为研究区，结合实测样本点数据、中华人民共和国植被图(1：1000000)数据(以下简称：植被图数据)、国家级自然保护区功能区划数据(以下简称：保护区数据)和ZY1‑02D高光谱数据共同构建样本标签数据，具体过程：
[0045] 步骤S1.1、在草原植被生长茂盛的7、8月份进行野外实地调查，并记录调查点的经纬度坐标和草原优势植被类型，继而选取云量小于20％且成像时间为7、8月份的4幅内蒙古锡林郭勒盟ZY1‑02D高光谱影像，在实测点周围做200m*200m的缓冲区，并将缓冲区内的像素作为样本标签数据；
[0046] 步骤S1.2、植被图数据是2008年由地质出版社出版，包含我国11个植被类型组，55个植被型的900多个群系和亚群系以及约2000多个群落优势种的地理空间分布，选取覆盖锡林郭勒盟的植被图数据；保护区数据是国家对有代表性的自然生态系统、珍稀濒危野生动植物物种的天然集中分布区，包含核心区、缓冲区和实验区，核心区是生态系统保护的核心区域，禁止任何单位和个人进入，也不允许从事科学研究活动，核心区外围可以划定一定面积的缓冲区，只准从事科学研究观测活动，缓冲区外围划为实验区，可以从事科学试验、教学实习、参观考察、旅游以及繁殖珍稀、濒危野生动植物等活动，选取锡林郭勒盟保护区数据；将植被图数据、保护区数据和影像数据在ArcGIS中进行叠加，核心区和缓冲区范围内的群落认为在多年内没有发生改变，但实验区内的群落需要结合当前影像进行样本标签选择；
[0047] 步骤S1.3、在植被图数据、保护区数据和影像数据三者叠加的前提条件下，对于特殊的区域，如建成区(市行政区范围内经过征收土地和实际建设发展起来的非农业生产建设地段)，以市区为主并参考植被图数据，认为建成区适宜生长的群落为市区所对应的群落；
[0048] 步骤S1.4、通过步骤S1.1‑1.3，共选取了10种草原群落类型，分别为大针茅草原(13655)、芨芨草盐生草甸(7775)、农作物(1384)、丛生禾草(3822)、大果榆疏林(2249)、糙隐子草(2137)、锦鸡儿灌丛(1266)、贝加尔针茅(1168)、杂类草草甸草原(1354)和丛生隐子草(1146)，其中括号内为选取的高光谱影像像素数目，共计35956。
[0049] 步骤S2、利用线性预测(LP)方法从高光谱数据中选取适当的波段子集，组成高光谱数据的光谱特征。
[0050] 步骤S2.1、从原始高光谱数据中随机选择一组波段B1和B2初始化算法，并生成波段子集Φ＝{B1B2}；
[0051] 步骤S2.2、估计和线性预测的方法能够评估波段之间的相似性，因此能够从所有的波段中选出与当前Φ最不同的波段B3，并且更新波段子集Φ＝ΦU{B3}，[0052] a0+a1B1+a2B2＝B′ (1)
[0053] 其中，B′是B3的估计和线性预测结果，a0，a1，a2是使线性误差e＝||B3‑B′||最小的参数；
[0054] 步骤S2.3、重复步骤S2.2，直到Φ中波段子集数为20，并令yspe＝Φ得到光谱特征，此时能够尽可能的通过最少的波段数表示高光谱数据的全部信息，减少数据冗余，便于空间特征的提取。
[0055] 步骤S3、将步骤S2中光谱特征Φ利用扩展形态学方法(EMP)提取空间特征，降低噪声对高光谱图像分类的影响，其主要包括开变换和闭变换。
[0056] 步骤S3.1、对高光谱像素x执行先腐蚀后膨胀开变换，操作定义为：
[0057]
[0058] 其中，MPγ表示经过开变换的影像；γk(x)表示对影像x进行的开操作；l＝5开运算的次数。
[0059] 步骤S3.2、对高光谱像素x执行先膨胀后腐蚀闭变换，操作定义为：
[0060]
[0061] 其中，MPΦ表示经过闭变换的影像；Φk(x)表示对影像x进行的闭操作；p＝5为闭运算的次数；公式(2)和(3)中k为结构元素大小，当k＝0时，表示未对原始图像进行操作。
[0062] 步骤S3.3、基于“开变换”和“闭变换”后的形态特征，选定对应光谱波段构建形态spa学轮廓，此时，像素x的形态学轮廓为m×(2l+1)＝220，其中m为光谱波段数，即y ＝EMP[0063] EMP＝{MPγ，x，MPΦ} (4)
[0064] 利用扩展形态学(EMP)特征提取方法能够获取丰富的空间信息，并且能够对结构性物体或物质进行区分。
[0065] 步骤S4、将高光谱图像的光谱特征和空间特征进行堆栈，以获取高光谱图像空‑谱特征，为了避免高维特性导致的“维度灾难”问题，本发明利用主成分分析法(PCA)在低维空间中进行融合。
[0066] 步骤S4.1、将步骤S2和S3中的光谱特征和空间特征进行堆栈并归一化得到融合特征 归一化操作能够减小数据间的差异。
[0067] 步骤S4.2、对数据 利用PCA从n＝220维空间投影到m＝50维空间(m＜n)：
[0068]
[0069] 其中，P＝[P1，...，Pm](m＜n)表示根据数据特征值排列顺序选择前m＝50个特征值(即m维空间)对应的特征向量组成的变换矩阵。
[0070] 步骤S4.3、将得到的空‑谱融合特征G和样本标签数据输入到RF分类器中进行训练，利用训练得到的模型对影像进行测试得到草原群落预分类结果。由于RF是基于像素特征进行分类，因此会产生椒盐噪声，导致草原分类不准确，因此有必要对预分类结果进行分类后处理操作。
[0071] 步骤S5、利用标签相似概率过滤器(LSPF)方法对分类结果做进一步分类后处理，利用与周围像素的空间相关性估计所有类别中所有像素的标签概率，然后将概率数据和样本标签数据输入到分类器中，以细化预先获得的分类图，减少椒盐噪声对分类结果的影响。
[0072] 步骤S5.1、根据像素是否属于类别c，将初始分类结果映射转换为一系列的二值标签图，属于该类则设置为1，否则为0，例如，影像场景中有10类地物类型，则预分类结果图分为10层二值标签图。
[0073] 步骤S5.2、依据现实世界中二维高斯分布的简单性和有效性，在局部空间窗口中利用指数函数建模LSP(Label Similarity Probability)：
[0074]
[0075] 其中，g(x，y)为同一空间窗口中的中心像素(值为1)与其相邻像素之间的LSP，e为自然常数，σ为标准偏差。
[0076] 步骤S5.3、进一步的在每个空间窗口内(5×5)，对于每个类c的LSPF(Label Similarity Probability Filter)可以表示为：
[0077]
[0078] 其中，fc(i，j)为中心像素(i，j)类别c的累积概率，⊙为哈达玛积(Hadamard product)，Gij＝[g(x，y)]为在(i，j)_th空间窗口中的LSP矩阵， 为预分类结果图中的布尔矩阵(Boolean matrix)。
[0079] 步骤S5.4、进一步的遍历所有类，得到所有像素的累积概率数据为可以表示为一个1×C的向量，继而将所有像素的累积概率和样本标签数据输入到RF分类器中，利用训练得到的模型对Fc(i，j)进行测试得到草原群落高光谱数据的分类图。
[0080] 综上可知，采用本发明系统性的高光谱草原群落划分方法，能够充分的挖掘高光谱数据的空‑谱特征，解决深度学习方法模型调参难、训练代价高等问题，在此基础上提高草原群落的识别精度，并且能够完全通过计算机软件一步实现，实用性较强。
[0081] 本发明着重强调了一种系统性的高光谱草原群落划分方法，也可以被认为是针对特定场景下一套完整的操作流程，包括样本合理性选取、光谱特征提取、空间特征提取、空‑谱特征融合以及分类后处理。基于本发明的上述具体实施例，本技术领域的技术人员在不脱离本发明原理的前提下，对本发明所作出的任何改进和修饰(包括更换应用场景)，皆应纳入本发明的保护范围。
