一种水体总磷浓度的反演方法

技术领域

本发明属于水质监测领域技术领域，尤其涉及一种水体总磷浓度的反演方法。

背景技术

水是生命之源，水域环境的健康与人类的命运息息相关。近些年来，随着经济的快速发展和人类活动的不断增多，我国内陆水体水质进一步恶化，水体富营养化问题日益严重，一些水体出现水华暴发和水体缺氧等现象，水质已沦为Ⅴ类或劣Ⅴ类，总磷(TP)是主要污染指标之一。每年大量的城市生活污水和工农业废水使得大量的氮、磷的污染元素直接或间接的进入内陆水体，以致水体中营养盐(氮(N)、磷(P))浓度普遍偏高。水域污染事件的频频发生，已严重威胁到人民的日常活动和生命健康，水体修复与水环境保护迫在眉睫。总磷(TP)浓度是水体污染程度的主要评价指标，总磷浓度的准确获取，能够为磷类污染物精准防治和水体修复提供重要的技术支撑。

目前，水体总磷浓度监测有人工移动监测、固定站点监测和遥感监测等方式。人工移动监测，依赖野外采样和实验室分析，工作量大、效率低，只能对特定时间段某一有限区域进行监测，难以反映全区域水质情况；固定站点监测，通过特定波长处吸光度的换算，可以实现定点实时在线监测，无法满足大范围同步动态监测的需求；遥感监测，包括可见光、红外光谱等监测手段，具有获取数据快、覆盖范围广和同步观测等优势。光谱遥感能够获取丰富精细的光谱信息，对水体总磷、总氮等元素感知度更灵敏、估算精度更高，已成为当前最有效的遥感监测手段。卫星光谱遥感成本较高、分辨率较低、很难推广应用。地面光谱遥感具有成本低、操作简单、受天气影响较小、无需大气校正等特点，是水体参数定量化研究的有效获取手段。

研究表明风速的变化会打破底泥与水之间的平衡关系，引起底泥再悬浮，进而导致底泥中的磷释放进入水中。不同类型的风浪会导致不同程度的沉积物再悬浮和磷释放，强风和中风条件下，悬浮固体和总磷均显著增加，且从沉积物中释放的磷主要以颗粒形式存在。短期强风条件导致更多的总磷释放，但风停止后总磷随悬浮固体迅速沉降到底部，并且没有促进藻类的生长。长期中风条件下，悬浮颗粒磷被磷酸酶水解为可溶性反应性磷，并促进藻类生长。藻类的生长反过来会分泌更多的磷酸酶，并促进颗粒物磷的再生。有研究认为气象因素可以解释太湖总磷浓度的58％年际变化，其中平均风速为主导气象因素，较低的风速和较长的低风持续时间导致湖底附近的溶解氧浓度较低，通过反硝化作用可以增加沉积物中磷的释放，增加氮损失，导致表层总磷浓度变高。

温度的变化会响水体中藻类、水生植物及微生物的生长，此类物质的生长消亡会影响总磷浓度的变化，进而对总磷的季节变化也有一定的影响。水体中各形态磷浓度在水平空间上沿水流方向递减，在垂直方向上各形态磷浓度向下递增，并受水温分层的影响。当湖水温度升高时，藻类植物繁殖加快，减少了上覆水体中磷的浓度，增加了对磷的需求，从而使平衡向着有利于磷释放的一方移动，促进沉积物中磷的释放，同时水生生物活动加剧，也促进了沉积物磷的释放。研究表明黑臭河沉积物中可溶性磷的浓度随着温度的升高而显着增加，与15℃时的浓度相比，在25℃和35℃时分别增加28％和87％，羟基氧化铁的还原和溶解性有机物竞争性吸收矿物质是沉积物中磷释放原因，而温度与这两者呈现正相关。

但是，当前总磷浓度遥感反演使用的数据多局限于水体反射光谱，结合环境气象要素的研究甚少，而且水体总磷浓度受到环境气象要素的影响，具有明显的区域性差异。因此，迫切需要一种考虑环境气象要素影响的水体总磷浓度遥感监测方法。

发明内容

基于上述存在的问题，本发明的目的是提出一种水体总磷浓度的反演方法，可作为传统水质监测手段的辅助方法，节省人力、财力，提高水体总磷浓度监测的广度、质量和效率。

本发明提供的一种水体总磷浓度的反演方法，包括：

S1、采集测试点的水体反射光谱、气象要素数据；

S2、对所述水体反射光谱和气象要素数据进行预处理；

S3、将所述预处理后的数据输入改进的随机森林，获得测试点的水体总磷浓度。

进一步地，所述气象数据为温度和风速。

进一步地，所述步骤S2中，所述预处理的方法包括：

S21、对同一地点采集的水体反射光谱构成的曲线进行平滑降噪，并获取A个光谱波段特征；

S22、对所述水体反射光谱进行归一化，公式如下：

式中，X

进一步地，所述步骤S2中，对同一地点采集多次数据，使用平均值作为最终数据。

进一步地，通过Savitzky-Golay平滑法进行所述平滑降噪。

进一步地，所述步骤S3中，改进的随机森林，利用相关系数对决策树结果进行加权以获取最终总磷预测结果并评价。

进一步地，所述步骤S3中，所述改进的随机森林的训练方法包括：

S31、训练样本数据集X，样本总数为p，样本包含光谱波段特征和气象要素特征，样本对应总磷浓度Y，A为光谱波段特征的数量，B为气象要素特征的数量，计算每个特征的皮尔逊相关系数，得到皮尔逊相关系数集合r＝{r

S32、使用训练样本数据集X与总磷浓度Y作为输入，采用bootstrap抽样技术从原始训练集随机产生t个训练子集，利用每个训练子集，生成对应的决策树{h

S33、从光谱波段特征中随机选择q(q≤A)个特征与B个气象要素特征共同作为候选特征，在决策树的每个节点按照预定的规则选择最优属性进行分裂,每个决策树都得到最大限度的生长，过程中完全分裂不剪枝；

S34、计算每棵决策树h

其中，q+B为构建该决策树使用到的特征数量，R

然后对p

S35、对于特征向量X，总磷浓度预测的结果表示为：

其中，

进一步地，计算每个特征的皮尔逊相关系数的方法为：

训练样本数据集X，样本总数为p，特征数量为m，样本对应总磷浓度Y，根据公式

进一步地，对反演模型性能进行评价的指标为决定系数R

式中：p代表样本数；y

进一步地，所述训练样本数据集通过步骤S1、S2构建。

本发明所述的一种水体总磷浓度的反演方法，具有如下有益效果：

(1)本发明对光谱进行了平滑、归一化，降低了背景噪声对总磷反演的影响。

(2)本发明针对当前水质监测局限于水体反射光谱，而忽略环境要素对反演精度的影响的问题，在建模过程中考虑环境要素，解决了环境要素对总磷浓度反演效果的影响，所构建模型具有预测精度高，对环境差异具有较好适应能力等特点。

(3)本发明针对传统随机森林中不同泛化能力决策树具有相同权重的问题，将提出的改进随机森林方法应用于光谱数据处理领域，实现了总磷浓度的预测，提高了总磷浓度反演精度和效率。

(4)本发明的方法适用性较高，具有高精度、高灵活性和强稳健性，为水体总磷浓度监测提供了新思路，提高总磷监测的广度、质量和效率。

附图说明

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

图1为本发明一个实施例的反演方法的流程示意图图；

图2为本发明一个实施例的改进的随机森林的训练流程示意图图。

具体实施方式

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

本发明的思路是：首先采集样本点的水体反射光谱、环境气象要素和水样。其中水样需要及时送至实验室进行化学分析，获取水体总磷浓度。然后对样本的反射光谱和环境气象要素进行预处理，使用预处理结果和总磷浓度作为输入，基于改进的随机森林算法建立总磷反演模型。模型建立完成后，按照相同流程对未知水体反射光谱和环境气象要素进行预处理，输入改进随机森林反演模型即可获取水体总磷浓度。

本发明提出的水体总磷浓度的反演方法，如图1所示，包括：

S1、采集测试点的水体反射光谱、气象要素数据；

S2、将所述水体反射光谱和气象数据进行预处理；

S3、将所述预处理后的数据输入改进的随机森林，获得测试点的水体总磷浓度。

在步骤S1中，光谱中包括多个波段，如270个波段，气象数据包括风速、温度等，一般来说波段数量远远大于气象因素，因为对总磷浓度有影响的气象要素是有限的。

训练改进的随机森林所需的数据集可以在水体处进行采集。比如，采集多个样点的水体反射光谱、环境气象要素和水样，其中同一样本采样点获取多条光谱曲线，水样当天送至实验室内进行分析获取总磷浓度。

在步骤S2中，要在对获取的数据(包括测试数据和待分析数据)进行预处理和归一化。首先，将同一样本点获取的多条光谱曲线，使用Savitzky-Golay平滑方法(当然，还有很多方法可以平滑去噪)对每条光谱曲线进行平滑，取这些曲线的平均值作为该样本的最终水体反射光谱曲线。

归一化处理使用线性函数归一化方法，该方法实现对原始数据的等比例缩放,线性函数将原始数据线性化的方法转换到[0,1]的范围，归一化公式如下：

式中X

在步骤S3中，将水体反射光谱P、风速U、温度T的归一化特征P

传统随机森林由多个决策树组成，对于随机森林{h

S31、计算得到每个特征变量与总磷浓度之间的相关关系，作为特征权重的分配依据，因为相关系数取值在-1到1之间，当相关系数小于0时，说明两者呈现出负相关性，两者之间的影响仍然存在，因此在利用特征的相关系数计算特征权重的时候应该对相关系数求绝对值，利用相关系数的绝对值进行对应的特征权重的计算。具体来说，训练样本数据集X，样本总数为p，样本包含水体反射光谱波段特征和气象要素特征，样本对应总磷浓度Y，A为光谱波段特征的数量，B为气象要素特征的数量，计算每个特征的皮尔逊相关系数，得到皮尔逊相关系数集合r＝{r

S32、使用训练样本数据集X与总磷浓度Y作为输入，采用bootstrap抽样技术从原始训练集随机产生t个训练子集，利用每个训练子集，生成对应的决策树{h

S33、从光谱波段特征中随机选择q(q≤A)个特征与B个气象要素特征共同作为候选特征，在决策树的每个节点按照预定的规则选择最优属性进行分裂,每个决策树都得到最大限度的生长，过程中完全分裂不剪枝；

步骤S32、S33与传统的随机森林的方法一样，此处不再赘述。

S34、对于每棵决策树h

上式中p

S35、对于特征向量X，总磷浓度预测的结果表示为：

上式中，各决策树的预测结果

在训练时，对随机森林算法中较为重要的两个超参数：决策树棵数t和候选特征数q进行调节优化。在所有候选的参数选择中，通过循环遍历，尝试每一种可能性，将输出的总磷浓度预测值与真实值进行比较，选择出表现最好的参数作为最终的结果(调节优化属于现有技术)。

对模型性能进行评价，评价指标为决定系数R

/>

式中：p代表样本数；y

模型训练建立后，在应用时，采集目标水体数据，对目标水体反射光谱以及环境气象要素作为输入，即可进行预处理和归一化，经改进的随机森林预测后，即可得到该水体的总磷浓度。

实施例1

依据本发明对某一未知总磷浓度的水体进行监测。光谱数据获取时间为北京时间上午9点-11点或下午14点-16点，波长范围350～1000nm，通道数为270，光谱分辨率优于1nm。

(1)数据采集。采集水体反射光谱曲线、环境气象要素(风速、温度)以及水样。

步骤1：采集水体反射光谱。采集水体反射光谱时需要选择晴朗无云的天气，采集时间在上午9点-11点或下午14点-16点，采集时要确保水体、太阳和设备探头之间无遮挡，选择合适的几何观测角度，每个样点重复采集10条光谱。

步骤2：采集环境气象要素。采样点处风速和温度通过手持风速仪获取，采集时需确保没有其他物体挡住风口。

步骤3：总磷浓度获取。通过塑料水瓶获取采样点处水样，水样需要及时送往实验室进行化学分析，继而获取水体的总磷浓度。

(2)数据预处理。

步骤1：对每个样点重复采集的10条光谱使用Savitzky-Golay平滑法进行平滑处理，降低水体反射光谱中的噪声，使光谱更加平滑，水体反射光谱记为P，每条水体反射光谱包含270个光谱特征。

步骤2：取步骤1获取的10条平滑后的平均曲线作为该样点的反射光谱曲线，以获取可靠的水体反射光谱。

步骤3：重复步骤2，对所有样点的最终光谱特征进行归一化处理，水体反射光谱数据集为：

归一化处理使用线性函数归一化方法，该方法实现对原始数据的等比例缩放,线性函数将原始数据线性化的方法转换到[0,1]的范围，归一化公式如下：

式中，X

将p个样点获取的所有归一化特征P

步骤4：将步骤3获取的所有样点的归一化结果X，与(1)中获取的总磷浓度(Y＝{y

(3)将(2)中的归一化数据X和总磷浓度Y输入到改进随机森林中进行训练，利用相关系数对决策树结果进行加权以获取最终总磷预测结果并评价。

改进的随机森林的训练方法包括：

步骤1：训练样本数据集X，样本总数为p，特征数量为m，样本对应总磷浓度Y，根据公式

步骤2：使用训练样本集X与总磷浓度Y作为输入，采用bootstrap抽样技术从原始训练集随机产生t个训练子集，利用每个训练子集，生成对应的决策树{h

步骤3：从光谱特征中随机选择q(q≤n)个光谱特征，与两个环境特征一起作为候选特征，在决策树的每个节点按照预设的规则(基尼指数、信息增益率等)选择最优属性进行分裂,每个决策树都得到最大限度的生长，过程中完全分裂不剪枝。

步骤4：传统随机森林回归，一般通过计算所有决策树给出的预测结果的平均数作为最终的预测结果。改进随机森林考虑到不同决策树的不同泛化能力，结合特征的相关系数赋予每个决策树不同的权重，将所有决策树的加权和作为最终的预测结果。

对于生成的随机森林模型为{h

上式中p

随机森林模型{h

上式中，各决策树的预测结果

对模型性能进行评价，评价指标为决定系数R

式中：p代表样本数；y

(4)目标水体总磷浓度反演。按照(1)中的标准采集未知水体的相关数据，采用(2)中相同的流程对数据进行预处理，最终得到该水体的归一化数据，将归一化结果输入到(3)构建的改进随机森林模型中得到该水体的总磷浓度。

将本发明结果与不考虑环境特征的传统随机森林总磷反演模型以及考虑环境特征的传统随机森林总磷反演模型进行对比，模型精度如表1所示：

表1 模型精度对比

本领域普通技术人员可以理解：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明权利要求所限定的范围。