一种基于多因素的海岸线变化识别方法与流程

本发明归属于海岸线变化识别技术领域，具体涉及一种基于多因素的海岸线变化识别方法。

背景技术：

海岸是人类社会经济活动最活跃、最集中的地区，与海运交通、领土主权、工程建设、资源开发、空间利用等诸多方面密切相关。海岸天然的地理优势，孕育了一个又一个自“大航海”时代以来的世界经济中心。进入21世纪，沿海城市在资源、环境、交通上优势更加明显，世界沿海国家经济重心向滨海地区转移；新时代我国加快推进海洋强国建设，海岸作为探索海洋的出发点与归宿点，正确理解、合理认知、科学探究海岸空间为实现海洋强国战略目标提供重要保障。

海岸线是海洋和陆地的分界线，是海岸乃至地球上最重要的边界线。全球海岸线总长44万千米，分布在100多个国家和南北极地区，全球已有超过一半的人口居住在距海岸线100千米内的沿海地区；我国拥有超过1.8万千米的大陆岸线和1.4万千米的岛屿岸线，沿岸滩涂超过2万平方千米，探明海滨矿藏储量15.3亿吨，沿海地区集中了40％的人口和60％的产值。海岸线的位置、走向和形态不仅体现了全球及海岸带环境过程，更是社会经济活动、人类综合作用的结果与反映。此外，海岸线作为陆海交界，除受到一般陆地和地表过程影响外，还受到潮汐、波浪、环流、信风等海洋动力作用以及生物活动影响，具有极大的复杂性和时空不确定性。

技术实现要素：

本发明提供一种基于多因素的海岸线变化识别方法，提出了用于海岸区人类影响分类的模糊模型，其集成了海岸线变化、归一化差分植被指数和沉降影响，以通过地理信息系统(gis)的图形可视化过程来增强数值语言模糊分类。显示了局部动力学的复杂行为，从而为环境问题和海岸带综合管理增加了有用的、实质性的信息。

一种基于多因素的海岸线变化识别方法，包括如下步骤：

s1.输入海岸线基线，采用遥感卫星图像覆盖所述海岸线区域的面积延伸；

s2.将来自每个所述遥感卫星场景的海岸线分割为多个长条多边形形状；

s3.采用模糊模型使用所述遥感卫星图像来提取关于时间变化的信息：(i)海岸线变化；(ii)归一化差分植被指数；和(iii)沉降演变；

s4.分别计算长条多边形的侵蚀、堆积和稳定的面积，所述侵蚀、堆积和稳定的面积被计算为相对于每个多边形区总面积的百分比，并将百分比用作第一输入(海岸线变化)的三个变量(x1、x2和x3)；

s5.计算归一化差分植被指数(ndvi)作为第四变量x4，ndvi＝(nir-red)/(nir+red)，nir和red分别表示近红外和红色反射率；

s6.分别建立针对每个多边形区和时间图像分类区域的建筑物累积变量及沉降演变，作为海岸线模糊分类模型的第五变量x5和第六变量y；

s7.采用模糊模型对每个变量x1、x2、x3、x4和y进行计算，

s8.以交互形式对每个变量x1、x2、x3、x4和y进行验证及微调，直到满足验证标准；

s9.对验证后的变量x1、x2、x3、x4和y值采用特定的语言进行标记：低、中、高，根据标记结果确定海岸带的人为影响分类等级。

进一步地，定义快照质量函数(snapshotquality)和历史代价函数(historycost)，快照质量函数衡量的是当前时间段变量x1、x2、x3、x4和y的聚类划分质量，表示聚类结果和当前时间段的数据的匹配程度；历史代价函数衡量的是相邻时间段变量x1、x2、x3、x4和y的聚类划分是否具有较好的时序平滑特性，表示的是历史模型与当前时间段的聚类结果的差异；对于包含t个时间序列数据s＝{x1,x2,x3,x4,x5}，且每个时间点的数据集x包含n个样本，xτ＝{x1,x2,x3,x4,x5}。定义如下的快照质量函数sq(cτ,mτ)以及历史代价函数hc(cτ,cτ-1)：

其中，cτ表示τ时间段数据的聚类结果，mτ是时间段数据的相似度矩阵。

进一步地，步骤s8中对每个变量x1、x2、x3、x4和y的验证方法包括如下步骤：

s11.对τ时间段的数据进行聚类时，将目标函数定义为：

将函数sq(cτ,mτ)及hc(cτ,cτ-1)带入目标函数j得到全局目标函数：

式中，c、n分别是聚类数和样本数，表示τ时间段的数据与τ时间段第i个分段(第i个聚类)的距离，表示τ时间段数据属于第i个聚类的隶属度，且满足每个样本属于所有类的隶属程度之和等于1，cτ指τ时间段的隶属矩阵，cτ∈{0,1}^c*k，而ci,k＝1表示τ时间段数据属于第i类，各数据在各时间段内只能属于某个确定的类，即∑τci,k＝1；

s12.最小化全局目标函数以求得最优分割参数，通过拉格朗日定理求解，定义下式为拉格朗日目标函数：

s13.通过求解对x，μ，λ的偏导数，并令这些偏导数分别为零。此时便可获得最优聚类中心及各个数据隶属于各类的隶属度；

s14.将数据的时间坐标作为一个额外变量并使其参与聚类过程，迭代计算方式如下：

计算τ时间段c个聚类中心vi：

计算聚类中心的模糊协方差矩阵fi：

计算距离函数值d²(xk,vi)：

利用距离函数值d²(xk,vi)更新模糊划分矩阵u：

当满足条件||u^(l+1)-u^(l)<ε||时，中止聚类算法，否则增加迭代次数，使l＝l+1，转到步骤s11，重复上述步骤，直至满足该条件为止。

进一步地，将所述遥感卫星图像进行图像预处理，进行大气校正，进而获得每个图像元素校正后的全光谱信息，且所述遥感卫星图像在所述海岸线区域上具有小于10％的云量覆盖。

进一步地，利用大气校正模型去除天空光和大气散射效应，使遥感图像更准确地反映特征的光谱值，大气校正模型由下式获得：

ρ为地表反射率，ρτ(θs,θv,φs,φv)为大气上界反射率，θs是太阳天顶角，θs是太阳方位角，θv是传感器方位角，t(θs)是大气透射率，tg(θs,θv)是太阳-目标大气路径透射率，t(θv)是目标-传感器大气路径透射率，ρr+α是分子散射和气溶胶散射所构成的路径辐射反射率，s是大气半球反射率。

进一步地，所述图像预处理采用暗对象减影方法，将由于太阳照度、传感器观察几何形状和季节性变化而引起的卫星图像辐射度差异进行校正和归一化。

进一步地，将所述遥感卫星图像进行图像预处理，进行几何校正，具体的操作步骤是地面控制点选择，像素坐标变换和像素亮度值的重新采样，像元坐标变换采用的是基于坐标变换的的方法，变换的公式为

x’＝a0+a1x+a2y+a3x²+a4xy+a5y²

y'＝b0+b1x+b2y+b3x²+b4xy+b5y²

其中，x’，y’为校正后输出的控制点坐标，x，y为地面控制点在原图像中的坐标，rmserror为每个控制点的均方根误差。

进一步地，模糊逻辑可以根据模型设计分类为低、低/中、中/高或高。

由于受潮汐、海岸地形等因素的影响较大，卫星过境时水边线刚好位于平均大潮高潮线的影像是很难获取的，海岸线自动解译算法的研宄大多是基于遥感影像中水边线的提取，即卫星过境时刻所记录的海陆分界线。本发明提出了用于海岸区人类影响分类的模糊模型，其集成了海岸线变化、归一化差分植被指数和沉降影响，以通过地理信息系统(gis)的图形可视化过程来增强数值语言模糊分类。诊断全世界沿海地区的人为影响(即与人类活动有关的影响)，例如沿海开发和规划、过度捕捞等。沿海环境的分类、风险和脆弱性评估是必不可少的环节，因此，本发明的算法最终显示了局部动力学的复杂行为，从而为环境问题和海岸带综合管理增加了有用的、实质性的信息。

本发明首次提出将海岸线变化(侵蚀、增长、稳定性)；ndvi(归一化差分植被指数)；和沉降演变作为随时间的变量影响因素，通过结合多个变量的替代方法，从而对与复杂环境系统相关的问题进行建模，并消除不精确和主观的概念，充分挖掘了利用地理信息中的遥感数据的潜力，并通过包括社会经济数据在内的其他影响因素以增强数值模糊分类的准确率。通过生成的空间地图实现图形可视化，这有利于地理特征的映射增强对每个部门演化模式识别的区分。

本发明的模型输出表示分数是从0到1的数字，其可转换为模糊语言分类变量；即低、中和高。同时，通过ndvi(归一化差分植被指数)使用gis通过图形增强了可视化程度。

本发明在对当前时间段的数据进行聚类时会考虑上一时间段的聚类结果。通过本发明的算法，可以有效地处理分析具有非线性及强耦合特性的海岸带数据，使得改进算法能够实时分析线性时间序列数据，加入了快照质量函数和历史代价函数，增强了聚类结果和当前时间段数据的匹配程度，且使得相邻时间段的聚类划分具有良好的时序平滑特性。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍。

图1为本发明海岸线变化识别方法的流程图。

具体实施方式

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图示出的结构获得其他的附图。

沿海岸带探测人类影响的方法包括海岸线侵蚀/增生模式评估，具体地，通过(a)地形剖面分析(考虑到跨岸形态以及作用在海滩上的破坏力和建设力之间的平衡)；(b)海岸线变化率；(c)土地利用/土地覆盖监测结合利用地理信息系统可视化的海岸线变化率。

本发明结合多个变量的替代方法，从而对与复杂环境系统相关的问题进行建模，并消除不精确和主观的概念，充分挖掘了利用地理信息中的遥感数据的潜力，通过包括社会经济数据在内的其他影响因素以增强数值模糊分类的准确率。通过生成的空间地图实现图形可视化，这有利于地理特征的映射增强对每个部门演化模式识别的区分。

一方面，模型输出表示分数是从0到1的数字，其可转换为模糊语言分类变量；即低、中和高。另一方面，通过ndvi(归一化差分植被指数)使用gis通过图形增强了可视化程度。

从遥感数据、海岸线变化、ndvi计算和结算影响中提取海岸线位置。模糊模型被设计为具有五个变量，即侵蚀、增长、稳定性、ndvi和累积。所有的语言学标签(模糊集)、隶属函数、模糊规则和去模糊化提供了输出，该输出是代表海岸带人为影响分类的明确数字。

第一步：数据处理

输入基线，使用遥感卫星图像来覆盖研究区域的面积延伸。卫星数据集是考虑到每年的相同/最近月份(例如八月和九月)而选择的，以寻求通过最小化季节变化来增加土地利用类别的分离。此外，所有选择的图像在研究区域上应当具有小于10％的云量覆盖。

由于遥感数据受到大气影响等诸多因素的影响，因此这些数据集无法用于进一步分析。卫星图像只能在执行若干图像预处理步骤(包括大气校正、几何校正)之后才能使用，以消除或最小化大气影响，进而获得每个图像元素(像素)的校正后的全光谱信息。本发明的大气校正方法可将由于太阳照度、传感器观察几何形状和季节性变化而引起的卫星图像辐射度差异进行校正和归一化。

作为输入，考虑三个方面，模糊模型使用前面描述的卫星图像来提取关于时间变化的信息：(i)海岸线变化；(ii)ndvi(归一化差分植被指数)；和(iii)沉降演变。

以海岸岛为例，由于研究区域在所有侧面(岛)都被水包围，因此在类似的缩放水平(1:5000的均匀比例)下，使用屏幕上手动数字化技术，将来自每个卫星场景的海岸线分割为多边形形状。这样可以有效描绘海岸线。

计算侵蚀和堆积面积(长条多边形)，对于每两个连续的多边形，使用gis环境中的空间联合工具，侵蚀、堆积和稳定性的面积被计算为相对于每个多边形区的总面积的百分比，并且那些值(％)被用作第一输入(海岸线变化)的三个变量(x1、x2和x3)。

ndvi(归一化差分植被指数)模糊模型中使用的第二个输入数据集是ndvi，该指数由图像中每个像素的新计算值组成，范围从-1到+1。ndvi由方程(1)和两个所需的输入波段，即近红外(nir)和红色(red)反射率计算。

ndvi＝(nir-red)/(nir+red)

在执行图像预处理(包括大气校正、几何校正)后，计算更具代表性的植被覆盖指数时需要反射值，得到每个多边形区域的ndvi平均值，并用作表示第四个变量(x4)的第二输入。

卫星传感器在获取地物光谱的过程中，来自太阳的辐射能量第一次穿越大气层会被吸收一部分，剩余部分到达地面后，会以透射和反射的方式与地物相互作用，这两种方式的辐射能量之和等于上述剩余部分的辐射能量，在上升通过大气层的过程中，地面物体反射的辐射能量将被大气吸收，最终传感器获得的辐射能量，一部分来自地面物体的反射辐射，另一部分来自大气的反射和散射的辐射能量。遥感影像记录的大气辐射能量，给影像造成一定的失真，为了正确评估遥感影像地物辐射值，就需要去除这些辐射噪声，从遥感图像中去除辐射噪声的过程称为大气校正，目前存在较多的大气校正理论模型，针对不同的遥感影像类型和状况选取合适的模型才可以获得较好的处理结果。

本发明遥感卫星图像在海岸线区域上具有小于10％的云量覆盖。本发明设计具有较好的辐射校正精度的大气校正模型，该模型考虑了地面非朗伯体情况。该模型使用的地面目标反射率与传感器入口处反射率关系式由式(1)所示：

在卫星飞行的过程中，会出现侧滚和距离地面高度的变化，这将会使获取的遥感影像发生几何形变。地面数据接收站会对获取的遥感影像做几何粗校正处理，在实际应用中，一般需要对遥感数据再进行几何精校正，在几何精细校正过程中，需要地面控制点gcps(地面控制点)。具体的操作步骤是地面控制点选择，像素坐标变换和像素亮度值的重新采样。像元坐标变换采用的是基于坐标变换的的方法，变换的公式为

x’＝a0+a1x+a2y+a3x²+a4xy+a5y²(2)

y’＝b0+b1x+b2y+b3x²+b4xy+b5y²(3)

其中，x’，y’--校正后输出的控制点坐标，x，y--地面控制点在原图像中的坐标，rmserror--每个控制点的均方根误差。

采用上述变换公式，将原始影像中的像元值坐标(x’,y’)转换为输出影像中的像元坐标值(x,y)，通过计算每个地面控制点的均方根误差(rmserror)，可以验证校正模型的有效性。

其中，ρ为地表反射率，ρτ(θs,θv,φs,φv)为大气上界反射率，θs是太阳天顶角，θs是太阳方位角，θv是传感器方位角，t(θs)是大气透射率，tg(θs,θv)是太阳-目标大气路径透射率，t(θv)是目标-传感器大气路径透射率，ρr+α是分子散射和气溶胶散射所构成的路径辐射反射率，s是大气半球反射率。

利用大气校正模型能够去除天空光和大气散射效应，实现遥感图像的大气校正。使遥感图像更准确地反映特征的光谱值。

沉降演变所使用的第三输入数据(y)是沉降影响(构建区域)。海岸线附近的基础设施和建筑物可直接影响海岸侵蚀以及洪水。在海岸线附近保持地形学方面的最小距离对于减少海岸带对人类的影响是非常重要的，然而，在本发明中，可以观察到，即随着时间的推移，在海岸线附近的沉降增加。采用数学函数(三角形或l-函数)根据特定范围和变量单位(x1、x2、x3、x4和y)对每个特定语言标记(低、中、高)。

而且，本发明基于对象的算法，通过分析集成邻域信息，将任何图像视为对象，这将增强分析并增加分类图像即lulc的准确性。因此，对于此，使用分割方法才用特征提取工具从每个卫星图像提取lulc类别(建筑物、植物和其它)。在该过程期间，测试缩放和合并级以获得三个类别的最佳结果，包括在所有卫星图像中的构建区域。由于模糊输入需要精确的结果，因此，将分割的光栅转换为矢量数据集，以更精确地描绘三个lulc类别在编辑会话期间gis环境中的每年的多边形区域。

建立针对每个扇区和时间图像分类区域的变量，然后用作海岸线模糊分类模型的第五变量(x5)。

模糊模型设计综合三个输入：海岸线变化、ndvi和沉降影响(建成区)，开发模糊模型设计。这三个输入由五个变量(x1、x2、x3、x4和x5)组成。从卫星图像中提取基线信息，提取所检测到的所有输入变量，根据具体的变量特征，输入不同的范围和单位。在这种情况下，考虑到多边形海岸线变化的总量(％)，x1、x2、x3(海岸线变化)在0到100(％)的范围内。ndvi(变量x4)，范围从-1.0到1.0，x5(建筑物)范围从0到100，然后通过时间变化进行评估。该模糊模型的输出是范围从0到1的数字，对海岸带人为影响进行分类分级。当输出数接近1时，其表明高的人为影响分类，而接近0是指低的人为影响分类，在该范围之间，模糊逻辑可以根据模型设计分类为低、低/中、中/高或高。所提出的模糊模型中使用的推理方法是基于模糊规则概念和表示输出的模型种类，通过由每个推理规则的聚合产生模糊集，在模糊模型中，第一输入(海岸线变化)被分成三个变量(基于海岸线的状态)、侵蚀(x1)、堆积(x2)和稳定(x3)，考虑到比较连续几年检测到的变化。为该变量考虑的语言标签被命名为低、中和高。根据表中给出的参数，所选择的隶属函数的类型是三角形或l函数。

最后，使用三个输入(海岸线变化、ndvi和沉降演变影响)、(建立区域)与五个变量(x1、x2、x3、x4和x5)构成模糊规则。规则由五个变量(x1，x2，x3，x4和x5)和它们各自的语言学标记(a1，a2，a3)，(b1，b2，b3)，(c1，c2，c3)，(d1，d2，d3)，(e1，e2，e3)组成。通过对五个变量进行积分来定义最终模糊规则输出y(f1，f2，f3)，其语言学标签使用"如果-则"规则格式。

例如，ifx1εa1andx2εb1andx3εc3andx4εd3andx5εe1thenyεf1，使用例如语言变量表达相同规则的方式：“如果侵蚀低且增长低且稳定高且ndvi高且累积低，则输出海岸带人类影响分类低”。所有这些变量和规则的集合需要被验证，否则其可能被分类为任意估计。在这种情况下，进行一些初步测试以便以交互形式微调这些功能的规则和参数，直到满足验证标准。验证步骤用于确定所实现的最终输出(模糊海岸带人为影响分类)的准确性和质量。

传统的数据验证问题中，数据集往往是静态且数量固定的。但在海岸线变化问题中，数据会具有两个重要特性，即演化特性和数量特性。演化特性是指数据的统计性质(例如均值、方差等指标)随着时间动态变化，具有潜在的内部变化规律以及时变性。这种内部变化规律既可能是由于该数据所描述的海岸线系统自身的自然且合理的长期变化，也有可能是由于数据噪声干扰导致的短期变动。数量特性则是指数据的规模庞大。由于在每一个固定的时间单位都会产生一批新的数据，系统的数据常常处于不断更新中，从而造成了总数目不断增加。因此，具有演化特性和数量特性的数据被称为演化数据。该类数据给数据验证带来了较大的难度。数据的演化特性要求算法能够处理数据随时间变化的趋势，合理的分析相邻时间段之间数据所发生的转变。数据的数量特性则要求算法能够在有限次遍历内完成分析，而非将整体数据全部存储于内存再运行算法聚类。

为了使聚类结果能够体现数据在时间域上的关联性并实时获得聚类的结果，本发明在对当前时间段的数据进行聚类时会考虑上一时间段的聚类结果。具体来说，每个新的时间段都有一批新的数据等待着被聚类算法归入某些聚类簇中。如果他们没有背离历史期望，而是与历史聚类结果相似，那么整个数据的聚类结果应该是与历史聚类结果相近的，聚类算法生成的整体聚类模型应该也是类似的。但如果他们的内部结构变化剧烈，那么算法应该修改聚类结果以反映出这种变化。演化聚类算法一方面要求在分析当前时间段的数据的同时能够结合之前时间段的模型，得出多个时间段之间数据统计信息的联系，另一方面要求聚类算法能够应对大规模的数据进行实时聚类得到结果。因此，为了刻画这两方面因素对聚类结果的影响程度，定义如下两个度量标准函数：快照质量函数(snapshotquality)和历史代价函数(historycost)。快照质量函数衡量的是当前时间段变量x1、x2、x3、x4和y的聚类划分质量，表示聚类结果和当前时间段的数据的匹配程度。历史代价函数衡量的是相邻时间段变量x1、x2、x3、x4和y的聚类划分是否具有较好的时序平滑特性，表示的是历史模型与当前时间段的聚类结果的差异。

所以，全局目标函数包含两个部分。一部分是对当前样本数据聚类质量的度量，通过快照质量函数sq(cτ,mτ)衡量。cτ表示τ时间段数据的聚类结果，mτ是时间段数据的相似度矩阵。该函数值越小，表示聚类后的模型和该时间段数据的潜在真实模型越符合。另一部分则是历史数据聚类结果与当前时间段的聚类结果之间的差异度，通过历史代价函数hc(cτ,cτ-1)衡量。该函数值反映了当前时间段的数据模型和上一时间段的数据模型之间的差异度。该值越小，表示相邻时间段的聚类结果越相似。故在对τ时间段的数据进行聚类时，为了同时考虑上述两个度量函数，将目标函数定义为：

对于包含t个时间序列数据s＝{x1,x2,x3,x4,x5}，且每个时间点的数据集x包含n个样本，xτ＝{x1,x2,x3,x4,x5}。定义如下的快照质量函数sq(cτ,mτ)以及历史代价函数hc(cτ,cτ-1)：

将上述函数(3)、(4)带入式(2)，得到全局目标函数：

式中，c、n分别是聚类数和样本数，表示τ时间段的数据与τ时间段第i个分段(第i个聚类)的距离，表示τ时间段数据属于第i个聚类的隶属度，且满足每个样本属于所有类的隶属程度之和等于1，cτ指τ时间段的隶属矩阵，cτ∈{0,1}^c*k，而ci,k＝1表示τ时间段数据属于第i类，各数据在各时间段内只能属于某个确定的类，即∑τci,k＝1。

最小化全局目标函数以求得最优分割参数，通过拉格朗日定理求解，定义函数(6)为拉格朗日目标函数：

通过求解对x，μ，λ的偏导数，并令这些偏导数分别为零。此时便可获得最优聚类中心及各个数据隶属于各类的隶属度。

在时序分割问题中，每个聚类的数据点都必须符合自身固有的时间顺序，因此将数据的时间坐标作为一个额外变量并使其参与聚类过程，所得到的聚类结果便是动态聚类结果的计算方法。迭代计算方式如下：

计算τ时间段c个聚类中心vi：

计算聚类中心的模糊协方差矩阵fi：

计算距离函数值d²(xk,vi)：

利用距离函数值d²(xk,vi)更新模糊划分矩阵u：

当满足条件||u^(l+1)-u^(l)<ε||时，中止聚类算法，否则增加迭代次数，使l＝l+1，转到步骤(2)，重复上述步骤，直至满足该条件为止。

通过本发明的算法，可以有效地处理分析具有非线性及强耦合特性的海岸带数据，使得改进算法能够实时分析线性时间序列数据，加入了快照质量函数和历史代价函数，增强了聚类结果和当前时间段数据的匹配程度，且使得相邻时间段的聚类划分具有良好的时序平滑特性。

以上所述仅是对本发明的较佳实施例而已，并非对本发明作任何形式上的限制，凡是依据本发明的技术实质对以上实施例所做的任何简单修改，等同变化与修饰，均属于本发明技术方案的范围内。