2002年—2018年太湖水体溶解二氧化碳浓度卫星遥感数据集

齐天赐，段洪涛，曹志刚，沈明，肖启涛，刘东，马金戈

1.中国科学院南京地理与湖泊研究所 中国科学院流域地理学重点实验室,南京210008;

2.中国科学院大学,北京100049;

3.西北大学 城市与环境学院,西安710027

1 引 言

近年来，大气中CO2浓度持续升高所带来的温室效应给社会和经济可持续发展带来严重威胁，已成为人类社会备受关注的热点问题。尽管全球湖泊总表面积仅占全球非冰川陆地面积的3.7%（Verpoorter 等，2014），但是湖泊水体往往因为其相对大气过饱和的溶解CO2浓度（cCO2）而成为显著的CO2排 放 源（Raymond 等，2013；Borges 等，2015），其排放潜力可能抵消约22%的陆地生态系统CO2净吸收（Ciais 等，2013），在全球碳收支中起到关键作用。然而，基于有限野外观测或间接计算的cCO2数据升尺度得到的全球CO2排放量的研究结果在数值上可以相差数倍（Cole 等，2007；Marotta 等，2009；Raymond 等，2013；Holgerson和Raymond，2016），准确估算全球湖泊CO2排放量仍是一项重大的挑战（Drake 等，2018），急需修订和改进。

受限于人力物力，目前已有的全球CO2排放估算研究中使用的湖泊cCO2数据多是短期的有限样点采样方法，这样的低频观测缺乏时空代表性，难以满足准确估算湖泊CO2排放的需求。由于湖泊内部物理和生物地球化学过程的影响（Balmer 和Downing，2011；Natchimuthu 等，2017；Xu 等，2019），cCO2表现出显著的时空异质性，需要较高的时空观测频率才能准确刻画cCO2变化的年均水平（Klaus等，2019；Loken 等，2019）。同时，目前在国内外只有极少数湖泊具有长时序、较高频率的连续观 测 记 录（Seekell 和Gudasz，2016；Xiao 等，2020；Yan 等，2021），而这些长期监测数据展示出cCO2高动态变化，突出了长期监测对CO2排放估算的重要性，也对新的高频率、大范围观测手段提出要求（杨平和仝川，2015；段魏岩和黄昌，2021）。

卫星传感器可以提供高频率、连续和大范围的观测数据，弥补野外观测方法的不足（张兵等，2021）。利用遥感可反演的环境变量（例如叶绿素a浓度、水温、黄色物质吸收等）来表征水体中控制cCO2的相关生物地球化学过程，从而实现水体cCO2的遥感估算，并重构其大范围、长时序的空间分布。卫星遥感技术在海洋及海岸带水体cCO2估算研究中已经取得了成功的应用（Bai 等，2015；Chen 等，2019），而在湖泊水体方面的研究也在稳步发展（Kutser 等，2015；Qi 等，2020）。对于瑞典湖泊，已有学者作了很多尝试工作，使用黄色物质吸收系数估算cCO2（Kutser等，2015）。对于富营养化湖泊，也有学者建立了以叶绿素a为主要变量的cCO2估算模型，实现了太湖长时序cCO2空间分布的重构（Qi等，2020）。

太湖作为中国第三大淡水湖，长期遭受富营养化问题困扰，频繁暴发蓝藻水华（Duan 等，2009），随之引发的各种环境问题也使得太湖成为生态环境和碳循环研究的热点（Lee 等，2014；秦伯强，2020），相关学者针对太湖进行了长期的野外观测。然而，虽然太湖的野外测数据时间序列长，但是其在空间和时间上的分布并不均匀（Xiao等，2020），因此基于现有数据得出的估计结果可能存在较大随机性和偏差。为了弥补实测数据的不足，本文以Qi等（2020）针对太湖构建的cCO2遥感估算模型为基础，先基于MODIS/Aqua 数据反演的水色参数等产品，通过经验模型逐像元估算太湖藻型湖区表层水体cCO2，再经统计平均得到2002-07—2018-12的月平均数据集。月平均cCO2数据集在时间尺度上对应常规生态环境监测频率，在此基础上利用卫星更高的观测频率来增强月均数据的时空代表性，减小采样随机性带来的偏差，为太湖生态环境研究提供支撑。

2 数据研发方法

基于覆盖太湖全湖（图1）的MODIS/Aqua 遥感影像数据，本文生产了2002年—2018年太湖月平均cCO2空间分布数据集。数据研发流程如图2 所示，总体思路为：首先，获取原始卫星数据及产品，对数据进行预处理得到遥感反射率数据并对数据进行标准化与质量控制；接着，利用已发表的遥感算法反演得到关键水环境参数，进而估算太湖逐像元cCO2；最后，对cCO2进行月尺度统计得到太湖月平均cCO2空间分布数据集。

2.1 MODIS遥感数据获取、预处理及筛选

本数据集使用美国国家航空航天局（NASA）地球观测系统星座中Aqua 卫星上搭载的中分辨率成像光谱仪MODIS （Moderate Resolution Imaging Spectrum-radiometer）传感器获取的太湖影像数据及其衍生产品作为模型输入变量估算太湖cCO2。从NASA 海洋水色网（https：//oceandata.sci.gsfc.nasa.gov/［2021-05-18］）下载太湖区域2002-07—2018-12 的MODIS/Aqua 的L1A 级数据供遥感反射率的计算使用，以及水面日平均光合有效辐射PAR（Einstein m-2·d-1）产品作为估算cCO2的输入变量。从NASA 地球观测系统数据和信息系统网站EOSDIS（https：//earthdata.nasa.gov［2021-05-18］）下载相同时间段的气溶胶光学厚度产品（MYD04_3K）用于大气校正，以及地表温度产品（MYD11A1）作为表层水温（LST）变量来估算cCO2（Liu等，2015）。

使用SeaDAS 7.5 软件对L1A 数据进行替代定标，并结合几何地理文件处理得到L1B 级数据。接着使用SeaDAS 7.5 中的l2gen工具，通过气象参数文件以及瑞利散射的查找表对L1B 数据进行处理，得到经过水汽、臭氧校正与瑞利散射去除后的反射率（Rrc）（Hu 等，2004）。另外，基于已有研究对主流大气校正算法在太湖的性能评估结果（Shen 等，2017），本文采用6SV 模型对太湖L1B数据进行精确大气校正以获取遥感反射率（Rrs）。具体步骤和参数选择如下：

（1）观测几何：从MODIS 辅助元数据中获取太阳天顶角、方位角、卫星天顶角和方位角。计算太湖区域平均观测几何输入6S模型。

（2）大气模型：选择中纬度大气模型，并根据影像日期选择冬季或夏季。

（3）气溶胶类型和光学厚度：选择大陆型气溶胶模型。气溶胶光学厚度使用MYD04_3K 产品在太湖区域取有效陆地像元平均值作为550 nm 处的光学厚度输入。该遥感产品在中国东部湖泊的适用性已得到验证评估（Shen等，2020）。

（4）地面和传感器高度：地面高度输入0 km；传感器高度设置-705 km。

（5）光谱响应函数：输入NASA EOSDIS 网站提供的MODIS传感器光谱响应函数。

（6）地面反射率特征：选择均一朗伯体和湖泊水体平均光谱。

将具有不同分辨率的MODIS 反射率与产品数据重采样采样到250 m的分辨率，以确保所有数据都具有相同的空间分辨率以支持进一步的计算。之后，为便于进行后续一致性检验和水质参数的反演，对数据进行以下处理：

（1） 将数据进行裁剪，获取太湖区范围（30.8°N—31.6°N，119.8°E—120.8°E）的数据；

（2）合成Rrc数据的真彩色图像（RGB），通过目视解译将存在云覆盖率高、湖体结冰、太阳耀斑的影像排除，保留云覆盖低的影像共1242 景，具体时间分布如表1所示；

表1 本数据集使用的MODIS/Aqua影像时间分布表Table 1 Temporal distribution of MODIS/Aqua images used in this dataset

（3）对于存在少量云像元的影像，将1240 nm波段的Rrc大于0.1的像元视作云像元剔除；

（4）使用浮藻指数（FAI）来剔除影像数据中的蓝藻水华像元（FAI＞-0.004）（Hu 等，2010），以避免水华浮渣对水色参数反演的影响；

（5）由于水生植被对遥感反射率（Luo 等，2020）和cCO2动态（Xiao 等，2017）都存在显著影响，因此将水生植被区连同陆地像元通过太湖矢量边界剔除。

2.2 关键水色参数反演与cCO2估算

叶绿素a 浓度（Chla）使用Shi 等（2017）提出的Chla 算法反演得到。该算法通过构建645 nm与859 nm 波段Rrc值的归一化差值指数与Chla 的经验关系实现对太湖水体Chla 的反演。该算法的决定系数R2=0.72，最低平均相对误差为24%，具有较高的反演精度。计算公式如下：

式中，Rrc（645）、Rrc（859）分别为MODIS 数据中645 nm及859 nm波段Rrc的值。

光合有效辐射漫衰减系数（Kd（PAR））通过490 nm 处的漫衰减系数（Kd（490））推算得到，两者的经验关系及Kd（490）的反演根据Huang等（2017）针对太湖开发的半分析算法反演得到。该算法通过645 nm 及531 nm 波段的遥感反射率Rrs（645）与Rrs（531）的比值计算645 nm的颗粒后向散射系数bbp（645）和吸收系数a（645）（式（2）—（5）），进而计算645 nm 的漫衰减系数Kd（645）（式（6）），最后由Kd（645）推算得到Kd（490）。由于该算法中各项经验系数是根据模拟与实测数据率定的，因此其遥感数据大气精校正的不确定性不会传递到参数中，故本文应用时不进行算法系数重新率定。计算公式如下：

式中，bbp（645）为645 nm 的颗粒后向散射系数，Rrs（645）与Rrs（531）分别为645 nm 及531 nm 波段的遥感反射率，a（645）为吸收系数，bb（645）为后向散射系数，U（645）为中间变量，rrs（645）为恰在水面以下的遥感反射率。

太湖水体逐像元的cCO2估算根据Qi 等（2020）针对太湖提出的经验算法进行计算。该模型以Chla、LST、Kd（PAR）及PAR 作为输入变量，构建了多元二次多项式经验模型实现对太湖水体cCO2的估算。计算方程如下：

式中，C，L，K，P为经过中心化的Chla，ln（LST），ln（Kd（PAR））和PAR 数据。以上处理计算均在MATLAB 2016a软件中进行。

2.3 月平均cCO2空间分布数据生成及产品展示

由于经数据质量控制后仍然可能存在少量异常像元，从而产生异常估算结果，故将日尺度的cCO2数据中数值小于0 及超过176 μmol·L-1（太湖野外采集数据中cCO2平均值+3×标准差（Qi等，2020））的数据视为异常值进行剔除。之后，将日尺度数据按月份分别进行平均，得到2002-07—2018-12的MODIS 月平均cCO2空间分布数据集，图3 展示的是由该数据集绘制的月平均cCO2空间分布。

基于数据集进行统计分析，发现太湖的月平均cCO2表现出较大的季节变化（图3）。在夏季和秋季（6—11 月），cCO2较低，而在冬季和春季（12—5月），cCO2较高。最低值发生在10月，最高值发生在1月。结合各湖区虚拟站点的cCO2时间序列变化，太湖中cCO2的空间分布在冬春季差异很大，而在夏季和秋季则有所减弱（图3、图4（a））。此外，1—5月，湖中cCO2的最大值在空间上有所变化，从西南地区向西北湖区移动。在年际尺度上（图4（b）），太湖年平均cCO2在2005年前较为平稳，在2005年达到峰值，此后逐年降低。2003年—2018年整个时段太湖年平均cCO2表现出显著下降趋势（0.80 μmol·L-1·a-1）。

图3 MODIS估算的太湖2002年—2018年月平均cCO2空间分布图Fig.3 Spatial distributions of monthly average MODIS-estimated cCO2 of Lake Taihu from 2002 to 2018

图4 MODIS估算的2002年—2018年太湖平均cCO2时间序列，虚拟站点位置如图1所示Fig.4 Time series of MODIS-estimated average cCO2 of Lake Taihu from 2002 to 2018，the locations of virtual sites were shown in Fig.1

3 cCO2遥感精度评估

3.1 验证数据获取及统计分析方法

本文中使用的野外观测验证数据来自太湖湖泊生态系统研究站（以下简称太湖站）2003年—2015年月度及季度湖泊水质调查的数据，通过水温、pH和总碱度计算得到表层水体cCO2。根据水文和营养状态将去除水生植被区的太湖分为西北湖区、梅梁湾、贡湖湾、湖心区及西南湖区5个子区域（Xiao 等，2017），分别位于5 个子区域的共21 个采样点分布如图1 所示。太湖站野外观测数据同MODIS 日尺度数据匹配的星地同步数据经过匹配准则（采样点采集时间与卫星过境时间差不超过±5 h，匹配像元周围3×3 像元窗口内数据值的变异系数小于10%）筛选，共得到的80 个匹配点对（Qi 等，2020）。匹配点对用于cCO2估算模型的精度评估与模型不确定性分析。精度评估所使用的均方根误差（RMSE）、平均偏差（MB）、平均比（MR）、无偏百分比偏差（UPD）和平均相对偏差（MRD）等统计指标的定义与计算公式，模型不确定性分析方法的细节，请参考Qi等（2020）。

3.2 模型性能评估

Qi 等（2020）提出的cCO2估算模型具有较高的精度，决定系数R2=0.84，均方根误差RMSE=11.81 μmol·L-1，无偏相对误差UPD=22.46%，模型训练cCO2变化范围为7.73—119.34 μmol·L-1。星地同步数据的对比结果显示（图5（a），表2），模型在整个藻型湖区的总体验证精度较好（RMSE=12.83 μmol·L-1，UPD=24.03%），其中西北湖区、梅梁湾、湖心区和西南湖区的验证精度优于总体精度（RMSE＜12.22 μmol·L-1，UPD＜24.98%），而贡湖湾的验证精度较总体精度要差（RMSE=16.48 μmol·L-1，UPD=32.97%）。此外，贡湖湾的MR接近1，说明模型估算结果在总体均值上无偏，通过平均得到的月均值能在一定程度上抵消像元尺度的误差。值得注意的是，西南湖区的匹配样点只有2 个，其统计指标可能缺乏统计意义。此外，太湖观测数据显示，MODIS/Aqua 传感器过境时刻的实测cCO2与24 h 连续观测得到的日均cCO2高度一致（R2=0.97，p＜0.01）（Qi 等，2020），说明卫星估算的cCO2具有全天均值的代表性，适用于长时间尺度cCO2动态变化的监测评估。

表2 野外观测与MODIS估算的太湖各湖区匹配样点—像元的cCO2差异统计结果Table 2 Statistic metrics of differences between cCO2 of pixel-sample matchups in different regions of Lake Taihu

不确定分析结果显示，cCO2估算模型的估算结果受输入变量的随机误差影响较小（RMSE=10.24 μmol·L-1和UPD=26.14%），仅占模型构建数据集范围（7.73—119.34 μmol·L-1）的约9%。同时，在随机误差的影响下模型的估算结果存在高估（MR＞1），相比真实值最大可能高估约30%，这一高估情况可能传递至月平均cCO2数据。更多关于模型性能评估的细节请参考Qi等（2020）。

3.3 数据集精度评估

由于野外观测的采样频率为月度或季度，其时空代表性不足以直接评估月平均cCO2遥感数据集的精度，故将野外观测数据与MODIS月均cCO2数据集统计至年尺度，即计算各个子区域中所有采样点或像元的年平均cCO2值，对两者进行比较以评估数据集的估算精度。结果显示两者保持着较好的一致性（图5（b））。其中，西北湖区和西南湖区的两个数据集年均值差异要大于其他3个湖区，西北湖区MODIS 估算结果偏低，西南湖区MODIS 估算结果偏高，而全湖年平均值的差异要小于各子区域年均值的差异（表3）。这主要是由于野外采样点空间分布不均匀及部分季度样点时间频率低导致的代表性不足对年平均估算结果的影响。西北湖区的月度采样点主要靠近河口区域，而开阔水域的样点区域均为季度采样点；西南湖区的采样点则全为季度采样点；而其他3个湖区中的采样点基本是月度采样点，进而导致了年均值估算的较大差异。

表3 野外观测与MODIS估算的太湖不同区域年平均cCO2差异统计结果Table 3 Statistic metrics of differences between annual average cCO2 in different regions of Lake Taihu

图5 野外观测与MODIS估算的太湖各区域cCO2验证散点图，各子湖区的空间位置如图1所示Fig.5 Validating catter plot of cCO2 estimated by field observation and MODIS in different regions of Lake Taihu，Locations of different regions were shown in Fig.1

综上所述，虽然MODIS 估算结果与野外观测数据仍有一定偏差，但是主要是由于数据本身时空代表性差异导致，本数据集对于深入理解太湖cCO2乃至碳循环过程的时空变化规律仍具有极高的应用潜力，值得推广使用。

4 数据集存储格式与命名

数据集中月平均cCO2以Geotiff 格式储存，数值单位为μmol·L-1，被掩膜的像元及被剔除的异常像元数值为NaN。数据文件的地理坐标系为GCS_WGS1984，数据命名方式为“Aqua_TH_yyyymm_CO2monthly.tif”，其中yyyy 为年份，mm 为月份。整个数据集共198 个文件（下载方式：https：//doi.org/10.5281/zenodo.4729048［2021-05-18］），从2002-07—2018-12，时间段内数据无缺失月份。

5 数据集优势与应用价值

本文提供的太湖月平均cCO2数据集相比于前人发表文献中的数据，在时空代表性上均有提升与改进，更加适合相关但非对应专业学者分析使用。首先，在统计月平均数据时cCO2表现出剧烈的时空变异（图6），即使是在cCO2较为均一的湖心区和西南湖区，以及夏秋季节（6—9月），其平均cCO2变异系数都大于15%，证实了传统野外观测方法在监测太湖cCO2动态变化方面的局限性与不足。本数据集生产过程中每个月可利用影像数量为6.27±3.04景（平均值±标准差，表1），基本满足cCO2动态监测的频率需要（Klaus等，2019）。

图6 月平均cCO2变异系数箱型图Fig.6 Box plots of coefficient of variation for monthly average cCO2 in Lake Taihu

其次，本数据集经过质量控制与统计处理，除了少数月份因有效影像数据不足导致的数据缺失外，每个月份都能展示出完整的太湖藻型湖区cCO2空间分布（图3）。这相比于以往基于野外采样数据插值得到的空间分布结果（Xiao 等，2020）具有更加合理且丰富的空间细节，有利于深度挖掘cCO2的时空动态变化。相比于日尺度的空间分布结果（Qi 等，2020），虽然时间分辨率有所下降，但是月平均数据极大程度的削减了由于云层、蓝藻水华覆盖所导致的数据缺失，降低了非遥感专业学者使用数据的门槛，而且月尺度的时间尺度可以更好的匹配传统生态环境野外调查的时间频率，支持其研究分析。

本数据集展示的太湖月尺度cCO2空间分布，对太湖CO2通量估算以及区域碳收支研究有着重要意义。相比于以往以发表或共享的基于野外调查采样点数据估算的多年平均CO2通量空间分布（Xiao等，2020）或是基于涡度相关方法测量的梅梁湾观测站点的2012年—2014年CO2通量数据（Xiao 等，2014），本数据集提供了太湖CO2通量更加细致的空间分异信息。以2018年的月平均cCO2数据为例，其与大气CO2浓度的饱和度比显示出剧烈的空间变化（图7）。例如，秋季（9—11月）饱和度比在西北湖区沿岸表现出类似于河口羽状流的CO2排放热点（饱和度比＞1），在开阔水域表现出大小不一的CO2吸收区域（饱和度比＜1）。这些空间信息不仅有助于减小太湖CO2通量估算的不确定性，还有助于深入理解湖泊碳循环的时空过程。

图7 MODIS估算的太湖2018年月平均cCO2相对大气的饱和度比空间分布图Fig.7 Spatial distributions of monthly average MODIS-estimated cCO2 saturation ratio in equilibrium with the atmosphere of Lake Taihu in 2018

6 结 论

本文使用基于MODIS 数据反演得到的水色参数以及NASA 发布的卫星产品，通过经验模型对太湖藻型湖区表层水体cCO2进行逐像元的估算，并对结果进行统计平均得到2002-07—2018-12 的长时序月平均cCO2数据集，以期与历史湖泊调查数据的时间频率一致，并通过统计平均以增强数据集的时空代表性，进而提升其对深入理解湖泊中cCO2时空变化规律与湖泊碳循环过程研究的应用潜力。本文给出了估算模型本身的性能评估以及不确定性分析结果，并进一步评估了遥感数据集与野外数据年均值的一致性，研究人员可以以此为依据判评估本产品对自身研究的适用性，并进一步分析其研究结果的不确定性大小及来源。此外，本文还展示了数据集的优势特点与应用价值，希望能给相关学者以启发。