基于高光谱的水体BOD含量模拟估算
王洪伟1, 王波2, 纪童3, 徐君4, 剧锋5, 王彩玲6,*
1.武警工程大学, 陕西 西安 710086
2.盐池县草原实验站, 宁夏 盐池 751506
3.甘肃农业大学草业学院, 甘肃 兰州 730070
4.西安航空学院, 陕西 西安 710077
5.中华人民共和国银川海关, 宁夏 银川 750000
6.西安石油大学, 陕西 西安 710065
*通讯作者 e-mail: azering@163.com

作者简介: 王洪伟, 1980年生, 武警工程大学副教授 e-mail: whwdyx@163.com

摘要

高光谱技术由于满足连续性与光谱可分性的要求, 具有能够区别同一种地物不同类别的能力, 且光谱数据获取速度快, 操作简易, 在监测水体分布状况、 水体指标上具有突出成就。 生化需氧量BOD是评价水污染的重要指标, 现行常规的测量方法为五日培养法, 这种方法消耗试剂、 操作复杂、 受干扰因素多、 测定时间长、 不能及时反映水质变化、 无法及时有效地预警突发水污染事件, 鉴于传统方法的缺点, 探索基于高光谱技术的水体BOD含量的估算和反演对水质评定具有重要意义。 以西安地区三处地表水为研究区, 共计60处试验点, 每处试验点重复测定10次光谱与BOD值, 取平均值作为原始光谱, 利用Person相关系数法筛选光谱与BOD值的敏感波段, 并使用主成分分析与最小二乘法消除光谱指标的多重共线性, 建立水质BOD指标的多元线性回归模型与偏最小二乘回归模型。 研究结果如下: (1)BOD敏感波段大体分布于600~900 nm, 共筛选出了35个显著相关的原始光谱指标, 其中758 nm相关系数绝对值最高(0.418); (2)经由主成分分析降维得出的 Z1 Z2与BOD指标的多元线性回归模型精度较好( R2=0.565, RMSE=0.007), 且主成分分析中可以明显区分0~0.2与0.4~0.6 mol·L-1 BOD浓度; (3)光谱指标与BOD指标构建偏最小二乘回归模型的精度 R2高达0.896, RMSEP=0.746 9(留一交叉法均方根误差); jack.test检验发现628 nm对反演水体BOD含量的影响极其显著, 889与893 nm波段对其影响较为显著; (4)根据模型拟合精度, 筛选的最优的BOD反演模型为偏最小二乘回归模型, 对偏最小二乘模型进行精度检验, 精度较好( R2=0.81)。 基于以上试验结果, 提出了一种基于偏最小二乘法高光谱水质BOD参数的反演方法, 为水质BOD参数动态检测提供了新方法。

关键词: 高光谱; BOD; 模型; 偏最小二乘法; 多元回归
中图分类号:P237 文献标志码:A
Simulation Estimation of BOD Content in Water Based on Hyperspectra
WANG Hong-wei1, WANG Bo2, JI Tong3, XU Jun4, JU Feng5, WANG Cai-ling6,*
1. Engineering University of CAPF, Xi’an 710086, China
2. Grassland Experiment Station of Yanchi, Yanchi 751506, China
3. College of Grass Industry, Gansu Agricultural University, Lanzhou 730070, China
4. Xi’an Aeronautical University, Xi’an 710077, China
5. Yinchuan Customs District P. R. China, Yinchuan 750000, China
6. Xi’an Shiyou University, Xi’an 710065, China
*Corresponding author
Abstract

Due to the requirement of continuity and spectral separability, hyperspectral technology has the ability to distinguish different types of the same ground object, and the spectral data acquisition speed is fast, and the operation is simple. Spectral analysis has made outstanding achievements in monitoring water distribution and water indicators. Biochemical oxygen demandis one of the important indicators to evaluate water pollution, the current conventional measuring method for 5 culture method, and this method consumes reagent, complicated operation, more interference factors, determination of time is long, can not reflect the water quality changes in time, can’t early warning of emergent water pollution events in a timely and effective manner, in view of the traditional methods of faults, explore the content of water, BOD estimation based on the technology of hyperspectral and inversion for water quality assessment is of great significance. This test three surface water in xi ’an area as the research area, a total of 60 sites, each site repeat 10 times spectra and the BOD value, average as an original spectrum and the BOD value, Person correlation coefficient method is used to filter the spectrum and the BOD value of sensitive wavebands, and principal component analysis and least square method are used to eliminate spectral index of multicollinearity, BOD water quality index of the multivariate linear regression model and partial least squares regression model. The results were as follows: (1) the BOD sensitive bands were generally distributed at 600~900 nm, and a total of 35 original spectral indicators with significant correlation were screened out, of which the absolute value of the correlation coefficient of 758 nm was the highest (0.418). (2) the accuracy of multiple linear regression model of Z1, Z2 and BOD indexes obtained by principal component analysis ( R2=0.565, RMSE=0.007) is good, and the BOD concentration of 0~0.2 and 0.4~0.6 mol·L-1 can be clearly distinguished in the principal component analysis. (3) partial least-squares regression between spectral index and BOD index shows that the model accuracy R2 of the partial least-squares regression model is up to 0.896, RMSEP=0.746 9 (root mean square error with one crossing method). By jack test, it is found that 628 nm has a very significant influence on the BOD content of inversion water body, and the bands of 889 and 893 nm have a significant influence on it. (4) according to the model fitting accuracy, the selected optimal BOD inversion model is the partial least squares regression model, and the accuracy of the partial least squares model is verified to be good ( R2=0.81). Based on the above test results, an inversion method based on partial least squares hyperspectral BOD parameters of water quality is proposed, which provides a new method for dynamic detection of water quality BOD parameters.

Keyword: Hyperspectral; BOD; Model; Partial least squares method; Multiple regression
引言

随着人类物质生活水平的提高和工业化的发展, 水污染已经成为当今社会普遍存在的问题, 其监测与治理也备受关注。 生化需氧量(biochemical oxygen demand, BOD)是水体中的好氧微生物在一定温度下将水中有机物分解成无机质, 这一特定时间内的氧化过程中所需要的溶解氧量, 是监测水中有机物染物的一个综合指标[1], 是地表水、 生活污水及绝大多数工业废水的必测指标之一。 BOD值越高表明水中溶解氧会被自身微生物消耗的数值越高, 造成许多的生态问题[2]。 “ 五日培养法” 为现下普遍的接受测定BOD的方法, 但测定时间长、 不能及时反映水质变化, 不适合现场监测。

自20世纪70年代以来, 随着遥感技术的快速发展, 高光谱技术已成为现代遥感技术的重要组成部分[3], 利用高光谱技术反演水质指数早有研究, 刘彦君等[4]利用多光谱数据, 进行线性与非线性模型反演研究, 对浙江农林大学东湖水体的总磷(TP)、 浊度(SS)、 悬浮物浓度(TUB)进行了反演。 林剑远等[5]利用水质化验数据和光谱反射率进行相关性分析, 建立了浙江省嘉兴市河网化学需氧量(CODcr)、 生化需氧量(BOD5)、 总磷(TP)、 总氮(TN)的反演模型。 周亚东等[6]利用GF-1号WFV遥感影像, 通过多元线性回归和RBF神经网络模型建立了武汉市周围水域综合营养状态指数模型。 这些成果有效解读了水体光谱特征规律, 为遥感监测水质, 生产生活提供了理论支撑与技术指导。

原始光谱反射数据有着数据量大, 指标彼此高度相关的特性; 原始指标高度相关的特性经常会导致多重共线性问题的产生, 从而导致模型失真[7], 因此如何对大量光谱数据进行处理和挑选一直是光谱反演模型的重点。 主成分分析法(PCA)与偏最小二乘法(PLS)作为常用降维方法在遥感上应用广泛[8], 许多研究结果也表明应用主成分分析与偏最小二乘法筛选的主成分参数可以更好的反演各自的指标。 杨国范等[9]利用比值线性回归模型与最小二乘支持向量机, 对铁岭清水河库叶绿素a浓度与Landast OLI卫星数据分析, 并建立了叶绿素浓度a的反演模型。 何金成等[10]利用近红外光谱数据结合偏最小二乘法回归建立了BOD预测模型。

现有文献报道中, 利用光谱估测水质参数BOD指标的报道较少, 基于此试验利用光谱数据进行水体指标BOD的反演, 测定水体样本光谱数据的同时收集水体样本并带回实验室测定BOD指标, 将采集到的光谱数据与BOD指标进行Person相关性分析, 挑选敏感光谱指标; 由于光谱指标之间的高度相关, 为避免模型失真, 在建立反演模型之前, 利用主成分分析和偏最小二乘法分别对光谱指标进行处理, 消除指标之间的多重共线性问题, 最终建立多元线性回归模型与偏最小二乘模型, 比较两种建模方法的建模精度与预测效果, 选出更加适合反演BOD指标的建模方法。 探索利用高光谱技术估测水体BOD值的可行性与最优方法, 为实时诊断水体状况提供理论基础和关键技术, 为实现对BOD指标实时监测提供可行的途径。

1 实验部分
1.1 试验地概况

于2018年对西安市地表水环境进行取样研究, 取样地点集中于渭河(林家村)、 浐河(田家湾)、 灞河(马渡王), 共计60处采样点, 每处采样点共计10次重复。

1.2 方法

1.2.1 光谱数据获取与校正

所用仪器为美国ASD (Analytica Spectra Devices., Inc)公司制造的适用于遥感测量、 农作物监测等方面的 FieldSpec® 4 Hi-ResASD便携式地物光谱仪, 其光谱范围为300~2 500 nm。

光谱采集选择干燥、 无风、 晴朗无云或少云的天气进行, 并根据天气条件及时进行标准白板校正, 采集时间尽量在10:00— 14:00之间, 此时光照条件良好。 进行地面水质采样和水体光谱数据等实验数据获取, 光谱采集参数设置时间为100 ms, 测量后及时进行白板校正[11]。 每块样本选择2~3个光谱采样点进行高光谱数据采集, 每个样点每次重复测量10次, 最后以该样点的光谱反射率均值制作光谱反射率曲线。

1.2.2 BOD指标的测定

采集水体样本时, 利用聚乙烯桶采集距离水面10~12 cm的水样, 不使漂浮于水面的物质混入, 每处试验点共取10次样本, 对水样加入保存剂, 以便将样本带回实验室, 利用标准稀释法[12]处理水样, 并在20 ℃培养箱中培养, 5 d后测出培养后的溶解氧含量, 取平均值作为BOD指标原始因变量。 样本BOD参数变化范围如表1所示。

表1 水质参数变化范围 Table 1 Variation range of water quality parameters
1.3 数据处理

普通的多元线性回归应用中有许多限制, 最典型的问题就是自变量之间的多重相关性。 为此, 利用主成分分析降维与消除指标间多重共线性的特性, 筛选多元线性模型的自变量, 已期解决多重共线性对参数估计的影响, 减小模型误差。 偏最小二乘回归中开辟了一种有效的技术途径, 通过对系统中的数据信息进行分解和筛选, 提取对因变量的解释性最强的综合变量, 辨识系统中的信息与噪声, 从而能够更好地克服变量多重相关性在系统建模中的不良作用。

2 结果与讨论
2.1 Person相关性

图1为原始光谱DN值与水体BOD含量的相关系数图, 因图中波段1 023~2 500 nm与水体BOD含量无显著相关性, 因此图中只展示了350~1 023 nm波段范围内的相关系数, 由图可知BOD指标与光谱在350~900 nm呈负相关, 960~100 nm为正相关, 350~490与920~1 000 nm与BOD指标无显著相关性, BOD敏感波段大体分布于600~900 nm, 其中758 nm处为相关系数绝对值最大值0.418, 根据相关系数大小与显著性原则, 筛选出了35个与BOD指标极显著相关的原始光谱指标, 作为多元线性回归模型与偏最小二乘模型的自变量, 筛选指标相关系数绝对值由大到小分别为: 758, 759, 853, 809, 1000, 810, 890, 813, 851, 1 012, 807, 893, 618, 864, 816, 806, 782, 787, 785, 888, 796, 808, 924, 845, 663, 530, 887, 724, 863, 889, 757, 683, 628, 909和689 nm。 主成分分析要求建模数据量高于变量数, 偏最小二乘法允许在样本点个数少于变量个数的条件下进行回归建模, Person相关系数法共筛选出35个光谱变量, 因此将60组样本数据分为建模组(40)与检验组(20)。

图1 相关系数图
注: 图中蓝色波长为极显著相关、 绿色波长为显著相关、 红色波长无相关性
Fig.1 Correlation coefficient diagram
Note: The blue wavelength is extremely significant and the green wavelength is significantly correlated. The red wavelength has no correlation

2.2 主成分分析

主成分分析结果如图2所示。

图2 主成分分析碎石图
注: 横坐标是主成分, 纵坐标为解释程度
Fig.2 Principal component analysis lithotripsy
Note: The abscissa is the principal component and the ordinate is the degree of interpretation

经分析共有10个主成分。 其中主成分1方差贡献率为94.9%, 主成分2方差贡献率为1%, 而主成分3~10累积方差贡献率不足10%, 且主成分2到主成分3, 斜率开始趋于平缓, 因此剔除主成分3~10, 只保留主成分1和2(Z1和Z2), 这2个主成分既能达到降维的目的, 又能反映原始数据95.9%的信息。

图3直观展现了各植被指数在主成分1和主成分2中的分布情况。 横纵坐标分别代表第一主成分与第二主成分以及各自的贡献率, Z1Z2累计贡献率高达95.9%, 可以解释原有变量中的大部分信息, 4个BOD含量分组中0~0.2与0.4~0.6 mol· L-1在4组中彼此独立, 可以明显区分, 0.2~0.4与0.6~0.8 mol· L-1彼此交叉分组不明显。

图3 主成分分析效果图
注: 图中横纵坐标代表主成分1与主成分2的方差贡献率, 箭头代表原始变量(x1-x35代表按波长大小排列的原始光谱变量), 其中方向代表原始变量与主成分的相关性, 长度代表原始数据对主成分的贡献度
Fig.3 Principal component analysis renderings
Note: In the figure, the horizontal and vertical coordinates represent the variance contribution rate of the principal components 1 and 2 and the arrows represent the original variables (x1-x35 represent the original spectral variables arranged by the wavelength), where the direction represents the correlation between the original variables and the principal components, and the length represents the contribution of raw data to the principal components

将特征向量代入主成分公式中, 得到主成分Z1Z2的表达式

Z1=0.168x1+0.168x2+0.17x3+…+0.166x34+0.169x35

Z2=0.382x1+0.353x2+0.287x3+0.28x4+0.151x5+0.122x7-0.311x26-0.264x27-0.208x29-0.364x30-0.196x32-0.271x34-0.13x35

将2个主成分分别代入多元线性回归中, 得到的方程

YBOD=-0.000004468z1+0.00005919z2+9.217(R2=0.656, RMSE=0.007)

多元回归模型中BOD与主成分拟合方程R2较大, RMSE值较小, 说明利用主成分Z1Z2通过多元线性回归, 可以很好的拟合水体BOD指标。

2.3 偏最小二乘模型的构建

由于自变量与因变量之间的量纲与数值都是不同的, 现将BOD值与筛选的光谱指标进行标准化处理, 利用R语言PLS偏最小二乘函数包建立水质BOD含量的估测模型, 各主成分贡献率结果见表2

表2 主成分贡献率 Table 2 Contribution rate of principal component

表2可知, 当主成分为3时, 解释率逐渐趋于平稳, 因此选取comps=3时建立模型。 y=0.015 703x1+0.124 092x2+0.423 545x3-0.181 04x4-0.255 47x5+…+0.331 165x34-0.189x35(R2=0.896, RMSEP=0.7469)。

使用函数包中jack.test函数对回归系数进行显著性检验结果见表3

表3 jack.test函数显著性检验 Table 3 Significance test of jack.test function

通过jack.test函数进行显著性检验, 表3中“ * ” 代表极显著影响, “ * * ” 代表显著影响, x1x35代表波长按由小到大排列的原始光谱变量, 由表3可知对水体BOD含量有显著影响的光谱指标有628, 889和893 nm, 其中对BOD有正向影响的光谱指标为628与889 nm, 对BOD有负向影响的光谱指标为893 nm。

2.4 最优模型筛选

比较多元线性回归模型与偏最小二乘法模型, 依据R2最大RMSE最小原则, 最终采用偏最小二乘法模型y=0.015 703x1+0.124 092x2+0.423 545x3-0.181 04x4-0.255 47x5+…+0.331 165x34-0.189x35(R2=0.896, RMSEP=0.746 9)。

2.5 模型精度检验

将检验组的20组BOD与光谱数据代入模型进行模型精度检验, 检验结果见表4

表4 模型精度检验 Table 4 Model accuracy test

表4可以看出, 偏最小二乘模型, 其均方根误差较低为0.12, 且估测精度R2较高。 说明利用偏最小二乘法可以建立精度较好的BOD反演模型。

3 结论

利用多元线性回归与偏最小二乘法建立水质BOD指标的反演模型。 在进行光谱单波段与BOD指标相关性分析时发现在350~500 nm波段相关系数偏低, 350~500 nm原始光谱反射曲线杂乱, 说明该波段可能受其他水质参数影响, 波段敏感性较差, 不能作为模型预测波段。 随着波长增加相关系数于758 nm达到最高值, 且光谱最优反演波段大多分布在600~900 nm处, 与林剑远等[5]得到的高光谱遥感数据与BOD指标敏感波段750~900 nm有一定不同但也有相似之处, 主要原因有以下几点:

(1)水体光谱受自然条件与人为干扰, 使边缘波段噪声很大, 导致350~500 nm波段敏感性较差, 与BOD指标的相关性较低。

(2)水体光谱易受时间空间等影响, 导致光谱区别较大, 但光谱趋势整体相似, 且受其他水质指标影响, 其光谱也会随之变化。

高光谱具有分辨率高, 波段连续性强的特点, 但光谱信息冗杂, 数据的筛选与模型的简化一直是光谱模型研究的重点[13], 通过主成分分析和偏最小二乘法综合筛选的光谱指标, 建立了多元线性回归模型与偏最小二乘法, 结果表明主成分分析与偏最小二乘法可以有效降低数据维度, 综合筛选指标特性, 提高光谱数据与BOD参数的相关性与模型精度, 其中偏最小二乘模型模型精度远高于多元线性回归模型, 因为偏最小二乘法是分别从因变量与自变量中提取成分因子, 保证成分因子尽可能多的反应变量的变异信息, 同时也保证了两者之间相关性最大[14], 且试验中样本个数与变量个数大致一致, 适用于偏最小二乘法模型。

在拟合偏最小二乘模型时, 利用jack.test函数得出对水体BOD含量有显著影响的光谱指标有628, 889与893 nm, 说明628, 889与893 nm可以作为反演BOD指标的敏感波长。 林剑远等[5]以高光谱数据研究是城市河网BOD指标中发现波段565 nm为单波段与BOD指标相关系数(0.44)最佳波段, 689/667 nm为组合波段与BOD指标相关系数(0.84)最佳波段, 与本工作筛选的敏感波段有所不同, 但有所相近。

以上试验结果为水质BOD指标的快速估算提供了依据, 也为水体质量评估提供更便利的方案。

参考文献
[1] HU Guo-qiang, YANG Hui-ling(胡国强, 杨惠玲). Chinese Journal of Environmental Engineering(环境工程学报), 1985, 12: 8. [本文引用:1]
[2] ZHANG Jun, YANG Wen-wu(张钧, 杨文武). Environmental Science and Management(环境科学与管理), 2011, 36(7): 139. [本文引用:1]
[3] WANG Li-wen, WEI Ya-xing(王莉雯, 卫亚星). Spectroscopy and Spectral Analysis(光谱学与光谱分析), 2013, 33(10): 2823. [本文引用:1]
[4] LIU Yan-jun, XIA Kai, FENG Hai-lin, et al(刘彦君, 夏凯, 冯海林, ). Research of Environmental Science(环境科学研究), 2019, 39(4): 1241. [本文引用:1]
[5] LIN Jian-yuan, ZHANG Chang-xing(林剑远, 张长兴). Remote Sensing Information(遥感信息), 2019, 34(2): 23. [本文引用:3]
[6] ZHOU Ya-dong, HE Bao-yin, KOU Jie-feng, et al(周亚东, 何报寅, 寇杰锋, ). Resources and Environment in the Yangtze River Basin(长江流域资源与环境), 2018, 27(6): 1307. [本文引用:1]
[7] MAO Han-ping, GAO Hong-yan, ZHANG Xiao-dong(毛罕平, 高洪燕, 张晓东). Transactions of the Chinese Society for Agricultural Machinery(农业机械学报), 2011, 42(5): 166. [本文引用:1]
[8] CHU Xiao-li, XU Yu-peng, LU Wan-zhen(褚小立, 许育鹏, 陆婉珍). Modern Instruments and Medical(现代仪器与医疗), 2007, 13(5): 13. [本文引用:1]
[9] YANG Guo-fan, YAN Meng-dong, YIN Fei(杨国范, 阎孟冬, 殷飞). Remote Sensing Information(遥感信息), 2016, 31(5): 77. [本文引用:1]
[10] HE Jin-cheng, YANG Xiang-long, WANG Li-ren, et al(何金成, 杨祥龙, 王立人, ). Acta Scientiae Circumstantiae(环境科学学报), 2007, 27(12): 2105. [本文引用:1]
[11] YANG Wei, LI Min-zan, ZHENG Li-hua, et al(杨玮, 李民赞, 郑立华, ). Transactions of the Chinese Society of Agricultural Engineering(农业工程学报), 2015, 31(s2): 164. [本文引用:1]
[12] GAO Yu-feng, WANG Xiao-na(高宇峰, 王晓娜). Jilin Water Resources(吉林水利), 2018, 1: 22. [本文引用:1]
[13] YANG Ke-ming, SUN Yang-yang, LIU Fei, et al(杨可明, 孙阳阳, 刘飞, ). Science Technology and Engineering(科学技术与工程), 2015, 15(15): 82. [本文引用:1]
[14] FENG Lei, CHEN Shuang-shuang, FENG Bin, et al(冯雷, 陈双双, 冯斌, ). Transactions of the Chinese Society of Agricultural Engineering(农业工程学报), 2012, 28(1): 139. [本文引用:1]