胶州湾水体悬浮物浓度遥感反演模式优化研究❋

侯琳琳， 马安青， 胡 娟， 单姽月， 邓洁敏， 韩嘉仪， 丁子彧

(中国海洋大学环境科学与工程学院，山东 青岛 266100)

胶州湾属于典型二类水体，承担着青岛及周边地区工农用水，也为其经济发展做支撑，因此对胶州湾水质状况的监测十分重要。悬浮物浓度是水质和水环境评价的重要参数之一，影响水体感观，降低水体透明度和富氧条件，直接影响水体的生态环境，是决定水质的一个主要因素。传统定点取样方法费用高，难度大，代表性差(以点代面)，且受水文、气候条件影响无法保证样点的连续性，很难反映整个海湾的实时动态变化[1]。卫星遥感技术具有超高分辨率、连续观测等优势，且监测面积大、周期性好、速度快、成本低，便于长期动态监测，可弥补传统水质监测的缺陷，有效解决水质监测中数据采集的局限性问题。有效地监测表面水质参数时间和空间上的动态变化对水质监测意义重大[2]。

目前对近岸二类水体水色物质浓度遥感反演主要利用经验、理论和半分析模型[3-4]。Pavelsk等[5]利用SPOT和ASTER影像数据，对Peace-Athabasca Delta水域悬浮泥沙浓度进行研究，揭示了三角洲流域形态和水源的强烈差异，并通过276个MODIS影像监测，追踪洪泛区湖泊的水文补给情况，提出不同遥感方法对悬浮沉积物的可见/近红外遥感的强大价值，以评估复杂流动环境中的水文过程。建立经验及半经验模型反演水体中悬浮物浓度是目前悬浮物浓度遥感反演估算的主要方法。但是大量研究表明，胶州湾水体浑浊，各种半经验模型的反演精度和模型适用性有较大差异，目前二类水体的遥感反演监测方法并不适用于这一独特的水域[6]，因此本文针对胶州湾地区探讨独特并适用于半封闭海域的水色遥感监测方法，促进遥感反演具有一定意义。

1 材料与方法

1.1 研究区域概况

胶州湾位于120°15′E～120°24′E，36°00′N～36°12′N，南北长约33 km，口窄内宽，口门最窄处仅3.1 km，东西宽约28 km，面积446 km2。近年来受近岸河口泥沙汇入以及人类活动等因素影响，胶州湾承纳了大量市区与郊区水体污染物，水环境质量遭到破坏，悬浮物、叶绿素a含量、黄色物质等水体成分物理量增加，水环境问题突出。

1.2 数据的获取与处理

1.2.1 胶州湾采样点数据的获取与处理 为获取胶州湾水域采样点数据，本研究在胶州湾选取32个采样点(见图1)。严格遵循2009年颁布实施的国家环境保护标准的水质采样技术要求进行水体采样，根据胶州湾水体功能、水文要素和污染源及污染物排放情况，确定采样频次为一年4次，采样时间分别在每年的2、5、8和11月9～14时，天气晴朗，云层较少，采样深度为水面下10 m以内，采集瞬时水样，在水样采入或装入容器中后，立即做好避光保存，并做好GPS定位、气温、光照等记录。悬浮物浓度的获得均采用烘干称重法(GB11901～89标准)[7]，利用孔径为0.45 μm醋酸滤膜放入恒重称量模盒，移入烘箱中放入温度为40～80 ℃的烘箱中烘干8 h后移入干燥中，使其冷却至室温，称重，反复烘干、冷却、称重直至两次重量差小于0.4 mg，并计算悬浮物浓度。

图1 胶州湾实测采样数据站点图Fig.1 Measured sampling data site in Jiaozhou Bay

1.2.2 遥感数据获取与处理 在获取了地面实测数据的同时，同步或准同步选取并下载与实测数据日期最接近的2月18—25日，5月1—10日，8月2—13日，11月6—14的遥感卫星数据。获取2001—2015年2、5和11月的Landsat TM/ETM遥感影像(免费下载地址http://glovis.usgs.gov/)和MERIS遥感图像(下载地址http：//lsatellite.cma.gov.cn/)，同时通过国家海洋一所获得2012—2013年HJ-1A/B CCD、HSI卫星数据。通过选取Landsat TM/ETM+卫星图像代表胶州湾4个季节的悬浮物的特点，卫星过境时间分别为2、5、8和11月。遥感影像受光谱、云层的影响与真实性之间存在误差，需对图像进行校正，校正前首先对TM/ETM图像进行辐射定标。目前相对完善的发展模型是6S大气辐射校正模型，进行3×3的低通滤波处理，对于可见光-近红外波段数据分析十分适用，该模型的计算精度高，且可靠、快速[8]。

图2 原始影像(左)与大气校正后(右)影像对比Fig.2 The original image (left) contrasts with the atmospheric corrected (right) image

1.2.3 实测光谱数据获取与归一化处理 光谱测量采用美国ASD公司生产的Field spec Hand Held手持式地物光谱仪，其波长范围为325～1 075 nm，共512个波段。采样日期应尽量选择当天天气晴朗、无云无云和水面平静。测量水体之前先用光谱仪测一下标准参考板，然后距离水面1 m处垂直水面进行测量，每个采样点至少测量10次光谱，剔除由于太阳直射对水面反射造成的数值较高的曲线，保留数值较低的曲线进行平均，作为该点的光谱反射值。在测量光谱时，面向太阳，避免仪器阴影影响水体的反射光谱[9]。为了更好地构建胶州湾的反演模型，选取敏感波段，本文将对实测悬浮物浓度和同步水体光谱曲线之间进行相关性分析，400～900 nm波长范围的归一化实测水体光谱曲线与悬浮物浓度的相关系数[10]如图3所示。

由图3可以看出，2月份，波长小于600 nm可见光区，水体反射率与悬浮物浓度呈负相关。波长大于600 nm短波红外波段水体反射率与悬浮物浓度呈现正相关，并且波长在670、710和830 nm中红外和近红外以及红外波段附分别出现峰值，且最高相关系数大于0.8， 这与经验研究的悬浮物光谱反射曲线两个峰值之一近红外波段(760～900 nm)相吻合，且峰值在830 nm处出现向长波段方向红移现象，说明在该波段附近水体反射率与悬浮物浓度存在很好的相关性。5月，在波长为575 nm黄光波段附近时出现最高负相关，这与悬浮物浓度较高情况下光谱峰值出现在560～590 nm较吻合。而较高正相关出现在760 nm之后，表明在这些波段水体的遥感反射率与悬浮物浓度的相关性较大。8月，水体反射率与悬浮物浓度在波长小于580 nm时呈现负相关，在500 nm附近有较高的负相关，在625和700 nm附近均出现了较高正相关，这一时期的属于夏季，胶州湾水体富营养化较严重，光谱与叶绿素浓度值相关性较好，对悬浮物敏感波段选取没有很大意义。11月，总体相关性与2、5和8月均不同，在570 nm波长之前呈正相关，在425 nm附近的正相关系数最大，而波长在570 nm之后水体反射率与悬浮物浓度开始呈负相关，虽负相关系数整体较小，最大负相关系数也不到0.6，但在600、770和860 nm波段处均出现小峰值，可以认为该光谱波段反射率与悬浮泥沙浓度相关性较好。针对每个月份最大正相关系数与负相关系数的排序结果见表1。

图3 2015年2、5、8和11月份悬浮物浓度与归一化实测光谱曲线的相关系数Fig.3 Correlation coefficient of suspended matter concentration and normalized measured spectral curves for February, May,August and November， 2015

表1 不同月份相关系数峰值分布Table 1 Correlation of peaks in different months

由表1可知，2月分别有两个较高正相关波段(704和657 nm)，5月在588 nm附近，水体反射率与悬浮物浓度呈现较高负相关，波长大于700 nm出现了较高的正相关，在波段890 nm左右出现最大正相关。11月，波长在700 nm附近水体遥感反射率与悬浮物浓度呈现最高负相关。因此，通过这些较高正相关和负相关波段可以选取反演胶州湾悬浮物的最优敏感波段，分别为黄光波段(560～590 nm)，近红外波段(760～900 nm)以及相对稳定的波段系数800～820 nm进行反演[11]。

2 模型建立与分析

2.1 基于HJ-1 的生物光学模型的建立与应用

基于环境一号卫星数据进行反演，本研究采用6S大气校正模型，在对获取的HJ-1A多光谱影像数据辐射定标之后，利用ENVI软件中的EFFORT Polishing的工具对表观反射率进行了适度修正，使光谱曲线更加接近于真实地物光谱曲线。通过传感器不同波段间比值分析，有效去除大气的部分影响，消除不同采样点在时空上的差异以及水面粗糙度变化的干扰，因此采用归一化处理的水体反射率和悬浮物浓度较大正负相关的波段之间做比值分析[12-14]。根据1.2分析，选取2012和2013年4个月份数据并基于HJ-1B卫星进行反演，为了更好地模拟悬浮物浓度值，2月份选择模型CCD3/CCD2，5月份选择模型CCD4/CCD2，8月份选择模型CCD3/CCD1，11月份选择模型CCD1/CCD4；从2013年的野外实测悬浮物浓度中随机选取22个样点建立模型，留出10个点用来检验模型的精度，把较大正相关与较大负相关的比值作为自变量，将实测悬浮物浓度值作为因变量，建立模型进行回归分析，针对每一种波段算法分别建立五种数学模型相关性分析。选取数学模型分别为指数函数、对数函数、幂函数、线性函数和一元二次方程算法分析(见表2)。

基于HJ-1B卫星数据，2月CCD3/CCD2的一元二次多项式能较好的估测悬浮物浓度，最高相关系数R2=0.62，RMSE=1.4 mg/m3；5月CCD4/CCD2的一元二次多项式模型相关系数较好，最好值为R2=0.68，RMSE=2.1 mg/m3；8月CCD4/CCD2的一元二次多项式模型相关系数和绝对均方根误差结果都较好，最高值R2=0.77，RMSE=1.2 mg/m3；11月CCD4/CCD2的一元二次多项式模型相关系数和绝对均方根误差都好，其中R2=0.74，RMSE=1.1 mg/m3。分析表中模型的均方根误差(RMSE)均远小于胶州湾实测平均悬浮物浓度，但是R2所达到的值都不太高，根据验证结果，这个模型不能够有效地定量反演胶州湾水体悬浮物浓度。

2.2 基于MERIS的半分析模型的建立与应用

目前国际上常用MERIS半分析算法反演水体水质[15-16]。水体中各组成成分浓度不同会导致水体反射率的不同变化，传感器接收到的信号也随之变化，即水体中组成成分浓度与遥感反射率相关，而遥感反射率与水体总吸收系数和后向散射系数相关。选取三波段模型和两波段比值模型应用于胶州湾悬浮物遥感反演研究。基于2014年11月、2015年2和5月MERIS遥感影像数据和同步实测悬浮物浓度数据及光谱数据进行分析，分别通过两波段比值算法和三波段比值算法，建立线性、指数、幂函数、一元二次多项式和对数函数分析模型，比较分析相关性系数R2和均方跟误差RMSE来确定模型精度，得到如表3的结果。

表2 基于实体反射率光谱曲线的悬浮物反演模型Table 2 Suspended matter inversion model based on solid reflectance spectral curve

如表3所示，基于MERIS卫星数据，2015年2月三波段算法的一元二次多项式能较好的估测悬浮物浓度，相关系数最高R2=0.53，RMSE=1.7 mg/m3；2015年5月两波段算法的幂函数模型相关系数较好，为最高R2=0.66，RMSE=1.9 mg/m3；2014年11月三波段算法的一元二次多项式相关系数和绝对均方根误差结果都较好，最高值R2=0.67，RMSE=2.0 mg/m3。但是相关系数R2整体精度值较低，不能很好适用于胶州湾悬浮物遥感反演，该模型精度不高因此不能达到较好的对胶州湾悬浮物进行遥感反演的预期结果，但为其他研究者做相似研究可以提供经验借鉴。

2.3 基于Landsat波段组合回归模型的建立与应用

上述半经验半分析算法适用于Landsat TM/ETM卫星数据，将半经验半分析模型的波段根据使用的Landsat进行调整，据此调整了二波段模型、三波段模型[17-20]。

由于到中红外波段之后，水体将大部分的光吸收，不能反映出水体组分信息，因此在本研究中对Landsat TM/ETM的前4个波段进行模型的构建，调整后的二波段模型(Band 4/ Band 1，Band 4/ Band 2，Band 4/ Band 3)、三波段模型(Band 4、Band 3、Band 2；Band 3、Band 2、Band 1)，与此同时又针对每一个模型采用多个函数分析法，分别为指数算法、线性算法、对数算法、一元二次算法。每一种算法使用22组样本点进行建模，预留10组样本点进行敏感性精度检验。用于模型的22组样本点和检验模型的10组样本点都均匀覆盖在整个胶州湾水域，具有代表性可推算模型的精度(见图4)。

表3 基于MERIS的反演模型和精度检验Table 3 Inversion model based on MERIS and accuracy test

图4 22组回归点和10组验证点散布整个区域Fig.4 22 sets of regression points and 10 sets of validation points spread the whole area

(1)TM数据反演波段模型：2001—2006年的遥感影像与实测数据的选取按照同步方法，5月14日Landsat遥感影像，悬浮物含量的实测数据使用的是5月6—20日的测量数据； 2月10日的Landsat遥感影像，悬浮物含量的实测数据使用的是2月8—18日的测量数据， 8月12日的Landsat遥感影像，悬浮物含量的实测数据使用的是8月2—13日的测量数据； 11月13日的Landsat图像，悬浮物含量的实测数据使用的是11月6—14日的测量数据。32个样点在胶州湾分布均匀，代表整个区域的悬浮物分布。在建立模型中，使用指数函数、对数函数、幂函数和一元二次方程来分析相关性(由y表示)，选取最优函数变换结果。为了找出比较合适的波段或波段组合以估测悬浮物含量，根据敏感波段分析，并参照前人的研究成果来确定了5种波段组合形式算法，并计算各个组合与实测悬浮物浓度的相关系数，结果如表4。

从表4中可以看出： Landsat TM图像经过校正处理后所建的波段模型相关性系数较好，且(1/RTM2-1/RTM3)·RTM4所代表的432波段与悬浮物浓度的相关系数最高，反演精度最高，通过y(多项式分析)达到的相关系数为0.79。除RTM4/RTM3使用y(指数函数)分析有着较好的相关性，其余波段模型均以多项式分析取得较好的相关性系数，这与光洁等[20]的研究比较一致。同时本文发现，RTM4/RTM1与RTM4/RTM3同样适合于不同水质情况下的悬浮物估测。实测水体光谱分析表明，估测悬浮物较好的波段位于TM的第2、3、4波段，表4中TM432波段组合相关性最好正和此相吻合，在其表达式中，Band 2、 Band 4(以下简称B)充分考虑不同浓度悬浮物的反射峰的特点，可有效反映出悬浮物的光谱特征；1/B4～1/B2是为了减少叶绿素对低浓度悬浮物遥感信息的干扰作用。基于悬浮物的低浓度和较高的叶绿素a浓度，因B2处于叶绿素的反射峰，故1/B2将减小，又因B4处于叶绿素的吸收峰值，故1/B4将增大，1/B4～1/B2相应减小，从而起到抑制叶绿素对悬浮物遥感信息干扰的作用较小。因此，可以认为(1/B2～1/B3)·B4是较好的估测悬浮物含量的因子。

表4 基于Landsat TM图像建立的模型形式和反演精度值Table 4 Model form and inversion precision values based on Landsat TM image

注：RTMi代表的是Landsat TM波段的中心波长反射率，y代表波段函数回归分析波段模型。

Note:RTMirepresents the central wavelength reflectance of the Landsat TM band.yrepresents the band function regression analysis band model.

(2)ETM反演波段模型：根据以上分析结果，由于B4与B3、B2有着很强的相关性，提取各期处理后的Landsat ETM+图像遥感数据，与同步实测悬浮物浓度值分别建立相应的悬浮物遥感估测模型。各期悬浮物遥感估测模型的形式和精度见表5，模型预测值和实测值的对比情况如图5。

表5 基于Landsat ETM图像建立的模型形式和反演精度值Table 5 Model form and inversion precision values based on Landsat ETM image

注：RETMi代表的是Landsat ETM波段的中心波长反射率，y代表波段函数回归分析波段模型。

Note:RETMirepresents the central wavelength reflectance of the Landsat TM band,yrepresents the band function regression analysis band model.

从图5中可以看出，校正后的ETM图像反映在432波段模型算法的相关性系数最好(高达0.88)，预测值与实测值相关性较好，总体来看，上述反演模型结果较好。以上建立的悬浮物遥感估测模型中，432波段模型算法是最好的，其它模型相关系数相对较低。一方面是由于叶绿素和黄色物质对悬浮物浓度估测的影响，时间上同步性也较差，2003年11月遥感图像与实测水质时间相差10 天，且水样主要集中于胶州湾北部，分布不均因而影响悬浮物估测模型精度值。另一方面，本研究采用的夏季遥感图像与水质实测时间的间隔较长，加上夏季春季雨水充沛，近岸水体交换强烈，对监测干扰大，对悬浮物的估测也存在一定的影响。

3 误差敏感性分析

考虑到现场观测仪器和测量数据的各种误差以及Landsat TM/ETM遥感数据误差对反演算法的影响，同时也为验证反演算法是否合理准确，本文对以上反演算法进行了误差敏感性分析。对预留10组实测点数据与对应算法估测出的悬浮物浓度值进行对比分析，以检验遥感反射率得出的估测值的可靠性，具体如图6所示。

图5 悬浮物浓度实测值与模型预测值散点图Fig.5 Scatter plot of measured values of suspended matter concentration and model prediction

图6 误差验证图Fig.6 Error verification diagram

如图6所示，针对所选的四种反演精度较高的算法，分别与对应的算法预留的10组实测数据和Landsat ETM 估测得出的数据进行相关性分析。其中估测数据是从图4反演数据总数中随机选取22组点与预留的10组实测验证点进行多项式线性比较分析。结果明显看出每一组的相关性达到0.83以上，尤其三波段算法的(1/RETM2-1/RETM3)·RETM4相关系数最高达到0.89。这和图4分析结果一致，表明Landsat ETM三波段反演算法较好，并且误差敏感性较小，用于胶州湾水体悬浮物浓度反演研究比较合理。

同时选取F检验法进行分析，通过计算实测反演值和误差分析值的方差之比来检验两组数据的差异。分别计算两组测定值的方差并进行统计量方差平方之比，结果F值为6.285。查F分布表判断出Fa(4，4)=6.388，即F

表6 胶州湾水体悬浮物浓度模型评价Table 6 Evaluation of the concentration model of suspended solids in Jiaozhou Bay

分析表6中模型的均方根误差(RMSE)小于胶州湾实测平均悬浮物浓度，R2也高达到0.89，根据验证结果，这个模型能够有效地定量反演胶州湾水体悬浮物浓度。由图6和表6分析可知，悬浮物预测模型的验证精度达到了水质参数的反演要求，对胶州湾悬浮物反演有一定的应用价值。

4 结论与讨论

本研究通过2001—2015年8年的数据并在分析前人的研究经验基础上，发现了一种较适合胶州湾悬浮物浓度估测的因子(1/RETM2-1/RETM3)·RETM4三波段反演算法。结果表明，悬浮物浓度的估测精度高，相关系数R2均在0.8以上，为胶州湾水体悬浮物遥感监测提供了一个较好的实用反演模型。

根据2000—2015年胶州湾水面实测数据分析，胶州湾水体因季节变化悬浮物有不同特征，因而对比分析了三波段模型、两波段模型两种方法选取最优模型算法估算胶州湾悬浮物浓度。基于MERIS、HJ-1和Landset TM/ETM+遥感数据与准同步地面采样数据，建立胶州湾悬浮物浓度反演模型并分析比较最优模型算法。结果表明利用Landsat ETM+遥感数据估算的三波段模型432算法与悬浮物浓度的相关性最好(R2>0.8)，反演精度较高，建立了悬浮物浓度的三波段半分析生物光学反演模型，并利用10组独立实验数据进行检验，结果发现算法反演值与实测值的R2均高于0.85，证明误差的敏感性差，该模型稳定性好。

由于本研究采用的实测样本数据并不是很多，模型的稳定性及代表性还需要进一步验证。在以后的工作中可以加大样点数量，提高样点的代表性。随着样本数据的增多，各模型中的参数趋于更加稳定，即可以在没有同步地面数据的情况下通过卫星遥感图像和模型来估测当时的悬浮物浓度。另外，虽然目前通过水样实验方法对单个水质参数进行研究，如叶绿素a、悬浮物、黄色物质等的光谱特性比较了解，但对水体组成成分之间的光学性质相互影响及各自光学特性分析了解相对不足，在今后仍需加强分析叶绿素a、黄色物质对悬浮物浓度遥感反演估测的影响。分离各水质参数物理量，明确它们之间的相互影响，才能更准确更有效的建立模型。

致谢：感谢国家海洋局第一海洋研究所崔廷伟、马毅和李晓敏老师提供的实测数据及测量光谱仪器。