1 引言
冰冻圈是全球气候变化的指示器， 也是受气候系统影响最直接、 最敏感的圈层（秦大河等， 2014）。近10年来， 随着全球气候变暖， 冰冻圈的分布范围在明显缩减（钟歆玥等， 2022）。湖冰是地球冰冻圈的重要组成部分， 是北半球湖泊变化的季节性现象之一。湖冰的出现影响着湖泊的物理、 化学、 生物过程（Xie et al， 2022）， 是目前地球环境和气候变化研究的一个敏感因素（康世昌等， 2020）。湖冰物候， 即湖冰冻结/消融的周期性变化， 会导致湖面性质的突变， 进而影响湖水与大气之间的能量交换（Woolway et al， 2020； 邰雪楠等， 2022）。湖冰物候又受到区域气候的控制， 被作为监测气候变化的良好指标（Ke et al， 2013）。随着全球气候变暖， 湖冰正发生着不同程度的变化。例如， Wang et al（2021）追踪了1985 -2000年北温带33333个湖泊冰盖的长期变化， 发现这些湖泊1 -3月的冰覆盖频率在观测时间段内显著降低。Lopez et al（2023）的研究表明在1951 -2014年北半球152个湖泊中97%的湖冰消融日期表现出提前趋势。湖冰变化将会对人类的生产生活造成影响， 例如， 湖冰阻塞会对下游生态和人类安全构成威胁（汪关信等， 2020）， 湖冰物候变化也可以影响基础生态系统的物理条件（如蓄热、 温度）。近年来， 全球变暖和不断升级的极端气候导致了全球和区域性湖冰范围减少以及湖冰物候改变， 预计该现象在未来会更加恶化， 因此开展湖冰物候监测研究工作十分必要。
积雪作为冰冻圈的天然固体水库（Barnett et al， 2005）， 是全球水循环的主要组成部分（Tarca et al， 2022）。北半球冬季的积雪覆盖面积占40%~50%（Tang et al， 2022）， 对地表水文和能量平衡起着至关重要的作用（Chen et al， 2015）。在全球气候变暖的背景下， 北半球的积雪覆盖范围在1922 -2010年、 1978 -2015年两个时间段内显著减少（Brown and Robinson， 2011； Hori et al， 2017）。适量积雪可以起到一定储热作用， 使农作物免受冷空气的侵袭（赵琴等， 2022）， 而过量或少量的积雪则会影响植被生长（Xie et al， 2021）， 间接对农业经济造成不可估量的影响（杨芳芳， 2021）。积雪具有较高的反照率， 其动态变化对全球或区域气候具有一定的反馈作用（Sood et al， 2020； Tang et al， 2022）， 而气候变化也必然会影响积雪积累及积雪物候， 进而引起极端降水增加（王卫国等， 2022）。积雪物候是积雪动态变化的重要指标（杨雅茹等， 2022）， 可以直接或间接地反映地区气候、 水文和生态条件， 对全球气候和生态系统等方面具有重要影响（Zhang et al， 2022； 郝晓华等， 2022）。近年来北半球积雪物候也发生了显著变化， 总体表现为积雪初日推迟、 积雪终日提前、 积雪日数缩短（Brown and Robinson， 2011； Peng et al， 2013）。积雪物候变化可能会导致春季融雪径流的峰值提前、 流量增大， 夏季径流减少， 从而增加春季洪水和夏季干旱等灾难性事件的发生（Barnett et al， 2005）。
对于湖泊生态系统而言， 湖冰和冰上积雪的存在可以在一定程度上增大湖面反照率（Henneman and Stefan， 1999； 牛瑞佳等， 2023）， 降低湖泊对太阳辐射的吸收， 减缓湖冰消融速度， 并影响湖泊的初级生产力（李兴东等， 2022； 邰雪楠等， 2022）。积雪以其较低的导热性， 在湖冰冻结初期起到了很好的保温作用， 一定程度上阻碍了湖冰的进一步冻结（李兴东等， 2022）。然而， 当积雪厚度超出湖冰承受能力时， 湖冰的下沉会使得湖水浸入积雪， 进而导致积雪向湖冰的转化。因此， 揭示流域内湖冰及积雪物候的变化规律对湖泊生态系统研究具有重要意义。
地表物候是描述地物周期性变化现象的常用指标， 可以用于量化陆地生态系统对气候变化的响应（杨雅茹等， 2022； Wang et al， 2022）。青海湖流域作为国内重要的生态敏感区， 研究湖冰物候、 积雪物候以及二者之间的相关关系有利于揭示流域内气候和生态系统的变化规律， 对当地水文、 生态和植被生产力等方面的研究具有重要意义。研究表明， 青海湖流域冬季云量较高， 2000 -2020年冬季平均云量可达70%左右（张敬书等， 2023）。目前无云积雪和无云湖冰数据的制作主要聚焦于全球或全国等大范围区域， 但仍缺乏对于单一流域内的描述。本研究基于MODIS V6积雪产品和反射率产品（李延等， 2023）， 采用青海湖流域内草地下垫面的最佳积雪判别阈值和云去除方法， 分别获取了该区域内的积雪和湖冰物候， 并分析二者的空间分布特征。同时， 采用Theil-Sen Median法和一元线性回归法分析了近20年流域内物候的变化趋势、 低海拔区域内积雪与湖冰物候的相关关系， 以及冬季负积温和气温对气候的重要影响， 为青海湖流域生态系统保护工作提供了理论依据。
2 研究区域
青海湖流域位于中国东部季风区、 西北干旱区和青藏高原高寒区三大自然区的交汇地带， 海拔分布在3142~5285 m［图1 （a）］， 干旱少雨、 太阳辐射强、 日夜温差大。流域内年平均气温在-1~1 ℃， 年降水量为270~400 mm， 属半干旱气候。布哈河、 沙柳河和倒淌河等河流的地表径流和湖面降水是青海湖水源补给的重要途径（张洪源等， 2018； 赵仪欣等， 2023）。流域内地表覆盖类型根据1∶100万中国植被图划分为草地、 湖泊、 人造地表、 裸地、 冰/雪、 湿地、 森林、 耕地， 其中草地是主要的下垫面类型［图1 （b）］， 以草原和草甸为主（李广泳等， 2014； 王琪等， 2019）。青海湖位于青海湖流域东南部（99°36′E -100°47′E， 36°32′N -37°15′N）， 是中国最大的内陆湖。近年来青海湖流域气候趋于暖湿， 流域内径流增大和降雨增多， 使得青海湖水位不断升高、 面积逐渐扩大（李林等， 2020； 赵仪欣等， 2023）。据统计， 2000 -2019年， 青海湖面积总体呈先减后增趋势， 截至2019年， 青海湖水域面积为4597.35 km2， 2017年是自2005年以来青海湖水位上升幅度最大的一年（祁苗苗等， 2020）。
3 数据来源及预处理
3.1 MODIS积雪数据
3.1.1 数据来源
青海湖流域积雪物候数据的生成主要采用了MODIS V6逐日积雪产品MOD10A1/MYD10A1， 该产品空间分辨率为500 m， 包括原始NDSI（Normalized Difference Snow Index， NDSI）数据波段、 NDSI_Snow_Cover_Class波段、 Snow_Albedo_Daily_Tile波段、 NDSI_Snow_Cover_Basic_QA质量波段、 NDSI_Snow_Cover_Algorithm_Flags_QA波段。本文主要使用该数据产品提供的“NDSI”和“NDSI_Snow_Cover_Class”两个波段。由于草地是青海湖流域主要的下垫面类型， 根据高扬等（2019）关于青藏高原不同下垫面NDSI阈值研究结果得知草地NDSI最佳阈值为0.33， 因此本文选用0.33作为NDSI阈值进行青海湖流域二值积雪制图。
3.1.2 数据预处理
对2000 -2020年每年10月至次年5月MOD10A1/MYD10A1积雪产品的“NDSI”和“NDSI_Snow_Cover_Class”波段数据按以下规则分类： （1）若MOD10A1/MYD10A1相同位置的像元至少有一个对应的NDSI值大于0.33， 且是无云像元， 则该像元被记做积雪像元（像元值赋值为1）； （2）若MOD10A1/MYD10A1相同位置的像元都是云像元， 则该像元被记做云像元（像元值赋值为0）； （3）若MOD10A1/MYD10A1相同位置像元对应的NDSI值都小于0.33， 且是无云像元， 则该像元被记做无雪像元（像元值赋值为-1）；
将像元分类后， 采用Chen et al（2020）提出的时空自适应去云算法对存在的云像元和缺失值进行填充， 得到研究区域内近20年逐日无云积雪二值数据。为消除瞬时降雪对积雪物候结果的影响， 对所得到积雪二值数据做以下处理： 若连续两日的影像对应位置的像元都是积雪像元， 则该像元最终确定为积雪像元， 否则该像元确定为无雪像元。根据以下积雪物候定义， 基于所生成的积雪二值数据， 提取研究区积雪物候参数。积雪物候的定义如下： （1）积雪初日（Start of snow cover days， SCS）定义为一年中某像元位置第一次连续3天被识别为积雪像元的第一天； （2）积雪终日（End of snow cover days， SCE）定义为一年中某像元位置最后一次连续3天被识别为积雪像元的最后一天； （3）积雪日数（Snow cover days， SCD）定义为一年中某像元位置被分类为积雪像元的总次数。
3.2 MODIS地表反射率数据
3.2.1 数据来源
研究中采用2000 -2020年每年10月至次年5月逐日地表反射率产品MOD09GA/MYD09GA进行湖冰的识别， 该产品包括空间分辨率为500m的MODIS 1~7通道逐日地表反射率和其他质量波段数据。本文基于其红光波段（620~670 nm）和近红外波段（841~876 nm）， 采用多波段阈值法区分湖冰和湖水像元（庞毓雯等， 2020）， 同时利用质量波段进行云掩膜， 将影像数据分类为湖冰、 湖水和云像元后， 再对数据做去云处理。
3.2.2 数据预处理
根据水、 冰、 植被在红光波段和近红外波段的反射率差异， 综合考虑漏分和错分误差， 将满足红光波段大于0.20、 红光波段和近红外波段差值大于0.04的像元分类为湖冰， 不满足上述条件的像元则分类为湖水， 并使用质量波段对分类结果进行云掩模， 随之采用去云算法获取逐日无云湖冰和湖水数据集。最终基于湖冰物候定义逐年获取湖冰物候， 湖冰物候的定义如下： （1）湖冰开始冻结日期（Freeze-up start， FUS）定义为一年中某像元位置第一次连续3天被识别为湖冰像元的第一天； （2）湖冰开始消融日期（Break-up start， BUS）定义为一年中某像元位置最后一次连续3天被识别为湖冰像元的最后一天； （3）湖冰存在期（Exist Duration， ED）定义为一年中某像元位置被分类为湖冰像元的总次数。
3.3 气象数据
本文气象数据来源于美国国家海洋和大气管理局（https： //www.ncei.noaa.gov/）， 包括2000 - 2020年青海湖流域刚察气象站点逐日气温等数据， 并使用ERA5数据获取青海湖和青海湖流域的年均气温。
4 研究方法
为定量化描述流域内积雪和湖冰物候的变化趋势， 本文采用Theil-Sen Median趋势分析方法逐像元计算2000 -2020年物候变化率， 同时使用一元回归分析方法计算二者在时间尺度上的变化趋势， 一元线性回归结果也可对Sen趋势分析结果进行检验。
4.1 Theil -Sen Median趋势分析
Theil-Sen Median趋势分析方法又称作Sen斜率估计方法， 用于长时间序列数据分析， 能降低或避免数据异常或数据缺失对分析结果的影响， 是一种稳健的趋势分析方法。本文采用Sen斜率估计方法分析积雪和湖冰物候的空间变化趋势， 具体公式为：
式中： 为物候变化趋势； 、 为物候所在年份； 、 分别表示时间序列中第 项和第 项的物候值； 函数代表取中值。当 大于0时， 物候呈现延长/推迟的趋势， 反之， 物候呈提前/缩短的趋势。
4.2 一元线性回归模型
一元线性回归模型是分析自变量 和因变量 线性相关性关系的方法。对于计算2000 -2020年积雪和湖冰物候参数在时间尺度上的变化趋势， 需要对二者分别求年平均， 然后计算年平均影像整体均值， 生成时间序列。以时间（年）为自变量 ， 物候参数值为因变量 ， 用回归拟合法拟合出物候随时间序列变化的曲线， 得到回归方程： ， 可取正负值。当 时， 说明物候随时间呈推迟/延长趋势； 当 时， 说明物候随时间呈提前/缩短的趋势。该结果也可以作为Theil-Sen Median趋势分析结果的验证。
5 结果与分析
5.1 湖冰物候时空分布格局及变化趋势
在湖泊深度、 风向、 湖泊透明度等因素的影响下（邰雪楠等， 2022）， 湖冰物候存在着明显的空间异质性。研究中基于已获取的青海湖流域湖冰物候， 分析了物候的空间分布特征和变化趋势。
空间分布特征结果如图2 所示。图中湖冰开始冻结日期和开始消融日期的单位为一年中第几天（day of year， doy）， 因研究时间段为每年10月至次年5月， 因此该数值会大于365， 而湖冰存在期的单位为一年中多少日（day）， 两者都简称为d。在湖冰冻结过程中， 湖面东侧的湖冰开始冻结日期早于西侧［图2 （a）］， 东侧最早开始冻结日期为321 d， 随着时间推移， 湖面气温逐渐降低， 西侧和中部湖泊表面加速冻结。湖冰存在期表现出由东侧外围-湖心-西侧逐渐缩短的空间分布特征［图2 （b）］， 在湖泊西侧， 湖冰存在期最短， 处于0~87 d， 湖泊中部湖冰存在期为87~92 d， 湖泊东侧湖冰存在期最长， 处于92~174 d。湖冰开始消融日期表现出明显的区域特征［图2 （c）］， 由西往东， 湖冰开始消融日期逐渐推后， 集中分布在425~464 d。
湖冰物候受多种非气象因子（经度、 海拔）的综合影响（庄立超， 2021）， 将青海湖按经度划分为5个部分， 湖冰物候在不同的经度范围内具有明显的空间异质性（图3 ）。湖冰开始冻结日期随经度的增加而提前， 但由于低经度范围（99.6°E -99.8°E）内湖泊水深较浅， 湖冰的开始冻结日期也较早［图3 （a）］， 这也使得二者的 仅为0.08［图3 （b）］。图3 （c）和（d）中表明湖冰开始消融日期和湖冰存在期随经度的增加分别呈推迟和延长的现象， 二者的 分别为0.80与0.74。
图3 湖冰物候在不同经度带的分布情况 （a） 不同经度带湖冰物候的时间节点， （b） 不同经度带湖冰开始冻结日期， （c） 不同经度带湖冰开始消融日期， （d） 不同经度带湖冰存在期Fig.3 The distribution of lake ice phenology in different longitude zones.（a） the time node of lake ice phenology in different longitude zones， （b） freeze-up start in different longitude zones， （c） break-up start in different longitude zones， （d） exist duration in different longitude zones |
利用Sen斜率估计方法得到的空间变化趋势结果可知， 湖冰开始冻结日期主要呈推迟趋势， 但青海湖东北侧的湖冰开始冻结日期呈提前趋势［图4 （a）］。研究中借助一元线性回归方法发现近20年湖冰开始冻结日期整体上呈推迟趋势， 推迟率0.3 d·a-1［图5 （a）］， 与Sen趋势分析结果保持一致。湖冰开始消融日期整体上呈推迟趋势， 主要分布在布哈河入口以及周围各个小岛附近［图4 （b）］。湖冰开始消融日期在一元线性回归分析结果中也呈推迟趋势［图5 （b）］， 推迟率为0.2 d·a-1。此外， Sen趋势分析结果表明湖冰存在期主要呈缩短趋势［图4 （c）］， 仅在中部湖面以及布哈河入口水域处呈延长趋势。值得注意的是， 一元线性回归分析结果也表明湖冰存在期呈缩短趋势［图5 （c）］， 缩短率为0.6 d·a-1。基于上述分析结果， Theil-Sen Median趋势分析与一元线性回归分析结果具有很好的一致性， 也从侧面验证了两个分析方法结果的准确性。
图4 湖冰物候空间变化趋势（单位： d·a-1）（a） 湖冰开始冻结日期， （b） 湖冰开始消融日期， （c） 湖冰存在期Fig.4 Spatial variation trends of lake ice phenology.Unit： d·a-1.（a） freeze-up start， （b） break-up start， （c） exist duration |
图5 湖冰物候线性变化趋势 （a） 湖冰开始冻结日期， （b） 湖冰开始消融日期， （c） 湖冰存在期， 圆形符号短划线表示湖冰面积比例所得， 方形符号直线表示湖冰像元分类所得Fig.5 Linear variation trends of lake ice phenology.（a） freeze-up start， （b） break-up start， （c） exist Duration.Circular symbol short line denotes lake ice area ratio.Square symbol straight line denotes lake ice pixel classification |
以上分析结果基于逐像元分类， 而湖冰物候也可以根据湖冰面积比例获取。研究中分别以10%、 90%为阈值， 即湖冰面积比例开始大于10%的日期为湖冰开始冻结日期、 湖冰面积比例开始小于90%的日期为湖冰开始消融日期、 湖冰面积比例开始小于10%的日期为湖冰完全消融日期。此处湖冰存在期定义为完全消融日期与开始冻结日期之间的天数。图5 表明， 基于逐像元和湖冰面积比例所得到的湖冰物候变化趋势相一致， 只在数值上略有差异。
5.2 积雪物候时空分布格局及变化趋势
积雪物候的空间分布也存在明显的空间异质性（图6 ）。积雪初日、 积雪终日、 积雪日数分别分布在275~404 d、 353~484 d、 3~209 d范围内［图6 （a）~（c）］。图中积雪初日和积雪终日的单位为一年中第几天（day of Year， doy）， 因研究时间段为每年10月至次年5月， 因此该数值会大于365， 而积雪的单位为一年中多少日（day）， 两者都简称为d。流域内海拔大于4500 m的高海拔地区气温较低、 降雪较早， 积雪初日明显早于低海拔地区， 大致分布在每年的10月上旬。同时， 该地区内具有常年稳定积雪， 积雪日数普遍大于2个月。随着早春温暖气候的到来， 高海拔地区升温缓慢， 积雪终日可达每年的4月下旬。相比之下， 低海拔地区（海拔小于4500 m）积雪初日较晚、 积雪日数较短、 积雪终日较早， 其中海拔小于3300 m的地区内积雪日数仅为10日左右［图6 （d）］。为了更好地说明积雪物候与海拔的相关关系， 研究中利用一元线性回归法对二者进行了相关性分析［图6 （e）方框所示］。分析结果表明， 积雪初日、 积雪终日分别与海拔呈负相关和正相关关系， R 2分别为0.98与0.96。高海拔低地区内温度普遍偏低， 有利于降雪的发生， 而较低的温度却不利于其消融， 最终导致积雪日数显著小于低海拔地区。［图6 （e）］也表明积雪日数与海拔高度呈显著正相关， R 2高达0.85。
图6 积雪物候空间分布特征 （a） 积雪初日 （单位： d）， （b） 积雪终日 （单位： d）， （c） 积雪日数 （单位： d）， （d）不同高度带积雪物候的时间节点， （e）不同高度带积雪物候Fig.6 The spatial distribution characteristics of snow phenology.（a） Start of snow cover days （unit： d）， （b） End of snow cover days （unit： d）， （c） Snow cover days （unit： d）， （d）the time node of snow phenology at different heights， （e）the lake ice phenology at different heights |
积雪初日整体上呈提前趋势， 占流域总面积的69.21%， 且主要分布在低海拔地区， 而推迟趋势则主要分布在高海拔地区， 占流域总面积的30.69%［图7 （a）］。研究中借助一元线性回归方法发现， 近20年积雪初日整体上呈提前趋势， 提前率为0.8 d·a-1［图8 （a）］。大部分流域范围内积雪终日呈推迟趋势， 占流域总面积的88.26%， 而呈提前趋势的地区仅占11.74%［图7 （b）］， 一元线性回归结果也表明积雪终日在整体上呈推迟趋势， 推迟率为1.1 d·a-1 ［图8 （b）］。此外， Sen趋势分析结果表明， 积雪日数在流域59.78%的面积范围上均呈延长趋势， 且主要分布在高海拔地区， 而缩短趋势则主要位于在低海拔地区， 占流域总面积的40.22%［图7 （c）］。［图8 （c）］表明， 流域内积雪日数整体上呈延长趋势， 延长率为0.6 d·a-1， 与Sen趋势分析结果保持一致。
图7 积雪物候空间变化趋势（单位： d·a-1）（a） 积雪初日， （b）积雪终日， （c）积雪日数Fig.7 Spatial variation trends of snow phenology.Unit： d·a-1.（a） start of snow cover days， （b） end of snow cover days， （c） snow cover days |
5.3 物候的影响因素
气温是地气相互作用研究中的重要变量（单帅等， 2020）， 很大程度上影响着湖冰和积雪物候变化。表1 表明， 湖冰开始冻结日期与冬季负积温和气温均具有正相关关系， 在0.01的置信水平下， 与气温的相关性系数可达0.553。而湖冰开始消融日期和湖冰存在期与冬季负积温和气温则分别均具有显著的负相关关系， 二者与冬季负积温的相关性系数分别为-0.611和-0.562（P<0.01）， 与气温的相关性系数也分别可达-0.492和-0.508（P<0.05）。冬季负积温和气温的升高会减缓湖冰的冻结、 加速湖冰的消融， 使得湖冰冻结日期推迟、 湖冰消融日期提前， 进而导致湖冰存在期缩短。
表1 物候与气温的相关性分析Table 1 Correlation analysis between phenology and temperature |
| 因素 | 湖冰物候 | 积雪物候 | ||||
|---|---|---|---|---|---|---|
| FUS | BUS | ED | SCS | SCE | SCD | |
| 冬季负积温 | 0.207 | -0.611** | -0.562** | 0.207 | -0.323 | -0.213 |
| 气温 | 0.553** | -0.492* | -0.508* | 0.365 | -0.378 | -0.692** |
**和*分别表示通过了0.01和0.05的显著性水平检验（** and * indicate passing the significance level test of 0.01 and 0.05， respectively） |
在王迎春等（2023）的研究中， 降雨和气温的不同组合会显著引发积雪物候变化， 而气温是引起变化的主要原因。表1 中， 积雪日数与气温之间呈显著负相关， 在0.01的置信水平下， 二者的相关性系数可达-0.692。2018 -2020年， 研究区内降温明显， 积雪日数逐年延长， 这表明气温变化会直接影响积雪的积累和持续时间， 积雪日数对气温具有良好的响应关系。值得注意的是， 冬季负积温和气温对积雪初日和积雪终日的影响并不显著， 后续的研究中需要进一步分析。
5.4 湖冰物候与积雪物候之间的相关关系
积雪物候随海拔变化显著， 而青海湖的平均海拔为3196 m， 位于青海湖流域内低海拔地区， 因此本研究仅针对低海拔区域分析积雪物候与湖冰物候的相关关系。研究中基于青海湖矢量边界设置了10 km缓冲区， 统计后发现该区域内的海拔集中在3520 m范围内， 因此本文将流域内海拔低于3520 m的区域作为低海拔流域， 分析此区域内积雪物候与湖冰物候的相关关系。低海拔流域内， 积雪初日、 积雪终日分别主要呈提前、 推迟趋势［图9 （a）， （b）］， 与整个青海湖流域的积雪物候变化趋势保持一致。而积雪日数则主要呈缩短趋势， 与青海湖流域积雪日数的整体变化趋势相反［图9 （c）］。
部分积雪与湖冰物候参数之间存有潜在联系。［图10 （a）］表明， 积雪初日与湖冰开始冻结日期之间具有较为显著的负相关关系， 相关性系数为-0.404。而积雪作为湖面保温层， 积雪日数的增加也会很大程度上减缓湖冰消融速度， 致使湖冰消融日期推迟， 因此二者呈一定正相关关系， 相关性系数为0.349。值得注意的是， 其余积雪与湖冰物候参数间的相关性较低［图10 （b）， （c）］， 主要原因有两点： （1）研究中为了避免湖面积雪和湖冰的错分， 对湖泊进行了掩摸处理， 仅获取低海拔流域内陆地表面的积雪物候， 从而导致陆地积雪物候和湖冰物候间的对应性较差。（2）低海拔区域内降雪频率较低、 积雪日数较短， 积雪无法对湖面冰层的冻结和融化过程产生长期影响， 最终导致二者之间的相关性相对较低。在后续研究中， 可直接利用气象站点的降水数据， 并采用双温度阈值法区分降雨和降雪， 从而有针对性地分析湖冰与积雪之间的关联。
图10 湖冰物候及积雪物候时间变化趋势 （a） 湖冰开始冻结日期及积雪初日， （b） 湖冰开始消融日期及积雪终日， （c） 湖冰存在期及积雪日数。方形符号直线是湖冰物候； 圆形符号短划线是积雪物候Fig.10 Temporal variation trend of lake ice phenology and snow phenology.（a） freeze-up start and start of snow cover days， （b） break-up start and end of snow cover days， （c） exist duration and snow cover days.square symbol is lake ice phenology， circular symbol is snow cover |
除以上传统气象要素之外， 湖泊深度、 湖泊面积、 湖水盐度等因素都会对湖冰物候产生影响（姚晓军等， 2015）。湖泊深度和湖泊面积会在一定程度上影响湖泊的冻结过程（王智颖， 2017）， 而湖水盐度则会直接影响湖冰的冰点温度， 进而影响湖冰物候（曹晓卫， 2021）。此外， 在全球变暖的大背景下， 冰川融水使得青海湖主要补给河流的地表径流增加， 这些河流的季节性变化是否会对湖冰物候产生影响也有待进一步探究。
6 结论
本研究基于近20年的MODIS V6积雪产品和反射率产品， 分别采用指数阈值法和多波段阈值法， 获取了青海湖流域每年10月至次年5月逐日无云积雪和湖冰数据集， 进而得到了青海湖流域的积雪和湖冰物候， 分析了物候空间分布特征及变化趋势， 并对气温与物候之间的相关性、 低海拔陆地积雪物候与湖冰物候之间的关联进行了深层次的分析。得出以下结论：
（1） 湖冰开始冻结日期和开始消融日期分别在321~389 d和425~464 d， 二者均呈推迟趋势； 湖冰存在期处于0~174 d之间， 呈缩短趋势。青海湖东侧的湖冰开始冻结日期早于西侧， 湖冰开始消融日期由西往东逐渐推迟， 湖冰存在期则表现出由东侧外围-湖心-西侧逐渐缩短的空间分布特征。
（2） 积雪物候在整体上呈现出积雪初日提前、 积雪终日推迟、 积雪日数延长的变化趋势， 同时积雪物候在不同海拔带上也表现出不同的分布特征， 如高海拔区域内积雪初日最早、 积雪终日最晚、 积雪日数最长， 显著区别于低海拔地区。
（3） 冬季负积温和气温会在一定程度上影响湖冰物候。冬季负积温和气温升高会导致湖冰开始冻结日期推迟、 湖冰开始消融日期提前、 湖冰存在期缩短。对于积雪物候， 积雪日数与气温之间呈显著的负相关关系， 气温下降， 积雪日数增加。
（4） 湖冰物候与积雪物候的相关性分析结果并不显著， 在后续的研究中可以借助气象站点的降水数据， 并采用双温度阈值法来区分降雪和降雨， 进而深入分析湖冰物候与降雪的响应关系。