近年来，受人类活动影响，内陆二类水体呈现不同程度的富营养化状态，严重影响了人类正常的生产生活^[1]. 叶绿素a浓度可以反映水中藻类物质含量，成为衡量湖泊富营养化的一个重要参数^{[2, 3]}. 国内外学者对复杂浑浊水体叶绿素a浓度反演已经展开了许多有价值的研究. 在具有代表性的3种算法^[4]：经验算法、 半经验算法及分析算法中，经验算法虽然形式简单，易于实现，但由于其物理基础薄弱，因此难以形成具有普适意义的模型，且精度一般较低^[5]； 分析模型具有坚实的物理基础，理论上一旦建立，即具有一定普适性^{[6, 7, 8]}，但其生物光学参数十分复杂，往往需要大量实验数据才能有效率定，不易实现； 半经验算法一般通过机理模型，或者水体反射率光谱特征，构建叶绿素a浓度的敏感指数，进而通过最小二乘拟合进行模型构建^{[9, 10, 11, 12, 13]}. 由于其模型形式简单，精度较高，又有较好的可解释性，成为应用较多的一种模型形式.

随着卫星遥感技术的不断发展，一些模型被成功应用于具有不同特性的卫星影像，实现了对内陆浑浊水体的遥感监测. 不同的卫星影像具有不同的特性和优势，如利用空间分辨率较高的TM数据^{[14, 15]}、 HJ 1-CCD数据^{[16, 17]}和RapidEye数据^[18]等构建估算模型，可以对叶绿素a浓度的空间细节信息进行更好的监测； 利用光谱信息较丰富的Hyperion数据^[19]，HICO数据^[20]和MERIS数据^{[21, 22]}，可以应用更加精确的半分析模型(如三波段算法^{[11, 12]}，以及据此发展而来的多种改进算法^{[13, 23, 24]})，实现叶绿素a浓度更加精准的观测； 由于内陆水体光学特性复杂，变化快，卫星影像的时间分辨率也逐渐成为研究内陆二类水体叶绿素a浓度高低变化，藻类生消周期的重要因素^[25].

GOCI作为世界上第一颗海洋水色观测静止卫星，其一天8次成像的时间分辨率为水体遥感动态监测提供了良好的契机. 其500 m的空间分辨率和平均20 nm的波宽，可以较好地捕捉水体的空间、 光谱特征. 然而由于其波段设置的缺陷，无法应用精度较高的三波段模型. 目前内陆二类水体叶绿素a浓度的GOCI估算模型，多为波段比值模型^{[26, 27]}，受建模数据集影响较大，精度难以保证. 针对GOCI数据的高精度叶绿素a浓度估算模型有待研究^[28].

因此，本研究基于DallOlmo等提出的标准三波段模型^{[11, 12]}，提出适用于GOCI影像数据的三波段算法，旨在拓展三波段算法的适用范围，提高GOCI数据在内陆二类水体中的叶绿素a浓度估算精度.

材料与数据 1 1.1 数据获取

本研究使用的野外实测数据采集自三峡水库、 巢湖、 洞庭湖和太湖. 包括三峡水库2009年8月1期数据，巢湖2009年6月1期数据，洞庭湖2013年8月1期数据，太湖2009年4月、 2010年8月、 2011年8月、 2013年5月和2013年8月5期数据. 共采集得到284组样本. 采样点分布如图 1所示.

图 1

Fig. 1

图 1 研究区示意 Fig. 1 Distribution of study areas and sampling stations

野外实验中对每个采样点进行了水体遥感反射率的测量，同时采集表层水样放入冷藏箱，并于当日送回实验室进行叶绿素a浓度及悬浮物浓度分析. 水体反射光谱的测量采用美国ASD公司生产的ASD FieldSpec Pro便携式光谱辐射计，其波段范围为350~1050 nm. 为减少水体镜面反射和船体自身阴影的影像，测量时采用唐军武等提出的 内陆二类水体水面以上光谱测量的方法. 具体方法见文献^[29]. 叶绿素a浓度的测定采用热乙醇法测量，即对采集的水样进行抽滤、 研磨、 离心、 定容后在分光光度计上分析得出该样点的叶绿素a浓度. 悬浮物浓度采用常规的干燥、 烘烧、 称重法(GB/T 11901-1989)测定. 数据统计信息如表 1所示. 其中统计了各期数据采集的地点、 样本数量，以及cChla和cs的均值、 方差、 最大最小值. 从中可以看出：首先，数据集样本覆盖范围较大，适合建模分析. 另外，三峡水库水质相对较好，叶绿素a浓度较低，太湖2009年4月叶绿素a浓度也较低. 除此之外，太湖8月样本和洞庭湖样本叶绿素a浓度动态范围较大； 从悬浮物浓度来看，三峡水库悬浮物浓度最高，其次为2009年4月太湖地区，其余日期的样本集动态范围和值域比较接近. 其中2013年5月的太湖数据用于GOCI影像的独立验证，不参与建模.

表 1(Table 1)

表 1 各个数据集采样点信息统计 Table 1 Statistical information of different datasets 获取时间 2009-04 2009-06 2009-08 2010-08 2011-08 2013-05 2013-08 2013-08 获取地点 太湖 巢湖 三峡 太湖 太湖 太湖 太湖 洞庭湖 点数 49 31 18 31 30 10 25 95 E(cChla)/μg·L-1 13.494 62.432 4.608 31.281 29.073 17.027 58.050 31.130 δ(cChla)/μg·L-1 14.804 43.166 5.906 25.155 8.190 7.639 45.516 30.125 max(cChla)/μg·L-1 79.796 192.882 18.860 110.592 52.109 34.317 229.303 124.248 min(cChla)/μg·L-1 2.975 23.815 0.039 2.745 14.046 5.122 3.720 1.860 E(cs)/mg·L-1 71.419 43.612 64.930 50.711 49.696 36.907 71.153 29.948 δ(cs)/mg·L-1 59.156 17.078 64.787 29.782 23.316 15.623 50.930 22.843 max(cs)/mg·L-1 244.900 83.933 253.000 132.333 109.266 79.200 215.750 101.000 min(cs)/mg·L-1 11.400 15.700 3.750 9.133 18.800 19.200 5.600 1.000 1)表中E代表均值，δ代表标准差

表 1 各个数据集采样点信息统计 Table 1 Statistical information of different datasets

除了野外实测数据，本研究另外使用了2013年5月13日的GOCI影像数据，用于模型的验证. 影像的几何校正通过自带的经纬度网格，在ENVI软件中进行校正. 大气校正使用6S模型(Second Simulation of a Satellite Signal in the Solar Spectrum)进行校正^{[30, 31]}.

Dall'Olmo等^{[11, 12]}提出了一种三波段概念模型来反演内陆水体的叶绿素a浓度，如式(1)所示.

该模型以生物光学模型为基础，与传统经验模型相比，有更严谨的推导和清晰的物理含义； 同时，模型中不涉及任何固有光学参数，因此比生物光学模型更容易应用和推广. 因此该模型在一些内陆湖泊和海湾中得到了广泛的应用^{[21, 22, 32, 33]}. 其推导过程如下.

根据水体生物光学模型，水体的光谱反射率与总吸收系数[a(λ)]和总后向散射系数[b_b(λ)]之间存在以下定量关系：

由于叶绿素a浓度与色素颗粒物吸收系数成正比，即Chla∝a_ph，因此只要将a_ph从a中剥离出来，使得式(2)的右侧只保留a_ph和常数，就可以反演出叶绿素a的浓度. 为了实现该目的，算法首先选取两个波段λ₁和λ₂，计算其反射率倒数之差，得到：

根据标准三波段模型原理，λ₁应位于a_ph峰值，λ₂与λ₁位置接近，且λ₂≥700 nm，λ₃位于近红外波段，λ₃≥730 nm. 然而，满足标准三波段算法3个波段设置的传感器较少，如图 2所示，高光谱数据，如Hyperion，在3个波长范围都有较多波段设置，可以满足需求； 而多光谱数据中，目前只有MERIS可以满足标准三波段模型的波段设置. GOCI数据均缺少λ₂处的波段设置，因此理论上无法使用标准三波段模型. 为了便于讨论，这里用3个波段处的中心波长来表达该模型，即可将式(8)表达为：

图 2

Fig. 2

图 2 三波段模型波段需求 Fig. 2 Spectrum requirements of the three band models

然而根据文献^{[32, 34, 35]}的研究，临近波段之间的吸收系数具有一定的比例关系. 可以用下式表示：

从GOCI数据的波段设置发现，B5中心波长位于660 nm附近，与B6(680 nm)比较接近，又由于非色素颗物与CDOM在长波波段的吸收系数相对固定，因此认为以下关系式成立：

与标准三波段模型相比，本研究提出的三波段模型利用a_ph波段间的固定斜率，更改了λ₂的位置(如图 2中λ^*₂所示)，使三波段模型能够应用于GOCI影像数据.

本研究中的GOCI数据利用6s模型(Second Simulation of the Satellite Signal in the Solar Spectrum)进行大气校正. 其中模型输入的卫星、 太阳几何角度，均获取自影像本身； 气溶胶光学厚度采用太湖站AERONET站点的实测值作为输入.

为了更加全面地评价本研究提出的模型，本研究首先利用GOCI和MERIS的波段响应函数，将地面实测高光谱数据模拟成为GOCI数据和MERIS数据，进而基于MERIS三波段因子[(1/B681-1/B708)×B753^[21]、 GOCI三波段因子和GOCI比值算法(B745/B680^{[26, 27]})分别估算叶绿素a浓度并进行对比分析. 估算模型通过常用的统计拟合方法获取. 3个模型均采用现行拟合方式获取.

为了使模型值域范围一致，先将数据按照叶绿素a浓度大小排序，再顺序抽取数据，最后得到186组建模数据，83组验证数据. 另外有10个样点用于实测GOCI影像的验证.

模型验证中，估算精度评价指标分别选取常用的平均相对误差(MAPE)和均方根误差(RMSE). 表达式如下：

2 结果与分析 2.1 遥感反射率特征

7期数据的遥感反射率曲线及各个波段与实测叶绿素a浓度的相关系数如图 3所示. 首先，7期数据都体现出较复杂的光学特性：400~500 nm由于叶绿素、 类胡萝卜素及可溶性有机物的强吸收作用，使得遥感反射率在该范围内相对较低； 550~580 nm范围内出现一反射峰，主要是由于叶绿素和胡萝卜素在该波段范围内吸收较弱以及细胞的散射作用形成的； 630 nm处由于藻蓝蛋白的吸收而形成一个小的反射谷； 675 nm附近由于叶绿素a的强吸收而形成一个明显的谷值； 700 nm附近由于叶绿素a的反射而形成一个明显的反射峰，该峰值是判断水体中是否含有藻类叶绿素的重要依据^[36]； 700 nm后由于水的吸收变强导致遥感反射率迅速减小，至730 nm左右出现一个拐点，此后遥感反射率光谱变化比较平缓，直到808 nm附近出现一个小的反射峰，该反射峰主要由于悬浮物的后向散射作用形成； 808 nm后遥感反射率继续减小，到900 nm之后反射率光谱出现较大噪声，信噪比急剧下降，一般不参与定量遥感反演. 同时，可以发现各期数据的相关系数曲线有较大差异. 由于叶绿素a浓度较低，2009年4月太湖数据集中，各个波段与叶绿素a浓度的相关性都较低(图 3)，相关系数集中在-0.3左右； 而在叶绿素a浓度更低的三峡水库，710 nm左右几乎没有荧光峰信息，叶绿素a浓度敏感波段移动到600 nm左右，与其余数据集有较大差异； 相关性最高的数据集是2013年8月洞庭湖数据集，在700 nm左右，遥感反射率与叶绿素a浓度相关系数达到0.8左右，并在长波波段持续保持在较高水平； 其余太湖8月数据集的相关性曲线体现出较为统一的特征，高相关系数都出现在近红外波段，在710 nm左右有一个小峰值. 与刘忠华^[18]的研究结论相似.

图 3

Fig. 3

图 3 7个数据集的遥感反射率数据及各个波段与叶绿素a浓度的相关系数 Fig. 3 Rrs spectrums of 7 datasets and their correlation coefficients with chlorophyll-a concentrations

总的来看，7个数据集体现出比较复杂的光学性质，仅从相关系数确立统计模型，缺乏有力的支撑.

所有采样点的色素颗粒物吸收系数曲线如图 4所示. 色素颗粒物吸收系数在440 nm出现一个峰值，该峰值由叶绿素a、 叶绿素b、 叶绿素c、 类胡萝卜素和其它辅助色素共同引起的. 480 nm附近，由于β胡萝卜素的吸收作用，吸收系数呈现肩状； 550 nm附近由于各种色素在绿光波段的弱吸收使得吸收系数在该处出现一极小值； 550 nm之后吸收系数逐渐增大，至625 nm附近由于藻蓝蛋白色素的吸收而出现一个小的峰值(仅部分吸收系数曲线值较高的样点体现出该特征)； 最后，在675nm附近存在明显的吸收峰，该峰值主要归因于叶绿素a的吸收作用^{[18, 37, 38]}； 在700 nm之后，色素颗粒物吸收作用明显减弱，吸收系数曲线急速下降，趋近于0.

图 4

Fig. 4

图 4 所有样点各个波段处的色素颗粒物吸收系数及与680 nm吸收系数的相关性 Fig. 4 The aph spectrums of all sampling stations and their correlation coefficients with aph(680)

将所有样点每个波段的色素颗粒物吸收系数，与680 nm处进行相关分析，结果如图 4所示. 从中可以看出，相关系数从400 nm开始，随波长增加，开始从0.6左右逐步增加； 至500 nm左右，收到色素弱吸收的影响，相关系数开始下降，直至550 nm达到谷值，接近0.6； 之后逐步升高，在660~690 nm之间达到最大，相关系数均在0.9以上； 随着波长进一步增加，色素吸收作用迅速减弱，相关系数曲线也开始急剧下降，到780 nm左右达到最低点. 由此可见，a_ph(660)与a_ph(680)之间相关性极高.

基于以上结果，绘制了所有样点a_ph(660)与a_ph(680)的散点图，如图 5所示. 从中可以看出，在截距为0的情况下，660 nm与680 nm处的色素颗粒物吸收系数的线性决定系数达到0.967. 因此可以认为，在本研究中，ε_{a ph}(680,660)=1.351. GOCI三波段模型假设合理.

图 5

Fig. 5

图 5 aph(660)与aph(680)的统计关系 Fig. 5 Statistical relationship between aph(660) and aph(680)

为了体现GOCI三波段算法的精度，将其与标准MERIS三波段算法和GOCI波段比值算法进行对比. 3个算法的建模数据如图 6所示. 从中可见，从整体来看，GOCI三波段模型散点与MERIS三波段模型较为接近，二者拟合效果均较理想. 其中GOCI三波段模型线性决定系数为0.809，略低于MERIS三波段模型的0.820，但相对于GOCI比值算法(R2=0.450)有显著提高.

对于GOCI比值模型而言，最明显的离群值出现在叶绿素a浓度低值区，如图 6(c)中虚框区域所示，在该区域中，叶绿素a浓度变化很小，集中在0~10μg ·L^-1之间，而比值模型自变量动态范围很大，涵盖了0.2~0.8左右的自变量范围，大于模型整体自变量范围的一半. 在模型应用时难免造成大量高估. 这主要是因为模型依赖两个波段的比值，且两个波段波长距离比较远的情况下，难以有效去除后向散射的影响，进而影响自变量对叶绿素a浓度的解释. 从另外一个角度来看，低值区多数点位属于2009年4月太湖数据集和2009年8月三峡水库数据集，图 3表明这两期数据的反射率信息与叶绿素a浓度正相关性很低，其中三峡水库各个波段都与叶绿素a浓度呈负相关关系，这也在一定程度上解释了构建模型时的异常情况：从统计意义上来看，三峡水库的数据与其余样区数据体现出相反的特征，因此难以用统一的模型对其进行表达.

图 6

Fig. 6

图 6 模型散点图 Fig. 6 Scatter plots of different models

从MERIS三波段模型的散点图[图 6(b)]来看，新波段的引入能够有效避免GOCI波段比值模型中的问题，消除后向散射的影响，使自变量对叶绿素a浓度的解释能力显著提高. 图 6(c)中在叶绿素a浓度低值区自变量出现的“扩散”现象得到有效抑制.

GOCI三波段模型效果与MERIS十分接近，表明模型中对λ₂波段位置的修改是有效且成功的. 但其效果略差于MERIS三波段模型. 因为相对于MERIS三波段模型，GOCI三波段模型的假设会带来更显著的误差：从图 4可以看出，a_ph在大于710 nm波长之后，几乎为零(红色均值曲线)； 而图 5表明，虽然所有样点吸收系数的线性相关性极高，仍有一些数据体现出斜率的微弱差异(如太湖2011年5月数据)，这种差异会最终体现在模型的构建中.

3个模型的叶绿素a浓度估算误差如表 2所示. 其中将MAPE分为叶绿素a浓度高值区和低值区两个部分，以便于分析误差分布. 总体上，模型误差与建模精度有较好的一致性：MERIS三波段模型误差最小，GOCI三波段模型误差略大，GOCI比值模型误差最大. 从MAPE(cChla<10)来看，3个模型效果都较差，这是受叶绿素a浓度测量误差与模型适用范围^[39]的影响，但仍然可以看出，GOCI三波段模型的误差与MERIS三波段模型接近，都明显优于GOCI比值算法. 从MAPE(Cchla≥10)来看，GOCI比值模型误差为0.426，GOCI三波段模型为0.325，提高了23％左右，相对于MERIS三波段模型的0.265，GOCI三波段模型高了18％左右. 从RMSE来看，GOCI比值模型误差为19.138 μg ·L^-1，GOCI三波段模型为16.308μg ·L^-1，降低了14％左右，相对于MERIS三波段模型的15.171μg ·L^-1，GOCI三波段模型高了6％左右.

表 2(Table 2)

表 2 模型误差 Table 2 Model errors 模型 自变量 表达式 MAPE (cChla<10) MAPE (cChla≥10) RMSE /μg·L-1 GOCI三波段 (1/B680-1/B660)×B745 y=786.659x-4.864 1.668 0.325 16.308 MERIS三波段 (1/B681-1/B708)×B753 y=261.629x-26.580 1.560 0.265 15.171 GOCI比值 B745/B680 y=125.274x-31.016 2.409 0.426 19.138

表 2 模型误差 Table 2 Model errors

综上所述，GOCI三波段模型可以得到与MERIS三波段模型接近的效果，无论从建模的确定性还是验证误差的大小，均明显优于目前广泛使用的GOCI波段比值算法.

将2013年5月13日的八景GOCI影像大气校正的结果，与10个地面实测点的光谱数据进行对比. 每个采样点都与时间上最邻近的图像进行对比,结果如表 3所示.

表 3(Table 3)

表 3 大气校正误差 Table 3 Atmospheric correction errors 波段 MAPE RMSE/sr-1 B412 0.423 8.030×10-3 B443 0.239 6.146×10-3 B490 0.127 4.436×10-3 B555 0.127 5.550×10-3 B660 0.140 4.393×10-3 B680 0.133 3.862×10-3 B745 0.238 2.528×10-3 B865 0.301 1.677×10-3 平均 0.216 4.578×10-3

表 3 大气校正误差 Table 3 Atmospheric correction errors

从表 3可以看出，大气校正效果整体较好，其中B412波段误差较高，MAPE达到0.423，RMSE为8.030×10^-3 sr^-1，其余波段误差均较低. 本研究中涉及到的B660、 B680和B745这3个波段的MAPE都在0.3以下，RMSE都在5×10^-3 sr^-1以下，精度较高.

将GOCI三波段模型和波段比值模型分别应用于2013年5月13日获取的8景GOCI影像，计算太湖叶绿素a浓度，并将与采样点时间最接近的影像估算得到的叶绿素a浓度，与实测叶绿素a浓度进行对比. 结果如图 7所示，可以看出，GOCI三波段算法较为均匀的分布在1 ∶1线附近，而GOCI波段比值算法有较明显的低估现象. 通过精度指标的计算，发现GOCI三波段算法的MAPE和RMSE分别为0.380和6.092 μg ·L^-1，而GOCI波段比值算法则分别为0.423和9.070 μg ·L^-1. 实测数据验证表明GOCI三波段算法具有更高的精度.

图 7

Fig. 7

图 7 GOCI三波段模型和比值模型的影像验证结果 Fig. 7 image validation result of the three band model and band ratio model

由于每天8景影像中，两种算法得到的结果差别类似，因此仅列举出11:00的影像进行说明，结果如图 8所示. 其中图 8(a)为太湖地区2013年5月13日的GOCI标准假彩色合成图像，图 8(b)为GOCI三波段估算结，图 8(c)为GOCI三波段估算结果，图 8(d)为三波段算法估算结果与波段比值算法估算结果的比值. 从假彩色合成图像中可以看到太湖北部，竺山湾、 梅梁湾、 贡湖湾北部水体以及部分沿岸水体影像偏红. 从GOCI三波段估算结果和GOCI波段比值估算结果中，可以看出叶绿素a浓度高值区都集中在竺山湾、 梅梁湾和贡湖湾水域中. 两种算法在高值区的结果较接近，但GOCI波段比值算法在低值区，即湖心区至南部沿岸区域间的大面积水体中，得到的结果较为集中，空间变化不明显. 而GOCI三波段算法得到的结果空间差异显著增大. 为了更直观地观察两种算法的差异，对生产得到的叶绿素a浓度影像进行比值计算，GOCI三波段叶绿素a浓度与波段比值叶绿素a浓度的比值结果如图 8(d)所示. 可见，在竺山湾、 梅梁湾以及贡湖湾，叶绿素a浓度较高的水体中，两种算法结果较为接近，在湖心区大面积水域，GOCI比值算法得到的结果明显低于GOCI三波段算法，在湖区西湖最为明显，二者差异达到5倍以上. 结合图 7的估算结果，认为GOCI波段比值算法在该期影像的应用中，对太湖大面积水域的叶绿素a浓度产生严重低估，GOCI三波段算法的结果更可信，适用性更强.

图 8

Fig. 8

图 8 GOCI影像叶绿素a浓度估算结果 Fig. 8 Chlorophyll-a concentration maps that yielded by GOCI image data

利用GOCI三波段模型，基于2013年5月13日太湖8景GOCI影像，得到2013年5月13日太湖叶绿素a浓度时空分布情况，如图 9所示. 从中可以看出：整体上，太湖叶绿素a浓度在2013年5月13日当天，呈现逐渐增加的趋势，湖心区叶绿素a浓度从08:00的5 μg ·L^-1，逐渐升高到15:00的20~40 μg ·L^-1. 其中08:00~11:00之间的升高趋势更明显. 与此同时，太湖北部梅梁湾、 贡湖湾沿岸湖区的叶绿素a浓度也有逐渐增加的趋势，高值区(叶绿素a浓度>100 μg ·L^-1)面积有所增长. 但南部沿岸湖区的高值区面积却明显减少. 这是由于受东南风影响，南部沿岸湖区聚集的藻类被驱散，而更多的藻类向梅梁湾、 竺山湾湖区聚集.

图 9

Fig. 9

图 9 2013年5月13日太湖叶绿素a浓度分布 Fig. 9 The cChla distribution at May 13,2013 of Taihu Lake

(1)680 nm色素颗粒物吸收系数与660~690 nm波段间的色素颗粒物吸收系数相关性极高，且680 nm色素颗粒物吸收系数与660 nm色素颗粒物吸收系数比值固定.

(2)适用于GOCI数据的三波段模型无论从建模数据集还是验证数据集来看，都与MERIS三波段模型效果接近，明显优于GOCI波段比值模型.

(3)新建立的模型克服了GOCI影像波段设置的劣势，为高时间分辨率、 高精度的内陆湖泊水体叶绿素a浓度监测提供技术支持.

(4) GOCI影像的独立验证结果表明，GOCI波段比值算法会对太湖叶绿素a浓度低值区造成显著的低估现象，GOCI三波段算法得到的结果更可信.