基于PLSR的珠江口城市河流水质高光谱反演

黄 华，李茂亿，陈吟晖，陈 耿，刘海龙，邢前国，蔡建楠

(1.广东环境保护工程职业学院，广东 佛山 528216；2.中山市环境监测站，广东 中山 528403； 3.广东省环境监测中心，广东 广州 510308；4.中国科学院烟台海岸带研究所，山东 烟台 264003)

珠江口地区河网密布，河流水质状况直接影响城市的整体风貌和居民生活质量[1]。近年来随着国家《水污染防治行动计划》的深入推进，城市河流的综合整治和水质改善已成为珠江口各市的重要工作任务。为实现对河流水质状况监控及河流整治效果的评估，往往需要配套开展密集、长期的水质监测工作，而采用常规的水样现场采集和实验室分析的方法存在耗时长、花费高、采样点位有限等局限，无法满足新形势下水环境管理工作的需要[2]。

水质遥感技术基于水体水质与光谱特征之间的响应关系实现对水质指标的反演[3]，获得的水质信息具有时间和空间上的相对连续性，可以作为传统水质监测方法的有效补充，近年来得到迅速发展。其中，高光谱水质遥感技术由于具有纳米级的光谱分辨率和连续的光谱谱段，极大地增加了水体光谱信息量，有助于精准建立各种水质参数与光谱响应关系以提高反演模型的精度，成为水质遥感研究领域的热点。目前，相关研究主要面向的监测对象为水库和海洋水体，而珠江口河网水体受陆源物质大量输入和海洋潮汐因素的共同作用，有着与内陆二类水体和开阔大洋水体不同的水体光学特征，开展珠江口水体水质光学遥感研究具有重要的意义[4]。目前城市河流水质遥感有关研究中的目标水质参数主要为水体浊度、悬浮物、叶绿素a、有色可溶性有机物(CDOM)等光学活性参数[5-7]，对于河流水质综合指标的遥感反演研究较少。

本文以珠江口城市中山市为研究区，通过开展河流高光谱监测和同步水质分析，研究受测水体高光谱特征与城市河流水质指数(city river water quality index，CWQI)的关系，并基于偏最小二乘回归(partial least squares regression， PLSR)建立高光谱数据与河流水质指数的反演模型，同时探索水体高光谱监测的最佳光谱分辨率，以期为城市河流水质光谱监测及有关高光谱传感器的设计和应用提供参考。

1 研究方法

中山市地处珠江口西岸，北接广州，与深圳隔海相望，紧邻港澳，处在粤港澳大湾区的地理中心位置。中山市境内河流密布，不同规模的内河涌、排洪渠高达1 000余条。改革开放以来，中山市经济迅猛发展，但随着人口和产业聚集程度不断提高，河流水体受到了一定程度的污染。

1.1 河流光谱监测和数据处理

1.1.1地面高光谱监测

于2019年7—9月选择中山市境内具有代表性的河流水体进行水质高光谱监测，共获得151个点位的现场实测数据。高光谱监测使用美国Ocean Optics公司USB4000光谱仪，光谱分辨率约为 0.20 nm，光谱范围为345.34～1 046.12 nm，视场角为15°。每次测量前均使用标准白板进行仪器定标。测量时将光谱仪探头垂直置于离河流水面 0.5 m 处测量河流水体上行辐射亮度Lu、天空辐射亮度Ls和太阳下行辐照度Ed，每个点位重复进行3次测量后取平均值作为测定结果。测点的光谱反射率Rrs计算公式[8]为

(1)

式中γ为天空光在水体表面的菲涅尔反射率，通常取经验值0.028。考虑到400～900 nm波段范围以外的光谱反射率更易受外界条件及仪器设计误差等因素的影响，可能对监测结果带来信号噪声[9]，因此仅选择400～900 nm范围内共2 534个波段的原始光谱反射率数据用于构建反演模型。

1.1.2无人机载高光谱监测

无人机遥感监测技术具有灵活性大、影像分辨率高、时效性强、成本低等优点，是水环境遥感研究的一个重要发展方向[10]。在利用地面实测高光谱构建CQWI反演模型后，将其推广应用于无人机高光谱影像中，以评价模型应用于无人机载高光谱监测的适用性。无人机载高光谱监测采用北京智科远达数据技术有限公司ZK-VNIR-FPG480型高光谱成像扫描仪，其有效谱段范围为400～1 000 nm，光谱通道数大于或等于270，光谱分辨率小于或等于 2.8 nm，空间分辨率为9 cm(100 m飞行高度)，视场角为26°。高光谱成像扫描仪搭载在大疆M600 PRO无人机平台上，在晴朗天气条件下，根据设定的航飞路线采集研究水体的光谱，并同步采集靶标的光谱数据，对水体数据进行定标运算得到水体反射率数据。无人机航拍区域主要为内河与外江交汇处、主要水工设施周边、码头作业区等潜在的水质突变点区域，以尽可能同时获取不同水质水体的光谱信息。

1.2 CWQI构建和水质监测

目前，国内外对于地表水环境质量评价方法主要有单因子水质类别评价法、水质指数法和不确定性方法三大类[11]，其中水质指数法可综合多项水质指标的监测结果，计算过程简单，评价结果便于比较。本文采用原环境保护部《城市地表水环境质量排名技术规定(试行)》中规定的河流水质指数作为反映河流水质状况的综合指标，单项指标(除溶解氧外)和溶解氧的CWQI计算公式分别为

(2)

(3)

式中：Ii、IDO分别为水质指标i和溶解氧的CWQI；ρi、ρ(DO)分别为水质指标i和溶解氧的质量浓度测定值；ρsi、ρs(DO)分别为水质指标i和溶解氧Ⅲ类地表水对应的质量浓度标准限值。所有水质指标的CWQI值相加即为该水样的CWQI。

《城市地表水环境质量排名技术规定(试行)》推荐CWQI计算的基本水质指标为GB 3838—2002《地表水环境质量标准》表1中除水温、粪大肠菌群和总氮以外的21项指标。为优化监测项目，对中山市各河流监测断面2018—2019年的数据(共计水质监测样本2 669个)进行了分析，结果显示中山市各河流的主要超标因子为化学需氧量(COD)、总磷(TP)、氨氮(NH3-N)和溶解氧(DO)，这4项水质监测指标计算得到的CWQI值(I4)与21项水质监测指标计算得到的CWQI值(I21)之比的均值和标准差分别为53.7%和0.08%，相关系数r=0.94(P<0.05)，两者具有显著的相关性(图1)，表明I4可较好地反映河流水质状况，故本文中选取上述4项指标作为与光谱测量同步进行的水质监测指标，并计算I4值作为河流水质综合评价结果。

表1 样品水质特征

图1 4项和21项水质监测指标CWQI值的关系

各河流监测点的水样采集与地面高光谱测量同步进行，水质监测指标为COD、NH3-N、TP和DO。其中DO的质量浓度使用美国YSI公司Pro 2030多参数水质分析仪现场测定，其他指标取河流表层 0～20 cm处混合样，装入洁净聚乙烯塑料瓶中并加入相应保存剂后置于低温保存箱及时送至实验室分析。COD质量浓度采用重铬酸钾法测定，NH3-N质量浓度采用纳氏试剂光度法测定，TP质量浓度采用钼锑抗分光光度法测定[12]，并使用空白样、平行样、标准样品和加标回收等措施对样品监测全过程进行了质量控制。水质样品各特征指标见表1，可见多数样品水质属于Ⅳ或Ⅴ类，主要定类水质指标为TP和NH3-N。

1.3 PLSR反演模型建立与评价

由于水体反射率数据波段数众多，相邻波段间隔紧密，而不同波段光谱反射率数据间最大方差膨胀因子大于80 000，表明不同波段光谱反射率存在极严重的多重相关性[13]。为提高反演模型精度，采用PLSR建立水体不同波段光谱反射率与CWQI的反演模型。PLSR基于高维投影思想，通过提取主成分的方法对系统中的数据信息进行分解和筛选，实现对系统解释性最强的综合变量的有效抽取[13-15]；PLSR常用于样本量小、变量间存在多重相关性的数据处理，适用于中高光谱遥感数据反演模型的建立[16-17]。

PLSR建模的基本步骤[18]为：①标准化处理水体光谱反射率数据矩阵X和水体CWQI数据矩阵Y，得到的标准化矩阵分别记为E0和F0；②分别提取E0和F0第一对成分T1和U1，T1和U1分别为E0和F0的线性组合并使之相关性最大；③分别求E0和F0在T1上的回归，并得到回归残差矩阵E1和F1；④用残差矩阵E1和F1代替E0和F0重复上述步骤；⑤检查收敛性以确定提取的主成分数。

将151个样本原始数据随机划分为建模集(120个样本)和验证集(31个样本)两部分，其中建模集用建立PLSR反演模型；验证集用于反演模型对CWQI反演效果的验证。模型验证效果分别采用相关系数(r)、均方根误差(RMSE)、平均相对误差(MRE)进行评价。

2 结果与分析

2.1 光谱反射率与CWQI的关系

根据CWQI计算结果，151个水质样品的CWQI值范围为1.36～10.75，平均值为3.60，标准差为1.91。将不同CWQI值范围的水样进行分类并计算其高光谱反射率均值，结果(图2)显示，在400～700 nm波段范围内光谱反射率对不同CWQI值的水样有较好的区分度，水体光谱反射率均值随CWQI的增大依次递减，其中I4≤2的水体反射率均值在整个波段范围内均高于其他类别的水体。但在 720～900 nm波段范围内，出现I4>5的较差水质水体的光谱分辨率高于2

图2 不同CWQI范围水样对应的高光谱反射率均值

2.2 不同光谱分辨率对PLSR反演模型的影响

为分析光谱分辨率对水体PLSR反演模型估算效果的影响，需找出水体CWQI高光谱遥感监测最佳光谱分辨率。采用分段三次Hermite插值多项式法对原始数据进行重采样，采样间隔依次为约 0.20 nm、1 nm、10 nm、20 nm、30 nm、40 nm、50 nm、60 nm、70 nm和80 nm，连同两端波段的总波段数分别为2 534、900、51、26、18、14、11、10、9和8。在此基础上建立不同光谱分辨率样本集的PLSR反演模型，并根据相应的建模与验证精度评价结果研究模型估算精度随光谱分辨率的变化。结果(表2)表明，在各光谱分辨率的数据集中，原始光谱数据集的PLSR反演模型的估算效果最差；对原始数据进行不同光谱分辨率的重采样后，均有助于提升PLSR反演模型的估算效果。PLSR反演模型的验证集RMSE随着光谱分辨率的降低呈现先减小后增大的趋势。光谱分辨率为50 nm、波段数为11个时，验证集RMSE最小达0.790，MRE为18.63%，模型估算效果最优。但总体上重采样后各验证集的PLSR反演模型估算效果均处于接近的水平。该现象表明，在CWQI的遥感监测中，过高的光谱分辨率可能会造成数据冗余和引入环境噪声，造成模型复杂化并影响模型稳定性。反之，在较低的光谱分辨率下，数量较少的波段可在基本涵盖水体主要污染物特征谱段的同时，增强光谱信噪比，从而提高估算精度。

表2 不同波段数的PLSR建模效果评价

2.3 提取的主成分数对PLSR反演模型的影响

确定提取的主成分数是PLSR反演模型建立和优化的另一重要环节。目前PLSR建模确定主成分数的方法主要有舍一交叉验证方法、分批交叉验证方法、分裂样本交叉验证方法、随机样本交叉验证方法等[17]，其中舍一交叉验证方法的应用较为广泛，但该方法得到的建模主成分数可能不是最优值，主成分数过多或过少可能会导致欠拟合或过拟合，影响模型的估算准确度[23]。本文采用50 nm分辨率的水体光谱反射率数据，分别提取2～10个主成分进行PLSR建模，以验证集RMSE最小时对应的主成分数为最佳的主成分数。图3表明，模型的RMSE随主成分数的增加呈先下降后增加的趋势，当主成分数为8时，模型的RMSE达到最低值。

图3 不同主成分数与验证集RMSE的关系

因此，确定建模最佳主成分数为8，并建立PLSR反演模型如下：

I=4.84-23.33R400+22.62R450+53.53R500+

4.69R550-24.19R600-102.32R650+

75.96R700+25.38R750-30.93R800+

46.15R850-19.99R900

(4)

式中Rj为不同波长处的水体光谱反射率。

2.4 PLSR反演模型估算效果评价

以31个验证集样本的CWQI实测值为纵坐标、CWQI高光谱反演模型(式(4))计算值为横坐标，绘制其实测值与计算值的对比散点图(图4)。结果表明，该反演模型的估算精度良好，RMSE、MRE分别为0.768和18.1%，说明该反演模型在有效缩减了高光谱输入参数的前提下，仍然保持了良好的泛化能力和反演精度。

图4 验证集样本CWQI实测值与反演模型计算值比较

2.5 反演结果

由于无人机航拍覆盖范围较大，本文未进行同步的水质分析，主要对基于无人机高光谱影像的CWQI反演效果做初步定性分析，并重点关注河流水质突变区域CWQI反演结果的合理性。从应用结果(图5)上看，反演得到的CWQI值可较合理地反映出河流水质的空间差异。如图5(a)反演结果突出了地理位置接近，但因水闸控制导致连通性差的水体水质的显著差异；图5(b)反演结果反映了河流交汇处水体混合过程中的局部水质差异；图5(c)反演结果突出了码头船舶作业引起的水体扰动和水质变化。该结果表明，基于地面高光谱测点数据建立的CWQI反演模型应用于空间影像数据反演CWQI具有较好的可行性。

(a) 某水闸两侧水体

受限于无人机航程的限制，每次作业时间约在20 min左右，对于河流水质大范围监测的能力较为不足，仅能重点开展局部区域的水质高光谱反演。同时，现场航飞采用多航带作业方式，但对于水体流动性强，涨落潮剧烈的水域，水体流动导致水质的时空变化不能被无人机同步采集；这种准同步的水体光谱数据与实际河流的水质分布会有一定差异。但整体上反演模型应用于无人机高光谱影像的CWQI反演可取得良好效果，表明该方法具有良好的业务化运行潜力。

3 结 语

2019年7—9月在珠三角地区中山市选择代表性河流进行了水体高光谱反射率监测和同步水质监测，并以基于COD、TP、NH3-N和DO 4项水质指标计算得到的CWQI来表征河流水质状况。CWQI与水体高光谱反射率的相关分析表明，不同波段高光谱反射率对水质状况有较好的光学敏感性；基于PLSR建立高光谱数据与CWQI的最优反演模型，其验证数据集RMSE和MRE分别为0.768和18.1%，具有良好的反演精度。将该模型应用于无人机高光谱影像的CWQI反演，可以较合理地反映出河流水质的空间差异。尽管无人机航程、水流变化复杂等因素可能对无人机高光谱水质反演产生影响，但水质高光谱反演模型和无人机结合，可为河流水质综合评价提供一种大范围覆盖、可视化表达的技术手段，具有良好的业务化运行潜力。