武汉大学水资源与水电工程科学国家重点实验室，湖北 武汉

收稿日期：2021年1月2日；录用日期：2021年2月1日；发布日期：2021年2月22日

摘要

陆气耦合过程复杂，研究土壤湿度对降水的影响机制，对气候预测和天气预报具有重要意义。目前土壤湿度对降水的反馈机制尚不明确，采用物理模型存在较大的不确定性，而通过数据直接反映非线性统计相关性的转移熵，是一种可解释水文系统的新范式，为研究此类问题提供了可能。本文提出并验证了归一化转移熵，使量纲不同的耦合强度具备可比性；使用长江流域2002~2018年168个气象站点观测的降水和蒸散发数据、全球陆面参数数据LPDR V2.0中的土壤湿度数据、以及MODIS中的植被指数NDVI数据，采用显著滞时、不同滞时的相对预测度和归一化信息熵3个信息指标识别土壤湿度对降水的反馈特征，并应用偏相关方法验证了长江流域土壤湿度对降水正反馈的合理性。结果表明：1) 土壤湿度对降水的平均显著滞时为4.3 d，与降水对土壤湿度的1.8 d相比滞时更长，说明土壤湿度对降水的反馈过程存在滞时，即土壤湿度响应前期降水后，延迟较长时段影响降水过程；2) 降水对土壤湿度的归一化转移熵为0.51，与土壤湿度对降水的0.13相比耦合强度更高，说明土壤湿度对降水的反馈显著弱于降水对土壤湿度的影响；3) 金沙江下游和洞庭湖、鄱阳湖、太湖水系等近湖区域由于水量交换频繁，土壤湿度对降水反馈较快且强度更高。

关键词

土壤湿度对降水反馈，长江流域，转移熵

Soil Moisture Feedbacks on Precipitation in the Yangtze River Catchment Based on Transfer Entropy

Wei Lou, Pan Liu^*, Lei Cheng, Kaijie Zou, Qian Xia

State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan Hubei

Received: Jan. 2^nd, 2021; accepted: Feb. 1^st, 2021; published: Feb. 22^nd, 2021

ABSTRACT

The land surface and atmosphere interact as a complexly linked system. The soil moisture-precipitation feedback mechanism is of great significance to weather forecast and climate prediction; however, it is uncertain in physical-based numerical simulation. Transfer entropy, which reflects non-linear statistical correlation directly through data, has been considered as a new paradigm to explain hydrological system and provides possibility to study the land-atmosphere coupling. Normalized transfer entropy was proposed in this study, and was verified to measure comparable coupling strength. Observation data on 168 weather stations during 2002~2018, MODIS-NDVI (MOD13A2) and soil moisture data from the global land parameters data record LPDR V2.0 were used. Characteristics of the feedback were explored by three information indexes. A positive soil moisture-precipitation feedback was also valid by the partial correlation. The results showed: 1) the average significant lag time of soil moisture-precipitation coupling was 4.3 d, which was longer than 1.8 d of precipitation-soil moisture coupling. It indicated that there was a delay effect during the feedback process. The soil moisture in response to precipitation affected the precipitation process after a longer period of time; 2) the normalized transfer entropy of precipitation-soil moisture coupling was 0.51, which was stronger than 0.13 of soil moisture-precipitation coupling; and 3) due to the frequent water exchange in the lower reaches of the Jinsha River, the Dongting Lake, the Poyang Lake and the Taihu Lake systems, soil moisture had a faster and stronger feedback on precipitation.

Keywords:Soil Moisture-Precipitation Feedback, The Yangtze River, Transfer Entropy

Copyright © 2021 by author(s) and Wuhan University.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

复杂的陆气反馈过程中，土壤湿度对大气降水的反馈机制是核心问题之一。一方面，土壤湿度作为陆面水文过程的重要状态变量，直接被大气降水驱动；另一方面，土壤湿度通过向大气补充可降水汽的直接方式，或改变大气条件的间接方式，从而触发降水事件发生。识别土壤湿度对降水的反馈关系特征，对搭建灵活易用的陆气耦合模型、提高水文模拟和气象预报精度有重要意义 [1]。

土壤湿度对降水反馈在机制上存在较多不确定性，指标较多且同质 [2] [3]，即使在同一研究区域，研究方法、数据尺度和评价指标不同，也会导致争议性结论。如伊利诺伊州，Findell和Eltahir [4] 基于插值的伊利诺伊气象站网观测数据，以滞后相关性统计方法识别出春末夏初存在土壤湿度对降水的反馈作用；而Salvucci、Saleem和Kaufmann [5] 认为插值存在人为正向偏差，并使用格兰杰因果检验证实伊利诺伊州土壤湿度对降水的反馈作用并不显著。

长江流域位于陆–气–海系统复杂耦合的东亚季风区，有58%的降水水汽来源于流域内陆地蒸散发，相较于海面，陆地土壤湿度为流域降水输送了更多的水汽 [6]。全球陆气耦合试验GLACE也表明，东亚季风区是土壤湿度对降水反馈的热点区域之一 [7]。K. C. Chow [8] 通过区域气候模式(RCM)做敏感性试验，证明青藏高原春季土壤湿度的增加会导致长江流域夏季降水量的增加；而王林、王磊等 [9] 用天气预测模型(WRF)进行数值模拟得到相反结论，春季高原土壤湿度偏湿，长江中下游地区初夏降水量大范围减小。

数值模拟和统计分析是研究土壤湿度对降水反馈的两种常用手段，数值模拟模型具有较大的不确定性且需庞大的建模工作，而统计分析中的指标难以度量非线性关系 [10]。基于香农熵的信息论可检验非线性因果，近年被认为是一种可解释水文系统的新范式 [11]。如用互信息确定水文气候变量的时空变异性 [12]、时间信息划分网络构建生态水文网络去推断系统的响应行为 [13] [14]、相对熵分析响应行为变化的趋势性 [15]、特定互信息表示降水对径流的预测程度 [16] 等。转移熵(TE)由Schreiber和Thomas [17] 提出，能够很好地分析两个变量的耦合强度和非对称驱动响应关系，但未应用于陆气耦合的研究。此外，TE仅以信息流大小衡量耦合强度，忽略了自信息量差异。

本文通过归一化的方式改进转移熵，统一量纲后使耦合强度具有可比性；基于长江流域实测气象站数据以及卫星反演的长序列土壤湿度数据，采用已有和改进的信息指标，探索从信息论的视角分析土壤湿度对降水的反馈特征及空间分布，为陆气耦合研究提供新途径。

2. 研究方法

先采用偏相关排除混杂因素影响，初步验证相互间的关系，再以信息论为基础的特征指标，分析得出变量关系的详细时空变化特征。

2.1. 偏相关分析

当变量y同时与变量x₁、x₂相关，可采用偏相关法剔除变量x₁的影响，只分析y和x₂的相关关系，被剔除的变量为控制变量。在包含3个自变量x₁、x₂、x₃和1个因变量y的系统中，偏相关系数的计算公式如下：

$r_{x_{1}y}=\frac{{\displaystyle {\sum }_{i=1}^n\left[\left(x_i-\bar{x}\right)\left(y_i-\bar{y}\right)\right]}}{\sqrt{{\displaystyle {\sum }_{i=1}^n{\left(x_i-\bar{x}\right)}^2}}\sqrt{{\displaystyle {\sum }_{i=1}^n{\left(y_i-\bar{y}\right)}^2}}}$ (1)

$r_{x_{1}y\cdot x_2}=\frac{r_{x_{1}y}-r_{x_1{x}_2}{r}_{x_{2}y}}{\sqrt{\left(1-r_{x_1{x}_2}^2\right)\left(1-r_{x_{2}y}^2\right)}}$ (2)

$r_{x_{1}y\cdot x_2{x}_3}=\frac{r_{x_{1}y\cdot x_2}-r_{x_1{x}_3\text{.}{x}_2}{r}_{y{x}_3\text{.}{x}_2}}{\sqrt{\left(1-r_{x_1{x}_3\text{.}{x}_2}^2\right)\left(1-r_{y{x}_3\text{.}{x}_2}^2\right)}}$ (3)

式中：x_i、y_i ( $i=1,2,3,\cdots ,n$ )是第i个数据中的x、y值， $\bar{x}$ 、 $\bar{y}$ 是 $x_1,x_2,\cdots ,x_n$ 和 $y_2,y_2,\cdots ,y_n$ 的均值；r表示变量间的相关关系，如 $r_{x_{1}y}$ 表示变量y和x₁的Pearson相关系数； $r_{x_{1}y\cdot x_2}$ 表示在控制变量x₂时，变量y和x₁的一阶偏相关系数； $r_{x_{1}y\cdot x_2{x}_3}$ 是控制变量为x₂、x₃时，变量y和x₁的二阶偏相关系数。

2.2. 转移熵方法

与模型模拟方法相比，信息论不受物理过程参数化和数据时空分辨率的影响；与现有统计方法相比，转移熵采用非参数估计的概率密度来描述，而不仅仅是数值的大小 [18]，可不带假设条件地衡量传统相关性指标等不能反映的非线性关系，从信息论的视角认识水文过程 [13] [14]。

2.2.1. 数据预处理

1) 归一化

每个时间序列变量x按下式将标准值映射至(0, 1)：

$x_{\text{Norm}}=\frac{x-x_{\min }}{x_{\max }-x_{\min }}$ (4)

式中： $x_{\max }$ 和 $x_{\min }$ 是样本数据的最大值和最小值，x和 $x_{\text{Norm}}$ 是标准值和映射结果值。

2) 零值效应

若变量x、y有大量零值，二维概率密度函数按下式计算：

$\{\begin{array}{l}p\left(x_i,y_i\right)=p\text{if}i=1,j=1\\ p\left(x_i,y_i\right)=p_{x0}{f}_{x0}\left(x_i\right)h_i\text{if}i>1,j=1\\ p\left(x_i,y_j\right)=p_{0y}{f}_{0y}\left(y_j\right)h_j\text{if}i=1,j>1\\ p\left(x_i,y_i\right)=p_{xy}{f}_{xy}\left(x_i,y_j\right)h_i{h}_j\text{if}i>1,j>1\end{array}$ (5)

式中： $f_{x0}\left(x_i\right)$ 和 $f_{0y}\left(y_j\right)$ 是边缘概率密度函数， $p_{xx}$ 是直接计数得到的各类数据点所占的比例。

3) 巴特沃斯高通滤波

对于变量x，巴特沃斯高通滤波可以保留短期波动以及忽略日或季节周期。

2.2.2. 信息熵和互信息

在信息论中，信息熵可定量描述不确定性，信息熵越大，表明不确定性越大。互信息是两个或多个变量共享的信息量，相关性越大，则互信息量越大。

$H\left(x\right)=-{\displaystyle \sum p\left(x_i\right){\log }_{2}p\left(x_i\right)}$ (6)

$H\left(x|y\right)=-{\displaystyle \sum {\displaystyle \sum p_{x|y}\left(x_i|y_i\right)p_{x|y}\left(x_i|y_i\right){\log }_2{p}_{x|y}\left(x_i|y_i\right)\text{d}x}}$ (7)

式中： $H\left(x\right)$ 是离散型随机变量x的信息熵， $H\left(x|y\right)$ 是随机变量x在得到变量y提供的信息后所剩余的不确定性，即条件互信息熵；x的可能取值x_i的概率为 $p\left(x_i\right)$，x_i 关于y_i 的条件分布函数为 $p_{x|y}\left(x_i|y_i\right)$。

互信息 $I\left(x;y\right)$ 为信息熵与条件互信息熵的差值。同理，偏互信息 $I\left(x;y|z\right)$ 定义为除去变量z的影响后，y对x提供的信息量。

$I\left(x;y\right)=H\left(x\right)-H\left(x|y\right)$ (8)

$I\left(x;y|z\right)=H\left(x|z\right)-H\left(x|y,z\right)$ (9)

信息熵计算的关键是估计联合概率密度，可用非参数方法求得，包括直方图法、核函数法和平均直方图法等。

2.2.3. 特征指标

在过程网络中，多个变量通过信息流相互关联。根据信息分离理论 [14]，U指仅由一个源变量单独提供目标变量的信息；R指多个源变量提供目标变量的重复信息源变量共同提供目标变量的信息。设每条信息流始于源变量x，终于目标变量y，本文用到以下特征指标：

1) 显著滞时

即x对y发生最大影响的时间 ${\tau }_0$，于 $I_{\max }\left(\tau ={\tau }_0\right)={\max }_{\tau =1}^{t}I\left(x\left(t-\tau \right);y\left(t\right)\right)$ 取得，其中 $\tau$ 是滞时( $\tau =0,1,\cdots ,t$ )。

2) 不同滞时的相对预测度

将总信息量 $I_{tot,\tau }$ 按式(10)分解，滞时1个单位时间跟 $\tau$ 个单位时间( $\tau =2,3,\cdots ,t$ )的x分别提供y的信息比例为 $\text{г}$ 和 $\eta$ ( $\text{г}+\eta =1$ )，即x的长短滞时历史对y现在状态的相对预测度(式(11)和式(12)) [19]。根据分离出的信息种类U、S、R，可指示链接冗余度， $\eta$ 值小则冗余度高。

$I_{tot,\tau }=I\left(x\left(t-1\right),x\left(t-\tau \right);y\left(t\right)\right)=I\left(x\left(t-1\right);y\left(t\right)\right)+I\left(x\left(t-\tau \right);y\left(t\right)|x\left(t-1\right)\right)$ (10)

$\frac{I\left(x\left(t-1\right);x\left(t\right)\right)}{I_{tot}}=U_{x\left(t-1\right)}+R_{x\left(t-1\right),x\left(t-\tau \right)}=\text{г}$ (11)

$\frac{I\left(x\left(t-\tau \right);y\left(t\right)|x\left(t-1\right)\right)}{I_{tot}}=U_{x\left(t-\tau \right)}+S_{x\left(t-1\right),x\left(t-\tau \right)}=\eta$ (12)

3) 归一化转移熵

当式(9)的z设置为x的滞时序列，得到的偏互信息称转移熵TE，可指示信息转移强度，但未考虑自信息差异。本文提出指标TE/H，以源变量的熵值为底进行归一化，使不同自信息量条件下的耦合强度可比。

$TE/H\left(x\left(t-\tau \right)\to y\right)=\frac{TE\left(x\left(t-\tau \right)\to y\right)}{H\left(x\left(t-\tau \right)\right)}=\frac{I\left(x\left(t-\tau \right);y|y\left(t-\tau \right)\right)}{-{\displaystyle \sum p\left(x\left(t-\tau \right)\right){\log }_{2}p\left(x\left(t-\tau \right)\right)}}$ (13)

以土壤湿度–降水耦合为例，降水时间序列为 $P\left(t\right)$，土壤湿度时间序列为 $SM\left(t\right)$ ( $t=\text{1},\text{2},\text{3},\cdots ,n$ )，设置滞时状态为 $t-\tau$ ( $\tau =0,1,\cdots ,t$ )。当 $I\left(P\left(t-\tau \right);SM\left(t\right)\right)$ 取最大值， ${\tau }_{01}$ 代表前期降水对土壤湿度发生显著影响的滞时，同理，当 $I\left(P\left(t\right);SM\left(t-\tau \right)\right)$ 取最大值， ${\tau }_{02}$ 代表前期土壤湿度对降水发生显著影响的滞时； $\text{г}$ 量化了仅通过了解滞时1个单位时间的土壤湿度 $SM\left(t-1\right)$ 的信息而减少 $P\left(t\right)$ 信息量的比例，包含 $SM\left(t-1\right)$ 单独提供的U信息，以及 $SM\left(t-1\right)$ 和 $SM\left(t-\tau \right)$ 提供的R信息。 $\eta$ 则表示 $SM\left(t-\tau \right)$ 提供的U信息以及 $SM\left(t-1\right)$ 和 $SM\left(t-\tau \right)$ 共同提供的S信息所减少 $P\left(t\right)$ 信息量的比例，该比例量化了较长滞时 $SM\left(t-\tau \right)$ 对 $P\left(t\right)$ 影响的相对重要程度。因此， $\eta =0$，代表 $P\left(t\right)$ 的不确定性通过了解 $SM\left(t-1\right)$ 获得了最大程度的降低， $SM\left(t-\tau \right)$ 提供的信息必定是冗余的。相反地， $\eta =1$，代表了解 $SM\left(t-1\right)$ 不能降低 $P\left(t\right)$ 任何不确定性，而要通过更长滞时 $SM\left(t-\tau \right)$ 的额外信息来降低；归一化转移熵指标 $TE\left(P\left(t-\tau \right)-SM\left(t\right)\right)/H\left(P\left(t-\tau \right)\right)$，即降水显著滞时下的历史序列 $P\left(t-\tau \right)$ 向当前土壤湿度状态 $SM\left(t\right)$ 转移的信息强度；而指标 $TE\left(SM\left(t-\tau \right)-P\left(t\right)\right)/H\left(SM\left(t-\tau \right)\right)$，即土壤湿度显著滞时下的历史序列 $SM\left(t-\tau \right)$ 向当前降水状态 $P\left(t\right)$ 转移的信息强度。

3. 试验验证

为验证指标适用性，采用蒙特卡罗试验模拟两变量x₁和x₂的互馈作用，随机产生1000个样本。模型1中，两变量以 $\tau =\text{3}$ 和4相互约束；模型2中，随机生成x₁再以 $\tau =\text{3}$ 约束x₂。设置 $f\left(X_i\right)=a{X}_i\left(1-X_i\right)$，令 $a=4$ 计算指标。此外，设定变化耦合系数a的9个情景，测试TE/H的强度指示作用。

模型1 $x_1\left(t\right)=f\left(x_2\left(t-4\right)\right),x_2\left(t\right)=f\left(x_1\left(t-3\right)\right)$ (14)

模型2 $x_1\left(t\right)=z~U\left(0,1\right),x_2\left(t\right)=f\left(x_1\left(t-3\right)\right)$ (15)

试验结果如图1，可得出以下结论：1) 改变变量的初始值，识别的显著滞时接近于理论值；2) $\eta$ 值可区分变量作用类型， $\eta =1$ 指示了x₁对x₂的单向作用，而 $\eta <1$ 指示了双向作用；3) TE、H曲线非单调，而TE/H随耦合系数a的增加而上升，验证其对耦合强度的指示作用。通过归一化，也解决了不同自信息量条件下的耦合强度相比时量纲不统一的问题。

4. 研究数据

4.1. 数据来源

长江流域(90˚33'~122˚25'E, 24°30'~35˚45'N)总面积180万km²，自西向东依次横跨19个省、市和自治区。流域流经青藏高原、四川盆地和长江中下游平原三大阶梯地貌，地形复杂，受西太平洋副高季风系统影响，降水量年内分配不均，空间分布差异明显，形成复杂的水文过程。

气象资料为长江流域2002~2018年168个气象站降水和蒸散发的逐日观测数据，取自中国气象数据网(http://data.cma.cn)的中国地面气候资料日值数据集(V3.0)，站点及水系分布如图2。

土壤湿度实测数据缺乏，微波遥感是实现土壤湿度大范围长时间连续观测的主要手段，不受大气影响且对地表具有一定穿透能力。先进微波扫描辐射计(AMSR-E)和第二代先进微波辐射成像仪(AMSR2)是搭乘在近年对地观测卫星EOS-PM1和GCOM-W1的传感器。Du、Kimball、Jones、Kim和Glassy [20] AMSR-E和AMSR2，进行了相似的校正和多频亮度温度(T_b)反演，得到新的全球陆地参数数据集(LPDR V2.0) (http://files.ntsg.umt.edu/data/LPDR_v2)，

图1. 蒙特卡罗试验结果：(a) 模型1跟2在不同源变量和目标变量下对应的τ、TE/H、η值；(b) 模型1在变化耦合系数下双向的H、TE、TE/H值

空间分辨率0.25^º。本研究调整坐标系为WGS-1984，并按气象站点位置，提取该数据集在点上的逐日土壤体积含水量。

为提供混杂因素，还使用了美国航天局NASA发布的MODIS NDVI (ID = MODIS/006/MOD13A2)数据，时间分辨率16日(https://ladsweb.modaps.eosdis.nasa.gov/)。数据在谷歌地球数据云平台(Google Earth Engine)中处理，利用数据集的质量控制字段(SummaryQA)，选取所有最优观测值(SummaryQA等于0)的最大值作为月NDVI。

4.2. 数据预处理

由于降水和土壤湿度在量级和变化幅度上存在较大差异，需要根据数据特征进行预处理，以还原它们在原始数据中所反映的相关性。本文采用以下预处理方式：1) 将所有数据按最值归一化至区间0~1；2) 连续的降雨数据中有大量缺测值，使概率密度偏移并影响信息熵计算，因此考虑了零值效应，以提升概率密度估算；3) 对单个变量，采用高通巴特沃斯滤波过滤日循环和季节循环，仅关注短期波动来挖掘降水和土壤湿度微小变化中的反馈机制 [13]。

图2. 长江流域水系及168个气象站分布

5. 结果与讨论

5.1. 基于偏相关的数据验证

常用的相关性方法可用来识别土壤湿度对降水的反馈，但两者的相互作用过程中常混杂其他因素，利用偏相关关系可排除干扰，更真实地反映变量间的相关性。

蒸散发和土壤湿度相互制约，且蒸散发与降水可能是正相关关系，也可能是负相关关系，使土壤湿度–降水机制具有较大不确定性；在植被覆盖的土地上，植被覆盖状况直接影响植被的蒸腾，从而影响土壤湿度的大小。因此选取降水、土壤湿度、蒸散发、NDVI四个变量，并将蒸散发及NDVI设置为控制变量。受NDVI数据尺度控制，所有变量数据均为2002~2018年的月尺度系列。

分前期降水–土壤湿度(P-SM)和土壤湿度–后期降水(SM-P)两组、设置0 d，5 d，10 d，15 d，20 d五种滞时开展研究，同时也检测控制变量与两个分析项的相关性，验证控制变量设置的合理性。结果见图3，可得以下结论：

1) 排除蒸散发及NDVI，P-SM在0~10 d滞时范围内存在显著(a = 0.01)偏相关关系，SM-P在1~8 d滞时范围也存在显著(a = 0.01)偏相关关系；

2) P-SM在0~20 d滞时范围内是正相关，SM-P在15~20 d滞时范围内是负相关；

3) 经检验，蒸散发及NDVI是强控制变量，P-SM和SM-P下均通过a = 0.01的验证。

综上，降水0~10 d内显著作用于土壤湿度，土壤湿度增多后又于0~8 d内显著作用于降水，均显示为正相关关系。即表明，排除混杂因素蒸散发及NDVI后，本研究所用数据序列存在土壤湿度–降水的正反馈。蒸散发及NDVI是强控制变量，可归因于存在于东亚季风区土壤湿度–降水的局地正反馈作用：较大的土壤湿度带来较大的蒸散发，输送水汽的同时，地表感热减少、潜热增加，波恩比变小使大气边界层变薄，增加了低层大气的湿熵，加大了对流不稳定性，水汽叠加边界条件改变导致降水发生。

图3. 长江流域土壤湿度对降水反馈的偏相关关系结果。红色空心折线表示前期降水–土壤湿度和土壤湿度–后期降水的偏相关关系，横坐标指代五种滞时情况；其他折线表示蒸散发–降水、蒸散发–土壤湿度、NDVI–降水、NDVI–土壤湿度的偏相关关系(见右侧图例)；深色和浅色虚线分别代表显著性水平a = 0.05和0.01下的相关性阈值

5.2. 基于转移熵的反馈特征

降水和土壤湿度的日时间序列 $P\left(t\right)$ 、 $SM\left(t\right)$ 可视为动态网络的节点，将序列按四季分为四组，设置滞时范围 $0,\text{1},\cdots ,\text{2}0\text{d}$，基于转移熵识别土壤湿度对降水反馈的时空特征。本文采用直方图法估计概率密度函数，代入可计算3个特征指标，分析如下。

图4是长江流域土壤湿度对降水的反馈空间特征。P-SM的显著滞时τ₀₁较小，空间差异性也较小，τ₀₁较长的区域集中于流域中部的汉江和乌江水系；SM-P的显著滞时τ₀₂较大，空间差异性也较大，τ₀₂较长的区域集中于流域西北部的金沙江、岷江、沱江和嘉陵江水系。整体上，P-SM的长滞时相对预测度η值大于SM-P，即降水滞时历史 $P\left(t-\tau \right)$ 对土壤湿度当前状态 $SM\left(t\right)$ 提供的信息较不冗余，提供的信息中，U信息与S信息比例居多，P-SM不冗余的区域位于金沙江下游和流域东部的洞庭湖、鄱阳湖、太湖水系，冗余的区域位于流域中部的嘉陵江水系、汉江水系、上游干流区间和乌江水系。而土壤湿度滞时历史 $SM\left(t-\tau \right)$ 对降水当前状态 $P\left(t\right)$ 提供的信息较冗余，即R信息比例居多，SM-P冗余的区域位于流域南部的金沙江下游和流域东部的洞庭湖、鄱阳湖、太湖水系，不冗余的区域位于流域北部的汉江水系、中下游干流区间。金沙江下游和流域东部的洞庭湖、鄱阳湖、太湖水系在P-SM和SM-P都较冗余，因此该区域内土壤湿度–降水反馈更受蒸散发、NDVI的混杂影响。

图4(a)包含指标 $TE\left(P\left(t-\tau \right)-SM\left(t\right)\right)/H\left(P\left(t-\tau \right)\right)$，图4(b)包含指标 $TE\left(SM\left(t-\tau \right)-P\left(t\right)\right)/H\left(SM\left(t-\tau \right)\right)$。整体上，P-SM的归一化转移熵TE/H明显高于SM-P，苏辉 [21] 以格兰杰因果关系和新型因果关系识别同深度下降水与土壤湿度的因果关系，得到相近结论。降水对土壤湿度具有直接的因果关系，因果影响较强，而土壤湿度对降水的因果混杂了间接关系，因果影响较弱，从因果机理上可验证信息熵识别耦合强度的适用性，可据此比较不同地区的耦合强弱。P-SM的信息转移强度在长江流域的干流区域较高，在其余支流区域上较低，推测是由于降水导致干流的河道水位上升，从而补充了周边土壤湿度，使P-SM信息转移强度较高。SM-P的信息转移强度与干支流关系不明显，强度更高的区域集中在流域东南部，尤其是洞庭湖、鄱阳湖地区，个别分布于流域西部的金沙江地区。

图5是长江流域土壤湿度对降水的反馈季节特征。SM-P的τ₀₂长度在秋季最长，τ₀₂长度较长的区域位于流域南部的乌江和洞庭湖水系；冬季的τ₀₂长度次之，集中区域向北部偏移至乌江水系和长江干流；春秋的τ₀₂长度相对短，分布也较为平均。整体上，春冬季的η值低于夏秋季；空间维度上，夏秋冬季的η值在中部的上游干流区间和嘉陵江水系较高，冬季的鄱阳湖和太湖具有较低的η值，而春季则在西部的金沙江下游及岷江水系较高，洞庭湖附近η较低。

图4. 长江流域土壤湿度–降水反馈空间特征。η以灰底色表示，颜色深度对应指标值；TE/H以橘色方块表示，方块边长对应指标值；显著滞时τ以黄色圆表示，圆半径对应指标值

图5. 长江流域土壤湿度–降水反馈季节特征。η以蓝底色表示，颜色深度对应指标值；显著滞时τ以黑色圆表示，圆半径对应指标值

图6和表1分别是168个站点土壤湿度对降水反馈显著滞时的频数分布图和统计值。可知，P-SM的区域平均显著滞时1.8 d，集中于1~2 d滞时范围；SM-P的区域平均显著滞时4.3 d，大部分站点显著滞时1 d，其余站点零散分布于2~20 d滞时范围。

图7和表2分别是168个站点土壤湿度对降水反馈归一化转移熵的频数分布图和统计值。可知，P-SM的区域平均归一化转移熵0.51，明显高于SM-P的0.13，且P-SM归一化转移熵的空间差异性略高。

表1. 长江流域土壤湿度对降水反馈显著滞时统计值

图6. 长江流域土壤湿度对降水反馈显著滞时频数分布图

表2. 长江流域土壤湿度对降水反馈归一化转移熵统计值

图7. 长江流域土壤湿度对降水反馈归一化转移熵频数分布图

综上，P-SM和SM-P均存在滞时。一方面，土壤湿度的记忆性使土壤湿度对降水补偿的响应存在滞后；另一方面，土壤湿度若要直接补充大气水汽，要以蒸散发为中介变量，若要通过间接反馈，则要改变成云致雨的大气条件，滞时较前者更长，并且受陆面状况影响很大，因此空间变异性较大。信息指标可反映序列间的时间依赖程度，据此分析得到P-SM和SM-P的滞时集中在1~2天。Lam、Bierkens和van den Hurk [22] 和李泽君 [1] 以广义线性模型得出前一天的土壤湿度会影响次日的降水概率，与信息指标的结果相近，且根据下渗率和土壤深度，距地表2~5 cm深度的土壤湿度应在较短的时间内对降水产生响应，因此从机理上看，信息指标的滞时识别更具合理性。显著滞时有明显的地区差异，因此可据此推断各地区陆气过程的滞时效应。

金沙江下游和洞庭湖、鄱阳湖、太湖水系等近湖区域土壤湿度对降水反馈较快且强度较大。基于陆气耦合模式LIS-WRF、CAM-CLM、GLACE等试验表明，位于大陆内部，远离海洋的干湿过渡地带是土壤湿度对降水反馈的热点地区 [7]，如北美大平原、美国中西部地区 [23]，同样为近湖区域，位于美国五大湖西部及南部。邹海波 [24] 指出鄱阳湖跟美国五大湖一样具有湖效应降水，同理，洞庭湖水域面积更大，识别出较大的反馈强度。近湖区域土壤湿度对降水可能的反馈机制可解释为：湖陆温差激发湖陆风环流，低层气流向湖中心辐合，抬升湖体上空的暖湿空气，进而触发对流形成湖效应降水，水量交换过程频繁，可能导致较低的η值，较高的г值和TE/H值，即近湖区域土壤湿度在短时间反馈于降雨且反馈强度较高。

6. 结语

本文基于长江流域2002~2018年的气象观测及陆面遥感反演数据，以显著滞时、不同滞时的相对预测度和归一化信息熵3个信息指标识别土壤湿度对降水的反馈特征，并以偏相关方法验证了数据系列中土壤湿度对降水的反馈，主要结论如下：

1) 长江流域前期降水对土壤湿度影响的平均显著滞时1.8 d，空间差异性小、归一化转移熵0.51，影响强度较大；土壤湿度对后期降水影响的平均显著滞时4.3 d，空间差异性大、归一化转移熵0.13，影响强度较小；

2) 金沙江下游和洞庭湖、鄱阳湖、太湖水系等近湖区域，由于水量交换频繁且存在湖效应降水，识别出土壤湿度对降水反馈较快且强度较大；

3) 排除蒸散发及NDVI后，偏相关分析得到降水0~10 d内作用于土壤湿度，土壤湿度增多后又在0~8 d内作用于降水，均显示为显著的正相关关系，可归因于土壤湿度对降水的局地正反馈，相比于信息论方法，偏相关识别出的滞时长度不合理。

应当指出，受限于高维概率密度函数估计，基于转移熵的因果推断未能排除多维变量影响；同时，今后可开展陆气耦合模式不同集合下的数值模拟试验，以期对前述的反馈特征予以验证。

基金项目

国家自然科学基金项目(51861125102)和清江水电开发有限责任公司(ENQJ-TD1-FW-2019185)资助。

文章引用

娄 巍,刘 攀,程 磊,邹铠杰,夏 倩. 基于转移熵的长江流域土壤湿度对降水反馈研究
Soil Moisture Feedbacks on Precipitation in the Yangtze River Catchment Based on Transfer Entropy[J]. 水资源研究, 2021, 10(01): 21-32. https://doi.org/10.12677/JWRR.2021.101003

参考文献