CN113255961A - 基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法 - Google Patents

Abstract

本发明提供一种基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，包括步骤一、使用长时期的中+高分辨率多光谱遥感观测数据和地面观测数据，基于机器学习理论建立河湖水环境参数反演模型；步骤二、基于所述河湖水环境参数反演模型进行水环境参数时空统计分析，得到水环境参数时空动态演变规律，并基于所述水环境参数时空动态演变规律设计地面水环境监测站网。本发明通过湖泊水环境历史变化规律分析来布设点位，以实现用尽可能少的地面监测点结合多光谱遥感卫星能获取满足尽可能高精度的湖泊水环境参数浓度分布。

Description

基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局 方法

技术领域

本发明涉及水环境监测技术领域，具体是一种基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法。

背景技术

当前，水环境监测方式主要采用布设若干个监测站网，定期采集水环境数据。监测站点的布设通常考虑周围的地理环境，在易于施工的区域布设监测站点。而定期的水环境采样则是基于等间隔的水平和垂直相交的格网，在所监测的湖泊区域中，以特定的空间距离间隔均匀获取采样点水样。采样点的密度和位置的确定以人的经验知识为主。然而实际情况中水环境参数浓度，如叶绿素、悬浮颗粒物等的空间变异并非完全独立，相邻的点之间具有一定的空间相关性，可以认为其浓度的变化是连续的。稀疏的采样点容易造成水环境分布模拟的精度达不到要求，而为了获得高精度的水环境分布数据，往往需要在一个湖泊区域内采集较多的水环境点，这不仅花费了大量的人力物力，同时可能造成采样区域的样点冗余。若在监测站点和定期采样布局设计时，能够充分考虑水体的空间相关性，了解该区域的水环境演变规律，依据水环境的空间变化特性来布设地面监测点，则能利用尽可能少的点来精确描述水体整体和局部的变化情况。

目前，通过公开文献可以查阅到关于湖泊地面水环境监测站点优化布局的专利申请及学术期刊论文，包括授权公告号为CN107194040B的专利“一种基于双向水流的水环境监测网络多目标优化部署方法”，该方法考虑水流方向，模拟计算水污染扩散模型，没有考虑到不同湖泊的水环境特性，未考虑历史水环境参数统计分布规律；

授权公告号为CN104270773B的专利“一种基于遗传算法多目标优化的流域传感器覆盖网优化方法”，只考虑使用遗传算法多目标优化将河流流域传感器覆盖网优化问题转化为0/1多目标规划问题，网点布设的经济效益，没有估计流域内部的真实水环境分布，布设的点位是否能够真实监测整个湖泊内的水环境变化具有不确定性。

论文“基于时序遥感数据的水环境时空采样布局优化方法研究”，提出了时空维度变化兼顾的湖泊水环境变异特征先验规律分析方法，引入模拟退火优化算法，实现了鄱阳湖的采样布局优化，然而该方法对仅考虑叶绿素a的反演结果，同时采用k均值聚类的方式无法自动实现分区的分类与合并，另外其使用的Sentinel-3 号卫星空间分辨率为300米，中分辨率卫星观测数据无法反映湖泊水环境空间变异的细节。

本发明以多光谱遥感卫星的长时期监测数据为依据，根据时序遥感影像反演水色参数，获得湖泊水色参数时序变化，建立湖泊水环境季节性和空间性变化规律，并以此为基础辅助水环境时空采样点优化。

发明内容

本发明的目的在于提供一种基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，在水环境监测点位设计过程中，顾及到特定湖泊的水环境的时空特异性，通过湖泊水环境历史变化规律分析来布设点位，以实现用尽可能少的地面监测点结合多光谱遥感卫星能获取满足尽可能高精度的湖泊水环境参数浓度分布。

本发明采用如下技术方案实现：

一种基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，包括如下步骤：

步骤一、使用长时期的中+高分辨率多光谱遥感观测数据和地面观测数据，基于机器学习理论建立河湖水环境参数反演模型；

步骤二、基于所述河湖水环境参数反演模型进行水环境参数时空统计分析，得到水环境参数时空动态演变规律，并基于所述水环境参数时空动态演变规律设计地面水环境监测站网。

进一步的，所述步骤一具体包括：基于中分辨率多光谱遥感卫星观测数据的时序序列，采用经验模型反演时序结果得到叶绿素a(Chl-a)与总悬浮颗粒物(TSS) 浓度分布，采用非监督自聚类对湖泊进行预分区，获得湖泊浓度分区结果；再根据分区结果在分区内部使用高分辨率多光谱遥感卫星和地面监测数据，基于机器学习理论中的回归方法分别建立分区水环境参数(Chl-a与TSS浓度)反演模型，从而得到河湖水环境参数反演模型。

进一步的，根据分区结果在分区内部使用高分辨率多光谱遥感卫星和地面监测数据，基于机器学习理论中的回归方法分别建立分区水环境参数(Chl-a与TSS 浓度)反演模型具体为：地面采样点与高分辨率多光谱遥感卫星进行同步观测，采集获取地面水环境参数Chl-a和TSS浓度作为真值，多光谱遥感卫星的波段组合作为特征向量，根据预分区结果，采用机器学习方法分区建立卫星遥感的水环境参数反演模型。

进一步的，所述步骤二具体包括：

(1)水环境参数时空统计分析：根据水环境参数反演模型训练结果与历史多光谱遥感卫星观测数据获得湖泊的历史月度Chl-a/TSS浓度分布序列，进而得到月度直方图；

(2)水环境参数时空动态演变规律：获取Chl-a/TSS年度直方图，计算最大值、最小值、均值点所在区域位置，得到位置时间序列；

(3)对位置时间序列进行统计直方图分析，统计年度均值Avg和方差SD，根据年度均值Avg以及标准差-3SD，-2SD，-SD，+SD，+2SD，+3SD进行湖泊分区，每个分区的中心点位置作为最终的监测点布设位。

进一步的，中分辨率多光谱遥感卫星为Sentinel-3号卫星，高分辨率多光谱遥感卫星为Sentinel-2号卫星。

与现有技术相比，通过本发明能够获得如下优异效果：

1)通过本发明，首先依据中分辨率遥感数据建立先验分区，再分区建立高分辨率遥感数据的水环境参数反演模型，在布设水环境监测站网时，在水污染风险高、中、低的区域分别布设站点，保证监测点位数据能够充分反映整个湖泊的水环境真实状况，分区反演模型精度相较于全局反演，其模型精度可提高20％；

2)通过本发明能够通过少数若干个地面监测站点辅助多源光谱遥感数据建立高精度的湖泊水环境参数反演模型，仅由少量的监测站点数据能准确获取水环境要素时空变化规律先验信息，实验结果表明，相对于等距间隔采样，通过本方法模型优化后，地面采样点位数可减少50％以上。

附图说明

图1是本发明基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法其中一个实施例的流程图；

图2是本发明基于水环境参数时空动态演变规律设计地面水环境监测站网的示意图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明提供一种基于长时间时序多源多光谱遥感影像分析的流域地面监测站网优化方法，其使用长时期的中+高分辨率多光谱遥感(如Sentinel-2/Sentinel-3 号卫星)观测数据和地面观测数据，基于机器学习理论建立河湖水环境参数反演模型，获取河湖水环境参数时空动态演变规律，并基于该规律设计地面水环境监测站网。本发明的算法使得在地面水环境监测站网布设时充分考虑的监测湖泊的水环境参数时空变化特性，通过自聚类算法对湖泊建立预分区，再根据分区分别建立水环境参数反演模型，获得多组水环境反演水环境参数反演模型，根据湖泊的水环境参数浓度时空变化自动提取监测站点位置，大大提高了监测站点布设的准确性和监测效率。

请参阅图1，本发明包括如下步骤：

步骤一、使用长时期的中+高分辨率多光谱遥感(如Sentinel-2/Sentinel-3号卫星)观测数据和地面观测数据，基于机器学习理论(如SVR支持向量回归模型理论)建立河湖水环境参数反演模型；

所述步骤一具体包括：

(1)用中+高分辨率遥感卫星观测数据(如实例Sentinel-3号数据与 Sentinel-2号数据)进行预处理：包括对Sentinel-2/Sentinel-3号多光谱卫星影像进行几何校正、裁剪配准、瑞利校正；

(2)湖泊预分区：基于中分辨率多光谱遥感卫星遥感(如Sentinel-3号卫星观测数据)时序序列，采用经验模型(如实例中的C2RCC方法)模型反演时序，获得水体中叶绿素a(Chl-a)和总悬浮物(TSS)浓度估算值，获得湖泊水体中叶绿素a和总悬浮物浓度时序向量SChl-a和TSS；基于SChl-a和STSS，采用非监督聚类方法(如实例中的ISODATA算法)进行自聚类，初始值按照浓度值高、中、低分为3类，获得湖泊浓度预分区结果。

(3)水环境参数反演模型构建：获得湖泊浓度分区结果后，再根据分区结果在分区内部使用高分辨率多光谱遥感卫星和地面监测数据分别建立水环境参数反演模型。地面采样点与高分辨率多光谱遥感卫星(如Sentinel-2号卫星观测数据)进行同步观测，采集获取地面水环境参数Chl-a和TSS浓度作为真值，多光谱遥感卫星的波段组合作为特征向量，根据步骤(2)中的预分区，采用机器学习方法(如实例中的SVR支持向量回归模型)分别建立分区水环境参数(Chl-a 与TSS浓度)反演模型，从而得到河湖水环境参数反演模型。

步骤二、基于所述河湖水环境参数反演模型获取河湖水环境参数时空动态演变规律，基于所述水环境参数时空动态演变规律设计地面水环境监测站网，如图 2所示；

步骤二具体包括：

(1)水环境参数时空统计分析：根据水环境参数反演模型训练结果与历史多光谱遥感卫星观测数据(如图2所示遥感影像1，遥感影像2，……数据N) 获得湖泊的历史月度Chl-a/TSS浓度分布序列(如图2所示Chl-a/TSS分布图1， Chl-a/TSS分布图2……Chl-a/TSS分布图N)，进而得到月度直方图。

(2)水环境参数时空动态演变规律：根据(1)得到的历史月度浓度分布序列获取，通过直方图统计分析获得Chl-a年度直方图及TSS年度直方图(年度直方图为全年月度直方图均值)，获得最大值、最小值、均值点所在区域位置，得到位置时间序列。

(3)统计中心位提取：对位置时间序列进行统计直方图分析，统计年度均值Avg和方差SD，根据年度均值及标准差(Avg、-3SD，-2SD，-SD，+SD，+2SD， +3SD)进行湖泊水环境风险等级分区，每个分区的中心点位置作为最终的监测点布设位。

本发明主要创新点为：

(1)采用“中分辨率数据预分区——高分辨率数据分区建模”双层水环境参数反演模型。实现首先用中分辨率多光谱遥感卫星观测数据经验模型算法得到Chl-a和TSS浓度分布，基于此结果采用非监督聚类对湖泊进行预分区，再根据分区结果在分区内部分别建立水环境参数反演模型，以提高水环境参数反演模型反演精度。

(2)基于时序统计分析候选监测站点的自优化方法。基于时序浓度分布的统计直方图分析，统计年度均值Avg和方差SD，根据度均值及标准差(Avg、-3SD， -2SD，-SD，+SD，+2SD，+3SD)进行湖泊风险等级分区进行湖泊分区，每个分区的中心点位置作为最终的监测点布设位。

下面以Sentinel-2/Sentinel-3号卫星遥感数据处理为例对步骤一进行具体说明：

1、多源多光谱卫星影像预处理：对多源遥感卫星观测数据进行常规几何校正、裁剪配准、瑞利大气校正，统一其空间坐标系统，去除大气影像。

2、湖泊预分区：

首先采用经验模型C2RCC方法对Sentinel-3号1级影像进行水色参数反演。以1c级数据为输入，选择适合研究区域及季节的大气校正参数，在算法中通过各有关的神经网络完成大气校正、水体光学量反演，得到各个组分吸收/散射系数等固有光学量参数，同时也计算其各个结果的数值范围不确定度等辅助信息。基于C2RCC的训练数据集中水体固有光学量与同步实测水环境参数的经验统计关系，可将富有植物色数吸收特性IOP_apig和一般颗粒物后向散射IOP_bp及钙质白色粒子后向散射IOP_bw转换为叶绿素chl-a浓度以及总悬浮物TSS干重。默认缺省转换公式为：

TSS＝1.72×IOP_bp+3.1×IOP_bw(g·m^-3)

基于Chl-a和TSS浓度分布图，采用非监督分类ISODATA算法进行自聚类，初始值为3类(浓度值为高、中、低)，获得湖泊浓度预分区结果。

以时序的Sentinel-3号数据叶绿素反演结果为依据，对于特定像元，在时间纬度组成特征向量，采用非监督聚类算法实现其时空特征聚类分区研究。本研究采用ISODATA算法，其基本工作原理是：

(1)输入N个模式样本{x_i,i＝1,2,…,N}，确定N_c个聚类中心

(2)将N个模式样本分给最近的聚类S_j；

(3)若S_j中的样本数目S_j<θ_N，则取消该样本子集，N_c-1.

(4)修正各聚类中心：

(5)计算各聚类域S_j中模式样本与各聚类中心间的平均距离

(6)计算全部模式样本和其对应聚类中心的总平均距离

(7)判断分裂、合并及迭代运算

若迭代次数达到阈值，则转至第(11)步。

若

即聚类中心的数目小于或等于规定值的一半，则转至第(8)步，对已有聚类进行分裂处理。

若迭代运算的次数是偶数次，或N_c≥2K，不进行分裂处理，转至第(11) 步，否则，转至第(8)步，进行分裂处理。

(8)计算每个聚类中样本距离的标准差向量

σ_j＝(σ_1j,σ_2j,…,σ_nj)^T

其中向量的各个分量为：

式中，i＝1,2,…,n为样本特征向量的维数，j＝1,2,…Nc为聚类数，N_j为S_j中的样本个数。

(9)求每一标准差向量{σ_j,j＝1,2,…,N_c}中的最大分量，记为{σ_jmax,j＝ 1,2,…,N_c}代表。

(10)在任一最大分量集{σ_jmax,j＝1,2,…,N_c}中，若σ_jmax>θ_s，同时又满足如下两个条件之一：

1)S_j中样本总数超过规定值一倍以上；

2)

则将z_j分裂为两个新的聚类中心，N_c+1。更新σ_jmax分量。

转至(2)。否则继续。

(11)计算全部聚类中心距离

D_ij＝‖z_i-z_j‖,i＝1,2,…N_c-1,j＝i+1,…,N_c

(12)比较D_ij与θ_c的值，按最小距离次序递增排序，即

(13)将聚类为

的两个聚类中心

和

合并，得到新中心：

式中，被合并的两个聚类中心向量分别以其聚类域内的样本数加权，使

为真正的平均向量。

(14)如果达到迭代次数阈值，则算法结束；否则，转至(2)，迭代次数+1。

ISODATA聚类的类间具有较大的差异性，类内具有良好的聚拢性。基于聚类结果，将整个湖泊区域分为M个聚类分区，记为(S₁,S₂…S_M)每个分区内进行采样点的布设优化，以保证每个分区的水环境特性都能较好的顾及。对于每一个分区，统计其叶绿素浓度的最大、最小、均值，并找到其所在位置。

3、地面采样点与高分辨率多光谱遥感卫星(如Sentinel-2号卫星)进行同步观测，采集获取水环境参数Chl-a和TSS浓度，作为观测真值，卫星观测数据的波段组合作为特征向量，根据步骤(2)中的预分区分别建立卫星遥感的水环境参数反演模型。本实例中水环境参数反演模型采用SVR支持向量退化模型。

对于每一个分区S_i,i＝1,2,…M,采用Sentinel-2号卫星校正后反射率及其波段组合作为输入，构建SVR支持向量回归模型。

表3.2 Sentinel-2号MSI载荷光谱参数

构建的特征空间为：

表3.3输入特征向量空间

特征序号	特征值
		X<sub>1</sub>	B<sub>blue</sub>
X<sub>2</sub>	B<sub>green</sub>
		X<sub>3</sub>	B<sub>red</sub>
X<sub>4</sub>	B<sub>NIR</sub>
		X<sub>5</sub>	B<sub>SWIR-1</sub>
X<sub>6</sub>	B<sub>SWIR-2</sub>
		X<sub>7</sub>	B<sub>green</sub>/B<sub>red</sub>
X<sub>8</sub>	B<sub>green</sub>/B<sub>blue</sub>
		X<sub>9</sub>	B<sub>red</sub>/B<sub>NIR</sub>
X<sub>10</sub>	B<sub>red</sub>/B<sub>blue</sub>
		X<sub>11</sub>	B<sub>red</sub>/B<sub>green</sub>
X<sub>12</sub>	(B<sub>green</sub>-B<sub>red</sub>)/(B<sub>green</sub>+B<sub>red</sub>)
		X<sub>13</sub>	(B<sub>NIR</sub>-B<sub>red</sub>)/(B<sub>NIR</sub>+B<sub>red</sub>)
X<sub>14</sub>	B<sub>green</sub>+B<sub>red</sub>/2
		X<sub>15</sub>	B<sub>red</sub>/B<sub>green</sub>+B<sub>NIR</sub>
X<sub>16</sub>	(B<sub>NIR</sub>-B<sub>SWIR-1</sub>)/(B<sub>NIR</sub>+B<sub>SWIR-1</sub>)
		X<sub>17</sub>	(B<sub>green</sub>-B<sub>NIR</sub>)/(B<sub>green</sub>+B<sub>NIR</sub>)
X<sub>18</sub>	(B<sub>SWIR-1</sub>-B<sub>SWIR-2</sub>)/(B<sub>SWIR-1</sub>+B<sub>SWIR-2</sub>)

采用SVR支持向量回归模型构建卫星遥感的水环境参数反演模型。地面采样点的水环境参数作为y_i，表3.2与3.3所描述的特征空间作为x_i：

问题转化为最小化损失函数：

其中，θ为权值向量，ξ_i和

是正松弛变量，C为惩罚因子。引入拉格朗日乘子法，则SVR回归模型可写为：

其中，

和α_i是朗格朗日乘子，K为核函数，

x为输入的特征向量。和函数采用径向基函数：

输入特征向量寻找局部最优模型，C的取值为1-20，参数γ的取值为1-500。

用80％的数据训练SVR模型，20％数据验证其精度。精度评定采用均方根误差RMSE，相关系数R²以衡量反演模型精度。

下面以Sentinel-2/Sentinel-3号卫星遥感数据处理为例对步骤二进行具体说明：

1、水环境参数时空统计分析：月度Chl-a/TSS浓度分布图计算，根据水环境参数反演模型训练结果与历史Sentinel-2号多光谱遥感卫星观测数据(如图2 所示遥感影像1，遥感影像2，……数据N)获得湖泊的历史月度Chl-a/TSS浓度分布序列(如图2所示Chl-a/TSS分布图1，Chl-a/TSS分布图2……Chl-a/TSS 分布图N)。

2、水环境参数时空动态演变规律：年度直方图分析，根据得到的浓度分布图序列获取月度直方图，通过直方图统计分析获得Chl-a年度直方图及TSS年度直方图(参数序列chla_annual，TSM_annual,统计直方图Hist_chla,Hist_TSM。)，获得最大值、最小值、均值点所在区域位置，得到位置时间序列。

3、统计中心位提取：根据得到的时间序列计算得到直方图均值位置 Hist_{chla_avg}，Hist_{TSM_avg}一倍标准差位置Hist_{chla_SD}，Hist_{TSM_SD}两倍标准差位置Hist_{chla_2SD},Hist_{TSM_2SD}，将直方图进行分区，例如，chl-a直方图分区表示为[-Hist_{chla_2SD},-Hist_{chla_SD}]，[-Hist_{chla_SD},Hist_{chla_avg}]， [Hist_{chla_avg},Hist_{chla_SD}]，

映射到地理区域，并求取地理区域中心位，中心位置序列则为水环境监测站点最终点位。

采用步骤二得到优化点位的地面观测数据+Sentinel-2卫星遥感数据，采用步骤一得到的分区水环境参数反演模型反演湖泊水环境参数Chl-a浓度与TSS浓度 (模型计算值)，并于传统地面点规则间隔采样得到的Chl-a浓度与TSS浓度分布(作为真值)相比较，计算相关系数R²>0.8，则证明通过本发明得出的水环境监测站点最终点位满足要求。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何属于本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

Claims (5)

1.一种基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，其特征在于：包括如下步骤：

步骤一、使用长时期的中+高分辨率多光谱遥感观测数据和地面观测数据，基于机器学习理论建立河湖水环境参数反演模型；

步骤二、基于所述河湖水环境参数反演模型进行水环境参数时空统计分析，得到水环境参数时空动态演变规律，并基于所述水环境参数时空动态演变规律设计地面水环境监测站网。

2.如权利要求1所述的基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，其特征在于：所述步骤一具体包括：基于中分辨率多光谱遥感卫星观测数据的时序序列，采用经验模型反演时序结果得到叶绿素a(Chl-a)与总悬浮颗粒物(TSS)浓度分布，采用非监督自聚类对湖泊进行预分区，获得湖泊浓度分区结果；再根据分区结果在分区内部使用高分辨率多光谱遥感卫星和地面监测数据，基于机器学习理论中的回归方法分别建立分区水环境参数(Chl-a与TSS浓度)反演模型，从而得到河湖水环境参数反演模型。

3.如权利要求2所述的基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，其特征在于：根据分区结果在分区内部使用高分辨率多光谱遥感卫星和地面监测数据，基于机器学习理论中的回归方法分别建立分区水环境参数(Chl-a与TSS浓度)反演模型具体为：地面采样点与高分辨率多光谱遥感卫星进行同步观测，采集获取地面水环境参数Chl-a和TSS浓度作为真值，多光谱遥感卫星的波段组合作为特征向量，根据预分区结果，采用机器学习方法分区建立卫星遥感的水环境参数反演模型。

4.如权利要求1所述的基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，其特征在于：所述步骤二具体包括：

(1)水环境参数时空统计分析：根据水环境参数反演模型训练结果与历史多光谱遥感卫星观测数据获得湖泊的历史月度Chl-a/TSS浓度分布序列，进而得到月度直方图；

(2)水环境参数时空动态演变规律：获取Chl-a/TSS年度直方图，计算最大值、最小值、均值点所在区域位置，得到位置时间序列；

(3)对位置时间序列进行统计直方图分析，统计年度均值Avg和方差SD，根据年度均值Avg以及标准差-3SD，-2SD，-SD，+SD，+2SD，+3SD进行湖泊分区，每个分区的中心点位置作为最终的监测点布设位。

5.如权利要求2所述的基于时序多源光谱遥感数据的湖泊水环境监测站点优化布局方法，其特征在于：中分辨率多光谱遥感卫星为Sentinel-3号卫星，高分辨率多光谱遥感卫星为Sentinel-2号卫星。

Images