-
气候变化及其对生态系统造成的影响已经成为全球变化研究的重点问题[1]。随着气候变化强度的不断加剧,生态系统的结构和功能必定会受到影响[2]。森林生态系统约占陆地生态系统总面积的30%以上,是陆地生态系统的主体[3],有降低大气二氧化碳(CO2)浓度、涵蓄和调节降水、减少地表侵蚀、调节区域小气候等重要作用[3−7]。因此,在气候变化背景下,研究森林生态系统过程,如生态系统潜热通量(FLE),对深入理解环境因子对生态系统的调控机制具有重要意义[8−10]。
涡动相关法已成为从小时到年际时间尺度监测生态系统碳、水和能量通量的主要方法,并成为自下而上估算全球生态系统碳、水平衡过程的支柱[11−13]。此外,涡度相关数据也越来越多地被用于生态系统模型的校准和验证[14−16]。为了定量估计生态系统过程,准确插补由于不利的气象状况或仪器故障等原因造成的缺失数据是非常重要的[17−18]。不同插补方法得到的结果可能有很大差别[19−21]。目前,大多数研究只关注插补结果,鲜有研究对插补方法进行比较分析,因此,对现有主要插补方法进行比较分析可以更好地理解生态系统对气候变化的响应。
如今在世界范围内利用涡动相关系统已经建立了许多通量数据集,尽管各大通量网系统(如FLUXNET、加拿大通量网等)建立数据集的数据收集方法得到了统一,但组织间数据处理方法仍有所不同。例如FLUXNET通量网通常采用边缘分布抽样法对缺失数据进行插补[22],而加拿大通量网通常利用线性与非线性回归法对缺失数据进行插补[23]。这导致相关数据可比性不强[24−25]。考虑到气候变化强度不断增强,增强通量数据的可比性是定量预测未来全球生态系统水文循环和改善森林管理决策的关键。
生态系统蒸散发(ET)作为全球水循环的第二大组成部分,在生态系统内部与外界进行能量和物质交换过程中发挥着重要作用,同时还影响着许多其他重要的生态过程[26−27]。本研究于北京松山自然保护区生态环境及生物多样性监测站开展,以天然落叶阔叶林生态系统为研究对象,通过涡度相关法连续监测数据,分析比较落叶阔叶林生态系统FLE缺失值插补的不同方法。本研究假设:插补结果会低估FLE,同时,不同插补方法得到的FLE模拟值有差异。本研究主要解决以下2个问题:①插补FLE与实测FLE是否存在偏差?②插补FLE与环境因子间关系是否和实测FLE与环境因子间关系相同?
-
松山自然保护区生态环境及生物多样性监测站位于北京市延庆区松山国家级自然保护区内(40°30′48″N,115°47′11″E),离延庆区张山营镇中心约10 km,离北京中心约104 km,北侧与大海坨自然保护区相邻,东南侧为佛峪口水库。
研究区属于北温带,气候为大陆性季风气候。研究区受所处地形的影响,气温低,湿度高,是典型的山地气候。降水季节分布不均,保证率低,多年平均降水量为450.0 mm,集中在6—9月;多年平均气温为8.9 ℃,极端最高气温为38.1 ℃,极端最低气温为−28.6 ℃,多年平均无霜期为153 d。2018年12月,以通量塔为中心,在冬奥延庆赛区外围松山国家级自然保护区内设立50 m×50 m的样方(40°53′N,116°18′E,海拔为1161 m)。样地内主要乔木树种为胡桃楸Juglans mandshurica,其他乔木树种为大果榆Ulmus macrocarpa、白蜡Fraxinus chinensis等。灌木以绣线菊Spiraea salicifolia、铁线莲Clematis florida为主。经样地调查,胡桃楸平均胸径为10.96 cm,平均树高为4.71 m;大果榆平均胸径为7.20 cm,平均树高为3.88 m;白蜡平均胸径为7.64 cm,平均树高为3.97 m。
-
通量监测采用开路式涡度相关系统,仪器初始安装高度为20 m。测量系统主要由三维超声分速仪(WindMaster,Campbell Scientific, Inc. Logan)、开路CO2/H2O红外气体分析仪(7500 DS, Campbell Scientific, Inc. Logan)和数据采集器(CR1000,Campbell Scientific, Inc. Logan)构成。同步测量的气象变量包括:净辐射[四分量辐射仪(CNR-4, Campbell Scientific, Inc. Logan)]、空气温度与湿度[空气温湿度传感器(ATMOS14, Campbell Scientific, Inc, Logan, )]、降雨量[雨量桶(ECRN-100, Campbell Scientific, Inc. Logan)],传感器安装高度均为20 m。采用CR1000数据采集器(CR1000, Campbell Scientific, Inc. Logan)采集气象数据,采样频率与涡度相关监测频率同步。
通量塔5 m半径内随机布设的5套热通量板(HFP01, Campbell Scientific, Inc. Logan)监测土壤热通量,传感器布设深度为10 cm。土壤温度和土壤含水量(CSW)分别由通量塔5 m半径内随机布设的5套土壤温度传感器(5TM, Campbell Scientific, Inc. Logan)和土壤水分传感器(5TM, Campbell Scientific, Inc. Logan)测定,布设深度为10 cm。采用CR1000数据采集器(CR1000, Campbell Scientific, Inc. Logan)采集土壤监测数据。
-
使用EddyPro软件(version 4, LI-COR, Lincoln)对原始10 Hz数据进行预处理,处理步骤包括峰值去除、倾斜校正(双轴旋转)、传感器滞后校正、光谱校正、去趋势化和0.5 h通量计算[20]。本研究采用PAPALE等[28]的方法剔除0.5 h通量数据异常值,利用数据相对绝对中位偏差(MAD)的偏离程度识别异常数据。具体如下:
其中:FLEi、FLEi−1、FLEi+1分别为第i时刻及其前后时刻的FLE数值,di为FLEi在数列中的相对位置,di如果满足以下条件则被定义为异常值:
其中:z为阈值,不同的z可以用来评价该方法对数据的影响,本研究中z取4[28]。Md为数列di的中位数。MAD的计算公式如下:
通过自举法确定摩擦风速阈值[28],夜间摩擦风速低于阈值的湍流均认为交换不均匀,剔除对应FLE值。本研究中摩擦风速阈值为0.204 m·s−1。另外,为了避免由于红外气体分析仪光路受影响导致的测量误差,光路信号值低于70的数据也被剔除。2019年由于仪表故障、系统维护和质量控制等原因剔除掉的FLE占全年数据的37.96%。
-
边缘分布抽样法(MDS)通过RStudio 1.1.463中的REddyProc程序包完成[29]。该算法主要分为2种不同的情况:①只缺失主要数据,其他气象数据可用;②部分或所有气象数据缺失。在第1种情况下,搜索7 d时间窗口内相似气象条件下是否存在数值,如果有,用满足条件数值的平均值插补缺失值。其中相似气象条件为总辐射(RG, W·m−2)、空气温度(Ta, ℃)和饱和水汽压差(DVP, Pa)相较于缺失值对应的RG、Ta、DVP变化幅度不超过50 W·m−2、2.5 ℃和500 Pa。如果在7 d时间窗口内没有相似气象条件,则将时间窗口增加到14 d。以此类推。
在第1种情况下,14 d时间窗口内相似条件下存在数值,插补数据质量为高等;21 d时间窗口内相似条件下存在数值,插补数据质量为中等;28 d时间窗口内相似条件下存在数值,插补数据质量为低等。在第2种条件下,FLE缺失数据的时间间隔小于1 h,插补数据质量为高等;FLE缺失数据的时间间隔小于1 d,插补数据质量为中等;FLE缺失数据的时间间隔大于1 d,插补数据质量为低等。
-
线性回归法主要采用AMIRO等[30]的方法。缺失值的时间间隔小于2 h采用线性内插法插补,其中夜间数据(即RG<20 W·m−2)设为0。缺失值的时间间隔超过2 h分为生长季和非生长季2种情况。生长季(Ta >0 ℃)采用240个连续数据时间窗口内的FLE与净辐射减去土壤热通量(即可利用能量)线性回归插补,时间窗口每次移动48个时间点(1个时间点为0.5 h)。非生长季(Ta<0 ℃)内采用缺失值前5 d至后5 d对应时间的平均值插补。
-
人工神经网络(ANN)是一种可以映射输入与输出数据间非线性关系的模型。本研究使用BP神经网络,模型为4层结构。输入层包含用于进行模型预测的自变量,输出层包含需要预测的变量,输入层和输出层由2个隐藏层连接。第1个隐藏层的激励函数为tansig函数,第2个隐藏层的激励函数为线性函数,训练方法为Levenberg-Marquardt。1层中的每个神经节点都与相邻层的所有神经节点相连。第1层隐藏层的神经节点数量范围由下式确定:
其中:m为第1层隐含层神经节点数,n和l为输入层和输出层变量个数,α为1~10的正整数。第2层隐藏层神经元个数为1。
本研究中输入变量为RG、Ta、DVP和CSW。为了强调昼夜造成的影响而加入RG,同时将RG<20 W·m−2的值设置为0。白天输入变量为RG、Ta、DVP、CSW;夜晚RG=0,输入变量为Ta、DVP、CSW。本研究采用均方误差(ERMSE)来衡量不同m构建的神经网络模型的优劣(表1)。m=9的神经网络模型的ERMSE最低,故本研究中m取9。
隐藏层神经元个数 均方误差 平均值 1 2 3 4 5 3 0.009 7 0.010 0 0.009 9 0.009 6 0.010 0 0.009 8 4 0.009 3 0.009 8 0.009 6 0.009 8 0.010 0 0.009 7 5 0.009 4 0.009 9 0.010 3 0.009 8 0.009 9 0.009 9 6 0.009 4 0.009 4 0.008 8 0.009 1 0.009 3 0.009 2 7 0.008 9 0.009 3 0.009 9 0.009 5 0.009 6 0.009 5 8 0.009 4 0.009 2 0.009 3 0.008 6 0.008 7 0.009 0 9 0.008 6 0.008 8 0.008 5 0.009 2 0.009 1 0.008 8 10 0.009 4 0.009 0 0.009 5 0.009 1 0.009 1 0.009 2 11 0.009 1 0.009 1 0.009 1 0.008 9 0.009 3 0.009 1 12 0.009 2 0.010 0 0.009 3 0.009 3 0.009 1 0.009 4 说明: 均方误差下的1、2、3、4、5分别指人工神经网络模拟的重复数(即不同的重复)。 Table 1. Result of the fitting of the number of different neurons in ANNs
-
本研究中涡度相关系统采用开路式红外气体分析仪。为了避免出现冬季由于仪器运作散发热量导致的错误通量[31],数据分析时只使用2019年4月11日至10月27日的数据。
为了评估3种数据插补方法的优良程度,采用随机取点的方法去除30%的FLE有效实测数据(2019年选定日期内数据共11 040个,其中有效实测数据7 576个,即随机去除其中2 272个数据),分别采用上述3种方法插补。插补完成后与实测数据进行比较。
由于线性回归法将生长季夜晚FLE的缺失值直接填补为0,为了更好地表现线性回归法的结果,本研究将其分为2种情况:第1种为正常线性回归插补(线性回归法Ⅰ),第2种为去除夜晚FLE数据插补(线性回归法Ⅱ)。
通过分段平均法分析人为去除的有效实测FLE数据及3种数据插补结果与环境因子间的关系,其中线性回归法采用去除夜间FLE插补数据进行分析。Ta按2.5 ℃分段,即每隔2.5 ℃取1次FLE均值; DVP按0.2 kPa分段,分段平均后再拟合。Ta与FLE使用指数方程y=aexp(bx)进行拟合,其中y为FLE,x为Ta,a、b为拟合参数;DVP与FLE使用指数-二次方程y=exp(a+bx+cx2)进行拟合,其中y为FLE,x为DVP,a、b、c为拟合参数。
指数-二次方程可以很好地模拟FLE随DVP增大先升高后降低的趋势,可以根据方程计算出最适DVP (DVPopt)为−b/2c。
除MDS方法使用RStudio软件外,其余所有分析均在MATLAB (Version 7.5.0, The MathWorks)中进行。
Comparative study of interpolation methods for missing latent heat flux data
doi: 10.11833/j.issn.2095-0756.20230526
- Received Date: 2023-10-24
- Accepted Date: 2024-04-07
- Rev Recd Date: 2024-04-04
- Available Online: 2024-07-12
- Publish Date: 2024-07-12
-
Key words:
- eddy covariance /
- latent heat flux /
- data interpolation /
- deciduous broad-leaved forest
Abstract:
Citation: | YANG Qiang, LI Xinhao, DU Tao. Comparative study of interpolation methods for missing latent heat flux data[J]. Journal of Zhejiang A&F University, 2024, 41(4): 810-819. doi: 10.11833/j.issn.2095-0756.20230526 |