CN115169263B - 基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法 - Google Patents

Abstract

本发明公开了一种基于TVGM‑Budyko耦合模型的多时间尺度径流变化归因分析方法。包括：步骤一、收集研究区域内水文、气象数据及土壤含水量数据；步骤二、数据前处理；步骤三、耦合TVGM模型和Budyko框架，将Budyko框架时间尺度由多年平均尺度扩展到月尺度及年际尺度；步骤四、对构建的耦合模型效果进行评估；步骤五、利用构建的耦合模型开展研究区多时间尺度径流变化归因分析。本发明通过耦合TVGM模型与Budyko框架，将Budyko框架的时间尺度由多年平均尺度拓展到月尺度及年际尺度，有助于准确分析不同时间尺度下的径流演变规律，进而开展多时间尺度径流变化归因分析，可为短期、中长期区域水资源管理提供科学依据。

Description

基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析 方法

技术领域

本发明涉及多时间尺度径流模拟与变化归因分析领域，具体涉及一种基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法。

背景技术

径流是水文循环的核心环节，受气候变化和高强度人类活动的影响，流域径流在不同时间尺度(如月、年和多年平均尺度)上呈现出不同的变化特征。如何准确辨识不同时间尺度下的径流演变规律并进行归因分析，是流域水循环研究的重点和难点之一。

当前径流模拟主要通过模型模拟实现，目前常用的模型分为数据驱动模型和物理过程驱动模型。数据驱动模型主要是指基于大量与径流相关的变量数据，根据人工智能或机器学习等方法对流域径流进行模拟预测，该方法具有运行速度快，效率高等特点，但其缺乏对流域水文过程物理机制的描述。物理过程驱动模型主要包括概念性集总水文模型和分布式水文模型，其具有坚实的物理基础和明确的物理过程，可以准确模拟不同时间尺度流域出口断面径流过程。Budyko框架是当前径流变化归因研究中常用的方法，该方法结构简单，参数少且易于获取，可结合弹性系数法准确推导出径流对多变量的弹性解析式，进而分离和量化多要素对径流的影响，但其适用于多年平均尺度，难以满足月或年际尺度流域径流变化归因研究。

发明内容

本发明的目的在于克服以上技术困难，提供一种基于时变增益模型(TVGM)和Budyko框架的多时间尺度径流变化归因分析方法，旨在解决不同时间尺度上的径流变化归因分析。

对于多年平均尺度而言，流域蓄水量变化可忽略不计，Budyko框架正是基于此假设。而对于月或年际尺度而言，流域蓄水变化不可忽略，因而拓展Budyko框架的重点在于对流域蓄水量变化的模拟。TVGM是一种结构简单且参数较少的流域水文模型，可模拟获取不同时间尺度流域土壤含水量。通过将TVGM与Budyko框架耦合，可以拓展Budyko框架的时间尺度，进而实现不同时间尺度径流变化归因分析，有助于流域水资源合理开发和科学管理。综上所述，耦合流域水文模型和Budyko框架是提升Budyko框架适用性，开展多时间尺度径流变化归因分析的重要方法。如何将流域水文模型与Budyko框架进行耦合是目前本领域技术人员亟需解决的关键问题。

本发明的目的具体是通过以下技术方案实现的：

一种基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，所述方法包括以下步骤：

步骤一、数据收集：收集研究区内径流数据、气象数据、土壤含水量数据；基于气象数据，利用FAO Penman-Monteith方程计算研究区在研究时段内的潜在蒸散发量；

步骤二、数据前处理：包括缺失数据的插补延长、数据空间插值和重采样、研究区数据面平均值计算中的一项或几项；

步骤三、耦合时变增益模型TVGM模型和Budyko框架：

3-1确定流域水量平衡方程并构建多时间尺度流域Budyko水热平衡统一框架形式；

3-2构建研究区域的TVGM模型；

3-3将构建好的TVGM模型代入到Budyko水热平衡方程，由此获得TVGM-Budyko耦合模型的解析式；

步骤四、对耦合模型的效果进行评估：基于研究区实测径流或蒸散发资料对TVGM-Budyko耦合模型的效果进行评估；

步骤五、多时间尺度径流变化归因分析：基于弹性系数法定量分析研究区域内不同时间尺度各影响因素对径流变化的影响。

进一步的，步骤一中所述气象数据包括气温、气压、相对湿度、风速风向、日照及辐射数据。

进一步的，步骤一中所述的FAO Penman-Monteith方程的具体计算步骤为：

式中，PE为蒸散发量，mm/d；Δ为饱和水汽压斜率，kPa/℃；R_n为地表净辐射，MJ/(m²*d)；G为土壤热通量，MJ/(m²*d)；γ为干湿表常数，kPa/℃；T_mean为日平均温度，℃；u₂为2米高处风速，m/s；e_s为饱和水汽压，kPa；e_a为实际水汽压，kPa。

进一步的，步骤二中所述的缺失数据的插补延长为线性内插法和线性外延法(优选利用窗口长度为5以内的已知值对缺失值进行线性内插或延长)，所述研究区数据面平均值计算采用泰森多边形法计算得到。泰森多边形法，其具体计算过程如下：

利用Gis软件将研究区内的各要素值提取到气象站点位置；将研究区域内所有相邻气象站连成三角形；依次做这些三角形各边的垂直平分线；将每个气象站点周围的若干垂直平分线围城一个多边形，该图形称作泰森多边形；用泰森多边形内包含的唯一气象站点处的要素值表示该多边形区域内的平均值；区域内所有的泰森多边形构成多边形网，将每个多边形在多边形网的权重系数与要素区域平均值相乘后累加可获得要素在整个研究区域的面平均值：

X＝f₁X₁+f₂X₂+L f_kX_k

其中，k为站点个数，X₁,X₂,…,X_j为各站点处的要素值，f₁,f₂,…,f_j为各多边形的权重系数。

进一步的，步骤三3-1中所述的流域水量平衡方程的具体形式如下：

P-ΔS＝Q_d+Q_b+E

P为降水/mm，△S为流域蓄水变化/mm，Q_d为地表径流/mm，Q_b为基流/mm，Q＝Q_d+Q_e，E为蒸散发/mm。

进一步的，步骤三3-1中所述的多时间尺度流域水热平衡统一框架具体形式如下：

Q＝(1-B_k(A,n))*P_e

其中，Q为径流量，P_e为有效降水量，单位均为mm；A为流域干旱指数，其值为潜在蒸散发量除以有效降水量，在多年平均尺度下，有效降水量约等于降水量，因此A＝Φ，Φ为多年平均干旱指数，A和Φ均为无量纲变量；n为Budyko参数，默认值为2.6。B_k为Budyko函数，其有多种形式，本发明具体实施方式中采用的是Choudhury-Yang型Budyko函数，其具体形式如下：

因此，本发明采用的Budyko水热耦合平衡统一框架可写为：

其中，Q、P_e和ET₀分别为径流、有效降水和潜在蒸散发，单位均为mm。n为Budyko参数，默认值为2.6。

进一步的，步骤三中所述的TVGM模型产流模块的具体形式如下：

Q_b＝k_r·AW

其中，g₁和g₂为时变增益因子的有关参数，C为与植被覆盖度有关的参数，均为无量纲参数；AW为土壤含水量，AWM为饱和土壤含水量，单位均为mm；k_r为土壤水出流系数，为无量纲参数。其余符合意义同上。

进一步的，步骤三中所述的耦合TVGM模型的多时间尺度Budyko水热平衡模型即TVGM-Budyko耦合模型的解析式形式如下：

对于月尺度水热平衡模型而言，径流被分为地表径流和基流，流域蓄水变化不可忽略：

Q_i＝Q_di+Q_bi

Q_i＝(1-B_k(A_i,n))*P_ei

在本发明的具体实施方式中：

ΔS_i＝AW_i-AW_i-1

其中，Q_i为第i月径流量/m，Q_di和Q_bi分别为第i月地表径流和基流量/mm；A_i为第i月流域干旱指数，其值为第i月潜在蒸散发量ET_0i除以第i月有效降水量P_ei，为无量纲变量；P_ei为第i月有效降水量，其值为第i月降水量P_i减去第i-1月至第i月的流域蓄水变化量△S_i；ET_0i为第i月潜在蒸散发量；AW_i为第i月土壤含水量，AW_i-1表示第i-1月土壤含水量/mm；AWM为流域饱和土壤含水量。上述变量单位均为mm。n和C均为与流域特性有关的无量纲参数，g₁和g₂为时变增益因子有关的无量纲参数。

对于年际尺度水热平衡模型而言，其Budyko方程形式与月尺度基本相同：

Q_j＝(1-B_k(A_j,n))*P_ej

在本发明的具体实施方式中：

ΔS_j＝AW_j-AW_j-1

其中，Q_j、Q_d,j和Q_b,j分别为第j年总径流量/mm、直接径流量和基流量/mm，由月径流量累加而得；P_ej为第j年有效降水量/mm，其值为第j年降水量P_j减去第j-1年至第j年的流域蓄水变化量△S_j；P_j,i为第j年第i月降水量/mm；ET_0j为第j年潜在蒸散发量/mm；AW_j,i为第j年第i月土壤水含量，AW_j-1为第j-1年土壤含水量/mm；AW_j,0表示上一年12月份的土壤水含量；AW_j为第j年土壤含水量；AWM为流域饱和土壤含水量。上述各变量单位均为mm。n和C均为与流域特性有关的无量纲参数，g₁和g₂为时变增益因子有关的无量纲参数。

对于多年平均尺度水热平衡模型而言，流域基流和流域蓄水变化均可忽略不计：

式中，

和

分别为多年平均径流量和降水量，单位mm；Φ为多年平均干旱指数，其值为多年平均潜在蒸散发量

除以多年平均降水量

为无量纲变量。n为与流域特性有关的Budyko参数。

在本发明具体实施方式为：

式中，

和

为多年平均径流量、降水量和潜在蒸散发量，单位mm。

进一步的，步骤四所述的模型效果评估指标为：相关系数(r)、相对误差(R_bias)、均方根误差(RMSE)和Nash效率系数(NSE)和Kling-Gupta效率系数(KGE)。其具体计算过程如下：

其中，x_i表示估计值，y_i表示实测值，

和

分别表示x和y的平均值，m为序列个数，r为相关性系数，

γ为x_i和y_i标准差的比值。

进一步的，步骤五所述的弹性系数法具体计算过程如下：

其中，E、P、ET₀和n的意义同前文所述；S_{Q_P}和

分别表示径流对降水和蒸散发的弹性。

本发明的优点和有益效果是：耦合时变增益模型与Budyko框架，将Budyko框架的时间尺度由多年平均尺度拓展到月尺度及年际尺度，扩大了Budyko框架的应用范围，可适用于不同时间尺度径流变化归因分析，对于流域水资源科学管理和综合利用具有重要意义。

附图说明

下面结合附图和实施例对本发明作进一步说明。

图1为本发明实施例的流程示意图；

图2为研究区径流变化趋势图；

图3为TVGM参数率定及验证结果；

图4～图6分别为月、年和多年平均尺度径流模拟结果；

图7～图9分布为月、年和多年平均尺度径流演变归因分析结果。

具体实施方式

一种基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，包括以下步骤：

步骤一、数据收集：收集目标研究区在研究时段内的气象数据、径流数据及土壤含水量数据；利用气象数据及FAO Penman-Monteith方程计算目标研究区在研究时段内的潜在蒸散发。其中气象要素数据包括降水、气温、气压、相对湿度、风速风向、日照及辐射数据；

FAO Penman-Monteith方程的具体计算步骤为：

假设作物植株高度为0.12m，固定的作物表面阻力为70m/s，反射率为0.23，即认为非常类似于表面开阔、高度一致、生长旺盛、完全遮盖地面而水分充分适宜的绿色草地的蒸散量为潜在蒸散量：

式中，PE为蒸散发量(mm/d)；Δ为饱和水汽压斜率(kPa/℃)；R_n为地表净辐射(MJ/(m²*d))；G为土壤热通量(MJ/(m²*d))；γ为干湿表常数(kPa/℃)；T_mean为日平均温度(℃)；u₂为2米高处风速(m/s)；e_s为饱和水汽压(kPa)；e_a为实际水汽压(kPa)。

式中关键参数的计算如下：

步骤1-1，计算Δ：

其中，T选用日平均气温；

步骤1-2，计算R_n：

R_n＝R_ns-R_nl (3)

即净辐射R_n为收入的短波辐射R_ns和支出的净长波辐射R_nl之差；

R_ns＝(1-α)R_s (4)

其中，α＝0.23，此处取绿色草地参考作物的反照率；

R_s为太阳辐射；a_s和b_s表示回归系数，推荐使用a_s＝0.25，b_s＝0.50；n为实际日照时数(h)，N为最大可能日照时数(h)；

w_s为日出时角；

表示纬度(rad)；

δ为太阳磁偏角；J为日序，取值范围为1到365或366，1月1日为1；

G_sc为太阳常数(MJ/(m²*min)),G_sc＝0.0820，

和δ意义同上，R_a为日地球外辐射；

d_r为日地平均距离，J的含义同上；

其中σ为斯蒂芬-玻尔兹曼常数(MK·K^-4·m^-2·d^-1)，σ＝4.903×10^-9；

T_max,K＝T_max+272.15 (12)

T_min,K＝T_min+272.15 (13)

T_max,K、T_min,K为最高、最低绝对温度(K)，T_max、T_min为日最高、最低温度(℃)；

e_a＝RH×e_s (14)

e_a为实际水汽压，RH为相对湿度，e_s为饱和水汽压；

T_max、T_min意义同上；

R_so＝(a_s+b_s)R_a (18)

R_a、a_s、b_s意义同上；

步骤1-3，计算G：

G≈0 (19)

步骤1-4，计算γ：

γ＝0.665×10³P (20)

其中，P为大气压；

步骤1-5，计算T_mean：

其中，T_max、T_min意义同上；

步骤1-6，计算u₂：

其中，u_z为10m高度处的风速，z＝10。

步骤二、数据前处理：数据前处理包括缺失数据的插补延长、数据空间插值和重采样及数据面平均值计算等；缺失数据利用线性内插或是线性外延方法进行插补延长；数据的空间插值采用反距离加权平均方法；数据面平均值计算方法为泰森多边形法，其具体步骤为：利用Gis软件将研究区内的各要素值提取到气象站点位置；将研究区域内所有相邻气象站连成三角形；依次做这些三角形各边的垂直平分线；将每个气象站点周围的若干垂直平分线围城一个多边形，该图形称作泰森多边形；用泰森多边形内包含的唯一气象站点处的要素值表示该多边形区域内的平均值；区域内所有的泰森多边形构成多边形网，将每个多边形在多边形网的权重系数与要素区域平均值相乘后累加可获得要素在整个研究区域的面平均值：

X＝f₁X₁+f₂X₂+L f_kX_k (23)

其中，k为站点个数，X₁,X₂,…,X_j为各站点处的要素值，f₁,f₂,…,f_j为各多边形的权重系数。

步骤三、耦合时变增益模型TVGM模型与Budyko框架：

步骤3-1：确定流域水量平衡方程并构建多时间尺度流域Budyko水热平衡统一框架形式，流域水量平衡方程具体形式为：

P-ΔS＝Q_d+Q_b+E (24)

其中，P为降水，△S为流域蓄水变化，Q_d为地表径流，Q_b为基流，Q＝Q_d+Q_e，Q为径流量，E为蒸散发，上述各变量单位均为mm。

多时间尺度流域Budyko水热平衡统一框架具体形式为：

其中，P_e为有效降水，P_e＝P-△S，单位mm。其余各符合意义同上文所述

步骤3-2：构建研究区域的TVGM模型产流模块；

Q_b＝k_r·AW (27)

其中，g₁和g₂为时变增益因子的有关参数，C为与植被覆盖度有关的参数；AW为表层土壤土壤湿度，AWM为表层土壤饱和含水量，单位mm；k_r为土壤水出流系数。

步骤3-3：将构建好的TVGM代入到Budyko方程，由此获得TVGM-Budyko耦合模型的解析式；

对于月尺度水热平衡模型而言，径流被分为地表径流和基流，流域蓄水变化不可忽略：

Q_i＝Q_di+Q_bi (28)

ΔS_j＝AW_j-AW_j-1 (32)

其中，Q_i为第i月径流量，Q_di和Q_bi分别为第i月地表径流和基流量；P_ei为第i月有效降水量，其值为第i月降水量P_i减去第i-1月至第i月的流域蓄水变化量△S_i；ET_0i为第i月潜在蒸散发量；AW_i为第i月土壤含水量，AWM为流域饱和土壤含水量。上述变量单位均为mm。n和C均为与流域特性有关的无量纲参数，g₁和g₂为时变增益因子有关的无量纲参数。

对于年际尺度水热平衡模型而言，其Budyko方程形式与月尺度基本相同：

ΔS_j＝AW_j-AW_j-1 (36)

其中，Q_j、Q_d,j和Q_b,j分别为第j年总径流量、直接径流量和基流量，由月径流量累加而得；P_ej为第j年有效降水量，其值为第j年降水量P_j减去第j-1年至第j年的流域蓄水变化量△S_j；P_j,i为第j年第i月降水量；ET_0j为第j年潜在蒸散发量；AW_j,i为第j年第i月土壤水含量，AW_j,0表示上一年12月份的土壤水含量；AW_j为第j年土壤含水量；AWM为流域饱和土壤含水量。上述各变量单位均为mm。n和C均为与流域特性有关的无量纲参数，g₁和g₂为时变增益因子有关的无量纲参数。

对于多年平均尺度水热平衡模型而言，流域基流和流域蓄水变化均可忽略不计：

式中，

和

为多年平均径流量、降水量和潜在蒸散发量，单位mm。n为与流域特性有关的Budyko参数。

步骤四、利用实测径流数据对多时间尺度Budyko框架效果进行评估：评估指标有相关系数(r)、相对误差(R_bias)、均方根误差(RMSE)和Nash效率系数(NSE)、Kling-Gupta效率系数(KGE)。其具体计算过程如下：

其中，x_i表示估计值，y_i表示实测值，

和

分别表示x和y的平均值，m为序列个数，r为相关性系数，

γ为x_i和y_i标准差的比值。

步骤五、定量归因研究区域不同时间尺度径流演变：基于弹性系数法，定量分析目标研究区域内不同时间尺度各影响因素对径流变化的影响。

弹性系数法具体计算过程如下：

其中，E、P、ET₀和n的意义同前文所述；S_{Q_P}和S_{Q_ET0}分别表示径流对降水和蒸散发的弹性。

实施例1：

以某流域为案例进行分析，对本发明技术方法的具体应用作进一步说明，具体应包括以下步骤：

步骤一、数据收集：

研究区内气象数据从中国气象局国家气象信息中心(NMIC)(http://cdc.cma.gov.cn)下载。利用气象数据根据FAO Penman-Monteith方程计算潜在蒸散发量，主要是利用降水、气温、辐射、气压、相对湿度、风速及日照时数。水文站点的径流数据从水文年鉴中获取。土壤含水量数据为遥感再分析数据，从国家青藏高原科学数据中心获取。数据具体介绍如表1所示。

表1数据准备

序号	数据类型	时间分辨率	空间分辨率	数据类型
					1	气象数据	逐日	-	站点
2	径流数据	逐月	-	站点
					3	土壤含水量数据	逐月	0.05°	栅格

步骤二、数据前处理：

采用线性内插法和线性外延法对缺失的降水、径流数据进行插补延长；采用反距离加权平均法对栅格数据进行空间插值，使各要素的空间分辨率一致，并计算栅格数据的面平均值；对于气象数据，利用泰森多边形计算研究区数据面平均值；

步骤三、耦合时变增益模型和Budyko框架：

对选取流域的月、年尺度径流序列进行分析，结果表明：多年平均月径流统计结果显示(图2)，该研究流域径流主要集中在每年的7月、8月、9月，年内变化明显；年径流趋势分析显示，自2005–2018年，该研究流域径流呈下降趋势，下降率约为5.34mm/年。

利用2005年1月–2008年12月逐月实测降水、径流及土壤含水量数据，对TVGM模型进行参数率定，2009年1月–2010年12月逐月实测降水、径流及土壤含水量数据对构建的TVGM进行验证，结果如图3所示：率定期R²和NSE分别为0.79和0.78，检验期R²和NSE分别为0.87和0.83。无论是率定期还是检验期，其R²和NSE均超过0.75，表明TVGM模型可用于该流域的径流模拟。

其次将具有固定参数的TVGM模型耦合到Budyko框架中，根据初始条件，迭代求解研究区逐月及逐年径流量。多年平均尺度Budyko框架选用5年滑动平均数据进行计算，2005-2018年共计10个时间段。

步骤四、利用实测数据对多时间尺度Budyko框架的效果进行评估：

利用耦合TVGM的Budyko框架对2005年1月–2018年12月的逐月径流进行模拟，结果如图4所示，R²＝0.76，NSE＝0.75，均大于(或等于)0.75，表明耦合TVGM的Budyko框架可以很好地模拟月径流过程；利用耦合TVGM的Budyko框架对2005–2018年的逐年径流进行模拟，结果如图5所示，R²＝0.83，NSE＝0.77，均大于0.75，表明耦合TVGM的Budyko框架可以很好地模拟年径流过程；利用Budyko框架对2005–2018年5年滑动平均径流进行模拟，结果如图6所示，R²＝0.88，NSE＝0.79，均大于0.75，表明Budyko框架可以很好模拟多年平均径流过程。综上，耦合TVGM模型的Budyko框架将Budyko框架的时间应用尺度从多年平均尺度拓展到月和年尺度，并能取得较好的模拟效果，说明耦合TVGM模型的Budyko框架具有较好的实际应用价值。

步骤五、多时间尺度径流演变归因分析：

利用弹性系数法开展不同时间尺度径流变化归因分析，结果如图7-9所示。结果表明，在月和年时间尺度，降水是导致径流变化的主要原因，潜在蒸散发对径流演变的贡献可以忽略不计：降水增加1mm将导致径流增加约0.93mm；潜在蒸散发增加1mm将导致径流减少约0.07mm。在多年平均尺度上，降水和潜在蒸散发均是导致径流变化的主要原因：降水增加1mm将导致径流增加约0.67mm，潜在蒸散发减少1mm将导致径流增加0.23mm。

最后应说明的是，以上仅用以说明本发明的技术方案而非限制，尽管参照较佳布置方案对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围。

Claims (10)

1.一种基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，所述方法包括以下步骤：

步骤一、数据收集：收集研究区内径流数据、气象数据、土壤含水量数据；基于气象数据，利用FAO Penman–Monteith方程计算研究区在研究时段内的潜在蒸散发量；

步骤二、数据前处理：包括缺失数据的插补延长、数据空间插值和重采样、研究区数据面平均值计算中的一项或几项；

步骤三、耦合TVGM模型和Budyko框架：

3-1确定流域水量平衡方程并构建多时间尺度流域Budyko水热平衡统一框架形式；

3-2构建研究区域的TVGM模型；

3-3将构建好的TVGM模型代入到Budyko水热平衡方程，由此获得TVGM-Budyko耦合模型的解析式；

步骤四、对耦合模型的效果进行评估：基于研究区实测径流或蒸散发资料对TVGM-Budyko耦合模型的效果进行评估；

步骤五、多时间尺度径流变化归因分析：基于弹性系数法定量分析研究区域内不同时间尺度各影响因素对径流变化的影响。

2.根据权利要求1所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤一中所述气象数据包括降水、气温、气压、相对湿度、风速、日照及辐射数据。

3.根据权利要求1所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤一中FAO Penman–Monteith方程的具体计算公式为：

式中，PE为蒸散发量，mm/d；Δ为饱和水汽压斜率，kPa/℃；R_n为地表净辐射，MJ/(m²*d)；G为土壤热通量，MJ/(m²*d)；γ为干湿表常数，kPa/℃；T_mean为日平均温度，℃；u₂为2米高处风速，m/s；e_s为饱和水汽压，kPa；e_a为实际水汽压，kPa。

4.根据权利要求1所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤二中所述的缺失数据的插补延长为线性内插法和线性外延法，所述研究区数据面平均值计算采用泰森多边形法计算得到。

5.根据权利要求1所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤三3-1中所述的流域水量平衡方程的具体形式如下：

P-ΔS＝Q_d+Q_b+E

式中，P为降水/mm，△S为流域蓄水变化/mm，Q_d为地表径流/mm，Q_b为基流/mm，E为蒸散发/mm。

6.根据权利要求5所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤三3-1中所述的多时间尺度流域Budyko水热平衡统一框架具体形式如下：

Q＝(1-B_k(A,n))*P_e

式中，Q为径流量，P_e为有效降水量，P_e＝P-△S，单位均mm；A为流域干旱指数，A＝ET₀/P_e；ET₀为潜在蒸散发，单位mm；B_k为Budyko函数；n为Budyko参数；在多年平均时间尺度上，流域蓄水变化忽略不计，流域有效降水量认为是降水量，即P_e＝P，A＝Φ＝ET₀/P，Φ为多年平均尺度干旱指数。

7.根据权利要求6所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤三3-2中TVGM模型产流模块的具体形式如下：

Q_b＝k_r·AW

式中，g₁和g₂为时变增益因子有关参数，无量纲；C为与植被覆盖度有关的参数，无量纲；AW为表层土壤的土壤湿度/mm；AWM为表层土壤饱和含水量/mm；k_r为土壤水出流系数。

8.根据权利要求7所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤三3-3中所述TVGM-Budyko耦合模型的解析式形式如下：

对于月尺度水热平衡模型，径流被分为地表径流和基流，流域蓄水变化不可忽略，具体形式如下所示：

Q_i＝Q_di+Q_bi

Q_i＝(1-B_k(A_i,n))*P_ei

ΔS_i＝AW_i-AW_i-1

其中，Q_i为第i月径流量/mm，Q_di和Q_bi分别为第i月地表径流和基流量/mm；P_ei为第i月有效降水量/mm，其值为第i月降水量P_i减去第i-1月至第i月的流域蓄水变化量△S_i；AW_i表示第i月土壤含水量/mm，AW_i-1表示第i-1月土壤含水量/mm；AWM为流域饱和土壤含水量/mm；A_i为第i月流域干旱指数，其值为第i月潜在蒸散发量ET_0i除以第i月有效降水量P_ei，为无量纲变量；

对于年际尺度水热平衡模型，具体形式如下式所示：

Q_j＝(1-B_k(A_j,n))*P_ej

ΔS_j＝AW_j-AW_j-1

其中，Q_j、Q_d,j和Q_b,j分别为第j年总径流量/mm、直接径流量和基流量/mm，由月径流量累加而得；P_ej为第j年有效降水量/mm，其值为第j年降水量P_j减去第j-1年至第j年的流域蓄水变化量△S_j；P_j,i为第j年第i月降水量/mm；AW_j,i为第j年第i月土壤水含量/mm；AW_j为第j年土壤含水量/mm；AW_j-1为第j-1年土壤含水量/mm；AWM为流域饱和土壤含水量/mm；A_j为第j年流域干旱指数，其值为第j年潜在蒸散发量ET_0j除以第j年有效降水量P_ej，为无量纲变量；

对于多年平均尺度水热平衡模型而言，流域基流和流域蓄水变化均可忽略不计，具体形式如下式所示：

式中，

和

分别为多年平均径流量和降水量，单位mm；Φ为多年平均干旱指数，其值为多年平均潜在蒸散发量

除以多年平均降水量

为无量纲变量。

9.根据权利要求1所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤四中所述的对耦合模型的效果进行评估的方法具体为：利用相关系数r、相对误差R_bias、均方根误差RMSE和Nash效率系数NSE和Kling-Gupta效率系数KGE作为评估指标进行评估。

10.根据权利要求5所述的基于TVGM-Budyko耦合模型的多时间尺度径流变化归因分析方法，其特征在于，步骤五所述的弹性系数法的具体计算过程如下：

其中，S_{Q_P}和

分别表示径流对降水和潜在蒸散发的弹性，Q为径流量，ET₀为潜在蒸散发。

Images