技术领域
[0001] 本
发明属于气象
水文测量技术领域，具体来说，涉及到一种基于多算法集成的全球陆表蒸散估算系统及方法。
背景技术
[0002] 当前，全球水资源日益匾乏，为了合理利用和分配水资源，迫切需要深入了解不同植被
覆盖和土地利用条件下的
蒸发耗水情况。
蒸散量包括植被截流蒸发量、植被蒸腾量、
土壤蒸发量和水面蒸发量，是区域水量平衡和
能量平衡的主要成分，不仅在
水循环和能量循环过程中具有极其重要的作用，而且是生态过程与水文过程的重要纽带。陆表蒸散是水圈、大气圈和
生物圈水分和能量交换的主要过程参量，是陆地表层能量循环、水循环和
碳循环中的关键变量。准确及时地获取特定区域的地表蒸散量数据，在农业、水文、森林和生态等领域都能够发挥十分重要的作用。因此，建立一种测定系统及方法，实现准确、快速的获取大区域陆表蒸散估算具有重要意义。
[0003] 然而，传统的测定方法对大范围内的地表蒸散状况代表性不强，且布设测点成本高，难以形成实用的观测网络。中国
专利(CN200810179250.0)就公开了一种基于遥感的区域蒸散量监测方法。该方法利用多时相的可见光/热红外遥感数据及每日常规气象数据建立了区域尺度上陆面蒸散量监测的
数据处理链，从而实现了区域蒸散量的运行性遥感监测。相对于传统通量观测站上测量蒸散的监测手段，遥感技术能够提供地表的多源多维多时相信息，具有大面积、宏观、实时和动态等优势，为陆表蒸散估算开辟了新的途径。遥感估算陆表蒸散的方法很多。主要有基于地表能量平衡的物理模型、经验统计算法、Penman-Monteith算法、遥感三
角形方法、数据同化方法等。这些模型方法是通过莫宁-奥布霍夫相似理论建立地表气象要素、遥感反演的关键参数(地表
温度、植被指数等)与陆表蒸散的定量联系来获取陆表蒸散，但是不同算法之间的估算的蒸散量差异很大，甚至差距达到2倍以上，总体
精度不高。而且，单一的蒸散算法存在很大的不确定性，少数采用简单的多算法平均值方法，简单平均值方法的不确定性较大，严重影响了陆表蒸散的计算精度。
发明内容
[0004] 为了解决上述技术问题，本发明提供了一种基于多算法集成的全球陆表蒸散估算系统及方法。
[0005] 本发明所述的一种基于多算法集成的全球陆表蒸散估算系统包括：MODIS蒸散算法模
块，订正的PM遥感蒸散算法模块，Priestly-Taylor蒸散算法模块，改进的Priestly-Taylor蒸散算法模块，基于经验方法的蒸散算法模块，观测数据提取模块，贝叶斯模型集成模块；首先利用MODIS遥感数据和MERRA再分析气象数据驱动MODIS蒸散算法模块、订正的PM遥感蒸散算法模块、Priestly-Taylor蒸散算法模块、改进的Priestly-Taylor蒸散算法模块、基于经验方法的蒸散算法模块，分别得出各自的蒸散量；然后根据
站点实测数据样本
质量由观测数据提取模块从全球通量站点观测数据提取地面蒸散观测值；最后通过贝叶斯模型方法，将提取的地面蒸散观测值和各个蒸散算法模块算出的蒸散量进行集成，形成贝叶斯模型集成模块。
[0006] 本发明所述的一种基于多算法集成的全球陆表蒸散估算方法，所述方法的具体步骤包括：
[0007] 1)构建MODIS蒸散算法模块：把冠层分为湿冠层和干冠层，把土壤分为水分饱和土壤和水分未饱和土壤，根据MODIS蒸散算法，计算与设计地表阻抗、空
气动力学阻抗以及
边界层参数；
[0008] 2)构建订正的PM遥感蒸散算法模块：根据传统的Penman-Monteith蒸散算法，对温度限制因子线性变化的
缺陷进行订正；
[0009] 3)构建Priestly-Taylor蒸散算法模块：根据现有的Priestly-Taylor模型，利用
生态系统胁迫因子扩展Priestly-Taylor系数，得到Priestly-Taylor蒸散算法；
[0010] 4)构建改进的Priestly-Taylor蒸散算法模块：采用表观热惯量参数化土壤水分蒸发因子修正土壤水分限制因子，得到改进的Priestly-Taylor蒸散算法；
[0011] 5)构建基于经验方法的蒸散算法模块：根据经验蒸散模型算法，设计蒸散输入、输出和计算模块；
[0012] 6)根据站点实测数据样本质量由观测数据提取模块从全球通量站点观测数据提取地面蒸散观测值；
[0013] 7)利用MODIS遥感数据和MERRA再分析气象数据驱动步骤1)的MODIS蒸散算法模块、步骤2)订正的PM遥感蒸散算法模块、步骤3)Priestly-Taylor蒸散算法模块、步骤4)改进的Priestly-Taylor蒸散算法模块、步骤5)基于经验方法的蒸散算法模块，分别得出各自的蒸散量；
[0014] 8)通过贝叶斯模型方法，将步骤6)提取的地面蒸散观测值和步骤7)的各个蒸散算法模块算出的蒸散量进行集成，形成贝叶斯模型集成模块。
[0015] 所述步骤1)中MODIS蒸散算法为：
[0016] LE＝LEwet_c+LEtrans+LEsoil
[0017]
[0018]
[0019]
[0020]
[0021]
[0022] 式中，LE为陆表蒸散，LEwet_c为湿冠层蒸发过程产生的蒸发，LEtrans为植被冠层蒸腾过程产生的蒸腾，LEwet_soil为饱和土壤水分蒸发产生的蒸发，LEsoil_pot为未饱和土壤潜在蒸发产生的蒸发，LEsoil为土壤蒸发，Δ为饱和水汽压与温度曲线的斜率，γ为干湿球常数，VPD为空气的饱和水气压与实际水汽压差，RH为
相对湿度，ra为
空气动力学阻抗，rs为地表阻抗，Fc为植被覆盖度，Fwet为湿度覆盖面积，ρ为空气
密度，Cp为空气的
比热容，Rnc为植被冠层所占的净
辐射能量，rhrc为冠层对应部分
空气动力学阻抗，rvc为地表热传输阻抗，Pa为
大气压，Rns为地表土壤所占的净辐射能量，rtot为整个空气动力学阻抗，ras为部分空气动力学阻抗，β，λ和ε为系数；
[0023] 所述步骤2)中温度限制因子订正方法为：
[0024]
[0025] 式中，mT为温度限制因子，Ta为空气温度，Topt为植被生长最适宜温度，设为25℃；
[0026] 所述步骤3)中Priestly-Taylor蒸散算法为：
[0027]
[0028]
[0029]
[0030] 式中，fwet为相对地表湿度，RH为相对湿度，fwet＝RH4，fg为绿色冠层覆盖度，fT为植被温度限制因子，,fM为植被水分限制因子，fSM为土壤水分限制因子；
[0031] 所述步骤4)中土壤水分限制因子为：
[0032]
[0033] 式中，DT为空气昼夜温差，DTmax为最大昼夜温差，取40℃；
[0034] 所述步骤5)中基于经验方法的蒸散算法为：
[0035] LE＝Rn(a0+a1NDVI+a2Ta+a3DTaR)
[0036] 式中，DTaR为日昼夜温差，NDVI为归一化植被指数，a0，a1，a2和a3为经验系数；
[0037] 所述步骤8)中贝叶斯模型集成模块计算公式为：
[0038]
[0039] 式中，r表示需要估算的蒸散，在某一时刻的该蒸散的观测值为rt，模型集合K，即{f1,f2,…,fk}，为估算r的所有单个模型的集合，p(r|fk)是单个模型fk估算蒸散r的概率密度函数，p(fk|rt)是单个模型fk的后验概率，能够反应单个模型fk与观测数据符合的程度。
[0040] 所述步骤8)中，贝叶斯模型集成模块得到的估算值的概率密度函数p(r|f1,f2,…,fk)就是单个模型概率密度函数的加权平均，所有单个模型的后验概率之和为1，即这样，可以把每个模型的后验概率看作是权重wk；
[0041] 假设p(r|fk)为高斯分布， 为该高斯分布的平均值，σ2K为方差，用参数和g(·)表示相关的高斯分布为：
[0042] p(r|fk)＝g(r|θk)
[0043]
[0044] 得到的估计值为r的条件期望，用期望-最大算法对代价函数进行优化得到权重wk。
[0045] 所述步骤1)通过植被覆盖度Fc来区分植被和土壤分配的全球陆表可利用能量，利用M ODIS遥感数据中植被吸收的光合有效辐射比例FPAR产品确定Fc；
[0046] 所述步骤2)订正的PM遥感蒸散算法中，统一的水分限制因子VPD的最大值定为2.79KPa；
[0047] 所述步骤3)Priestly-Taylor蒸散算法中的Rns可以用地表净辐射Rn和叶面积指数LAI计算，具体为Rns＝Rnexp(－kLAI)，k为经验系数；
[0048] 所述步骤4)改进的Priestly-Taylor蒸散算法的土壤水分限制因子中，DT是通过日昼夜温差计算得到或通过地面
温度计算得到；
[0049] 所述步骤5)中经验系数a0，a1，a2，和a3是利用全球通量站点观测数据作为
训练数据回归得到的；
[0050] 所述步骤6)的观测数据提取模块中，站点实测数据样本质量是通过质量标识判断而得到的；
[0051] 所述步骤7)从MODIS遥感数据中选取叶面积指数、植被指数，从MERRA再分析气象数据中选取基本气象参数，包括温度、相对湿度，驱动各个算法模块计算得出各自的蒸散量。
[0052] 所述全球陆表包括森林、灌木、农田、草地和湿地。
[0053] 与现有方法相比，本发明所述的一种基于多源数据和多种算法集成的全球陆表蒸散估算系统及方法的有益效果是：通过理论分析、实验测试和实际应用，充分利用地面实测资料、遥感数据和多算法相结合，以为贝叶斯模型依托，既有明确的物理
基础，又具有
稳定性好、简单、灵活、易于操作等优点，可对其中的MODIS蒸散算法、订正的PM遥感蒸散算法、P riestly-Taylor蒸散算法、改进的Priestly-Taylor蒸散算法以及经验方法的蒸散算法综合集成，比单一的蒸散算法具有更强的稳定性，适用面更加宽广，有广阔的应用前景。
附图说明
[0054] 图1为本发明采用的基于多算法集成的全球陆表蒸散估算系统的结构示意图；
[0055] 图2为本发明基于多算法集成的全球陆表蒸散估算方法的
流程图；
[0056] 图3为全球通量站点实测值验证。
具体实施方式
[0057] 下面结合具体的
实施例对本发明所述的一种基于多算法集成的全球陆表蒸散估算系统及方法做进一步说明，但是本发明的保护范围并不限于此。
[0058] 实施例1
[0059] 本发明设计的基于多算法集成的全球陆表蒸散估算系统包括：MODIS蒸散算法模块，订正的PM遥感蒸散算法模块，Priestly-Taylor蒸散算法模块，改进的Priestly-Taylor蒸散算法模块，基于经验方法的蒸散算法模块，观测数据提取模块，贝叶斯模型集成模块；首先利用MODIS遥感数据和MERRA再分析气象数据驱动MODIS蒸散算法模块、订正的PM遥感蒸散算法模块、Priestly-Taylor蒸散算法模块、改进的Priestly-Taylor蒸散算法模块、基于经验方法的蒸散算法模块，分别得出各自的蒸散量；然后根据站点实测数据样本质量由观测数据提取模块从全球通量站点观测数据提取地面蒸散观测值；最后通过贝叶斯模型方法，将提取的地面蒸散观测值和各个蒸散算法模块算出的蒸散量进行集成，形成贝叶斯模型集成模块。
[0060] 与之相应的基于多算法集成的全球陆表蒸散估算方法的具体步骤包括：
[0061] 1)构建MODIS蒸散算法模块：把冠层分为湿冠层和干冠层，把土壤分为水分饱和土壤和水分未饱和土壤，根据MODIS蒸散算法，计算与设计地表阻抗、空气动力学阻抗以及边界层参数；
[0062] 2)构建订正的PM遥感蒸散算法模块：根据传统的Penman-Monteith蒸散算法，对温度限制因子线性变化的缺陷进行订正；
[0063] 3)构建Priestly-Taylor蒸散算法模块：根据现有的Priestly-Taylor模型，利用生态系统胁迫因子扩展Priestly-Taylor系数，得到Priestly-Taylor蒸散算法；
[0064] 4)构建改进的Priestly-Taylor蒸散算法模块：采用表观热惯量参数化土壤水分蒸发因子修正土壤水分限制因子，得到改进的Priestly-Taylor蒸散算法；
[0065] 5)构建基于经验方法的蒸散算法模块：根据经验蒸散模型算法，设计蒸散输入、输出和计算模块；
[0066] 6)根据站点实测数据样本质量由观测数据提取模块从全球通量站点观测数据提取地面蒸散观测值；
[0067] 7)利用MODIS遥感数据和MERRA再分析气象数据驱动步骤1)的MODIS蒸散算法模块、步骤2)订正的PM遥感蒸散算法模块、步骤3)Priestly-Taylor蒸散算法模块、步骤4)改进的Priestly-Taylor蒸散算法模块、步骤5)基于经验方法的蒸散算法模块，分别得出各自的蒸散量；
[0068] 8)通过贝叶斯模型方法，将步骤6)提取的地面蒸散观测值和步骤7)的各个蒸散算法模块算出的蒸散量进行集成，形成贝叶斯模型集成模块。
[0069] 所述步骤1)中MODIS蒸散算法为：
[0070] LE＝LEwet_c+LEtrans+LEsoil
[0071]
[0072]
[0073]
[0074]
[0075]
[0076] 式中，LE为陆表蒸散，LEwet_c为湿冠层蒸发过程产生的蒸发，LEtrans为植被冠层蒸腾过程产生的蒸腾，LEwet_soil为饱和土壤水分蒸发产生的蒸发，LEsoil_pot为未饱和土壤潜在蒸发产生的蒸发，LEsoil为土壤蒸发，Δ为饱和水汽压与温度曲线的斜率，γ为干湿球常数，VPD为空气的饱和水气压与实际水汽压差，RH为相对湿度，ra为空气动力学阻抗，rs为地表阻抗，Fc为植被覆盖度，Fwet为湿度覆盖面积，ρ为空气密度，Cp为空气的比
热容，Rnc为植被冠层所占的净辐射能量，rhrc为冠层对应部分空气动力学阻抗，rvc为地表热传输阻抗，Pa为大气压，Rns为地表土壤所占的净辐射能量，rtot为整个空气动力学阻抗，ras为部分空气动力学阻抗，β，λ和ε为系数；
[0077] 所述步骤2)中温度限制因子订正方法为：
[0078]
[0079] 式中，mT为温度限制因子，Ta为空气温度，Topt为植被生长最适宜温度，设为25℃；
[0080] 所述步骤3)中Priestly-Taylor蒸散算法为：
[0081]
[0082]
[0083]
[0084] 式中，fwet为相对地表湿度，RH为相对湿度，fwet＝RH4，fg为绿色冠层覆盖度，fT为植被温度限制因子，,fM为植被水分限制因子，fSM为土壤水分限制因子；
[0085] 所述步骤4)中土壤水分限制因子为：
[0086]
[0087] 式中，DT为空气昼夜温差，DTmax为最大昼夜温差，取40℃；
[0088] 所述步骤5)中基于经验方法的蒸散算法为：
[0089] LE＝Rn(a0+a1NDVI+a2Ta+a3DTaR)
[0090] 式中，DTaR为日昼夜温差，NDVI为归一化植被指数，a0，a1，a2和a3为经验系数；
[0091] 所述步骤8)中贝叶斯模型集成模块计算公式为：
[0092]
[0093] 式中，r表示需要估算的蒸散，在某一时刻的该蒸散的观测值为rt，模型集合K，即{f1,f2,…,fk}，为估算r的所有单个模型的集合，p(r|fk)是单个模型fk估算蒸散r的概率密度函数，p(fk|rt)是单个模型fk的后验概率，能够反应单个模型fk与观测数据符合的程度。
[0094] 所述步骤8)中，贝叶斯模型集成模块得到的估算值的概率密度函数p(r|f1,f2,…,fk)就是单个模型概率密度函数的加权平均，所有单个模型的后验概率之和为1，即这样，可以把每个模型的后验概率看作是权重wk；
[0095] 假设p(r|fk)为高斯分布， 为该高斯分布的平均值，σ2K为方差，用参数和g(·)表示相关的高斯分布为：
[0096] p(r|fk)＝g(r|θk)
[0097]
[0098] 得到的估计值为r的条件期望，用期望-最大算法对代价函数进行优化得到权重wk。
[0099] 所述步骤1)通过植被覆盖度Fc来区分植被和土壤分配的全球陆表可利用能量，利用MODIS遥感数据中植被吸收的光合有效辐射比例FPAR产品确定Fc；
[0100] 所述步骤2)订正的PM遥感蒸散算法中，统一的水分限制因子VPD的最大值定为2.79KPa；
[0101] 所述步骤3)Priestly-Taylor蒸散算法中的Rns可以用地表净辐射Rn和叶面积指数LAI计算，具体为Rns＝Rnexp(－kLAI)，k为经验系数；
[0102] 所述步骤4)改进的Priestly-Taylor蒸散算法的土壤水分限制因子中，DT是通过日昼夜温差计算得到或通过地面温度计算得到；
[0103] 所述步骤5)中经验系数a0，a1，a2，和a3是利用全球通量站点观测数据作为训练数据回归得到的；
[0104] 所述步骤6)的观测数据提取模块中，站点实测数据样本质量是通过质量标识判断而得到的；
[0105] 所述步骤7)从MODIS遥感数据中选取叶面积指数、植被指数，从MERRA再分析气象数据中选取基本气象参数，包括温度、相对湿度，驱动各个算法模块计算得出各自的蒸散量。
[0106] 本发明所述一种基于多算法集成的全球陆表蒸散估算系统及方法配合使用时，其操作步骤具体如下：
[0107] 1、根据MODIS蒸散算法(Mu et al.2011)，利用MODIS的LAI、NDVI以及EVI(增强型植被指数)产品结合再分析MERRA气象资料，估算区域蒸散量；
[0108] 2、根据订正的PM蒸散算法(Yuan et al.2010)，利用MODIS的LAI、NDVI以及EVI产品结合再分析MERRA气象资料，估算区域蒸散量；
[0109] 3、根据Priestly-Taylor蒸散算法(Fisher et al.2008)，利用MODIS的LAI、NDVI以及EVI产品结合再分析ERRA气象资料，估算区域蒸散量；
[0110] 4、根据改进的Priestly-Taylor蒸散算法(Yao et al.2013)，利用MODIS的LAI、NDV I以及EVI产品结合再分析ERRA气象资料，估算区域蒸散量；
[0111] 5、根据改进的经验性遥感蒸散算法(Wang et al.2008)，利用MODIS的LAI、NDVI以及EVI产品结合再分析ERRA气象资料，估算区域蒸散量；
[0112] 6、通过全球通量观测站点数据质量提取地面蒸散观测值；
[0113] 7、建立贝叶斯模型方法，将不同来源的估算结果融合在一起，得到最优估计结果。r表示需要估算的变量，在某一时刻的该变量的观测值为rt，模型集合K，即{f1,f2,…,fk}，为估算r的所有单个模型的集合。根据全概率定律，需要被估算算的蒸散r的概率密度函数可以由每个模型的概率密度函数表示如下：
[0114]
[0115] 其中，p(r|fk)是单个模型fk估算蒸散r的概率密度函数，p(fk|rt)是单个模型fk的的后验概率，能够反应单个模型fk与观测数据符合的程度。所有单个模型的后验概率之和为1。这样，可以把每个模型的后验概率看作是权重wk,，这个权重就是需要估算的参数。可以表达为：
[0116]
[0117] 因此，得到的估算值的概率密度函数p(r|f1,f2,…,fk)就是单个模型概率密度函数的加权平均，权重值是各个函数的后验概率。假设p(r|fk)为高斯分布，fkA为该高斯分布2 2
的平均值，σK为方差，用参数θK＝{fkA，σK}和g(·)表示相关的高斯分布，如下：
[0118] p(r|fk)＝g(r|θk)
[0119]
[0120] 得到的估计值为可以表示为r的条件期望，如下：
[0121]
[0122] 很显然，fkA是单个模型的估算值，那么问题的关键则是如何估算模型的后验概率wk，使得到的r概率密度函数与观测值更加接近。根据贝叶斯理论，已知定观测数据，当似然函数公式(4)值最大时，贝叶斯方法得到的估算值为真值的概率最大(Duan and Phillips,2010)。为了方便计算，应对该似然公式的对数函数进行优化：
[0123]
[0124] 其中，∑(s,T)表示在所有空间点s和所有时间点T的观测之和，rs,t表示某一时间点T和空间s上，用期望-最大算法(expectation-maximization，EM)算法对代价函数进行优化(Raftery et al.2005)，单个模型fkA的估算值和相应的观测值rs,t作为训练数据计算每个模型的后验概率。
[0125] 8、选择测点
位置的测量蒸散量与实际蒸散量进行相关分析，得到研究区内两者的关系；
[0126] 9、评价模型在研究区的世界蒸散估算状况。
[0127] 验证数据集为全球FLUXNET200个通量站点观测数据。全球FLUXNET站点数据包括了200个站点，分布亚洲、欧洲、北美洲、南美洲、非洲和大洋洲。横跨热带、亚热带、温带、寒带各种
气候带，覆盖了沙漠、绿洲、草原、森林、萨瓦娜、农田和湿地等各种植被覆盖类型，满足至少有一年的观测数据，能够满足各种植被类型和气候带条件下蒸散的估算要求。以全球陆表为研究区，该研究区的通量观测数据是2000-2009年时间段收集的，相应的MODIS数据和MERRA再分析数据也是2000-2009年时间段的。为了验证模型的蒸散量估算效果，计算得到蒸散量量与全球通量观测站的观测数据进行相关分析。根据全球FLUXNET各个观测点的测量的蒸散量和贝叶斯算法集成估算的蒸散量建立散点图，其结果如图3，平均误差2 2 2
为2.3W/m，相关系数平方(R)为0.64，均方差为24.2W/m。结果表明，模型估计值和实地观测数据具有较高的相关性。通过分析可见，贝叶斯模型集成算法对估算全球陆表蒸散量很有效。
