在一定的气候和干湿条件下, 以水为介质的物质交换、能量流动或有机物运输或运动过程, 会受到地理单元组成及结构格局之间水文连续性的影响^[1], 其中这个过程为水文功能连通, 而影响这个过程的水文连续性则为水文结构连通^[2-3]。水文结构连通反应的是静态格局, 其在空间上表现出三维特征：沿河道或潮沟方向即纵向, 水平正交于河道或潮沟方向即横向, 和沿重力方向正交于地表切线方向即垂向。

水文连通最早是针对河流系统或流域提出的, 旨在分析研究区域的水文结构连通与产流等其他水文过程之间的关系, 因此相关研究多着眼于基于地形特征^[4-5], 结合景观生态学和图论理论开发结构连通指数并探究其对水文过程或生态过程的影响^[6-8]。滨海湿地水文连通的研究在初期借鉴了该思路：基于潮沟形态和分布量化纵向水文结构连通并寻求其与生态因子和生态过程之间的关系^[9-10]。然而, 纵向水文连通的计算方法的基础假设认为潮沟或河道分布是产生水文连通的结构基础。但是对于滨海湿地而言, 海水通过地表或浅层地下潮流迅速向陆相滨海湿地汇集；而在低潮时, 滨海湿地包气带和潜水层储存的水会向海相渗出^[11], 因此潮间带及潮上带部分生境斑块虽然在地理上与潮沟隔离, 但仍然存在横向水文连通, 其在静态结构上表现为土壤含水量的空间差异^[2]。由于地形和土壤结构的高度空间异质性, 即使在较小的空间尺度内, 滨海湿地横向水文功能连通仍然存在着高度的时空差异^[12-13]。

目前已有研究多基于水文模型量化滨海湿地横向水文连通, 但这一方法仅适用于较大的区域或流域尺度, 当其应用于生态系统尺度时, 会由于数据在空间上的一致化假设而忽略细节信息。因此, 对横向水文功能连通研究的基础和关键是运用现场观测的方法来对不同截面进行数据采集。近两年有学者基于黄河三角洲的野外表层土壤含水量观测数据, 结合微地形以及图论理论在样地(1 m²)和群落(100 m²)两个尺度量化了不同类型滨海湿地的横向水文结构连通并分析其对植物群落的影响^[14-16]。黄河三角洲植物根系分布较浅, 表层土壤与植物之间的关系密切^[17], 因此该方法建立的水文结构连通指数不仅能反映横向水文连通的空间异质性, 对湿地植被管理和恢复更有实用价值。目前的问题在于其如何推演到生态系统尺度(>100 hm²)上。

由于采样点分布不均导致的空间分辨率不一致, 在黄河三角洲所建立的基于野外数据监测的横向水文连通可能会存在尺度效应。已有研究表明, 滨海湿地潮间带的水流在40 m和100 m的观测尺度上就会出现显著差异^[18], 但是各个尺度之间变异性及其来源目前并不明确^[19-20]。此外, 水文过程与相关因子之间的统计学规律会随着尺度的改变而发生非线性变化^[21], 并显示出更加复杂化和多样性的相关关系^[22]：在某一尺度得出的研究结论很难以线性关系推演到另一尺度中。因此, 将样地(小)尺度方法应用于生态系统(大)尺度上时, 需要借助无标度网络模型进行重推演。此外, 寻找各个尺度之间的差异性来源也将更好的帮助我们揭示滨海湿地水文连通特征。

综上, 本文提出以下研究目标：1)改进目前样地和群落尺度的滨海湿地横向静态水文连通的参数化方法, 将其应用到较大的生态系统尺度并用以分析黄河三角洲滨海湿地的水文连通, 2)结合已有研究结果分析黄河三角洲滨海湿地横向静态水文连通在不同尺度的变异性及其变异性来源。

1 试验地概况

黄河三角洲位于山东省东营市河口区, 山东黄河三角洲国家级自然保护区内(118°32.98′—119°20.45′ E, 37°45.77′—38°12.31′ N)^[23], 本文的研究区域如图 1所示, 位于国家级自然保护区内黄河以北到渤海一侧的滨海湿地。

2 研究方法 2.1 实验设计

根据黄河三角洲的湿地类型^[24]和分布, 将研究区域中潮间带、潮上带、河滨带以及常年采取人工补水管理的芦苇湿地分别定义为潮间区(Inter-tidal Flat, IF)、盐沼区(Tidal Marsh, TM)、河滨区(Artificial Pound, AP)和生态补水区(River Side, RS), 四个区域面积分别为10.81 km², 14.5 km², 8.10 km²和1.70 km²。结合已有研究, 本文究将1 m², 100 m²和10⁶ m²定义为小尺度、中尺度和大尺度并如图 1所示, 在四个区域中分别设置12、19、12和11个采样点进行环刀和铝盒采样。样品采集和保存按照林业行业标准(LY/T 1215—1999)及农业行业标准(NY/T 1121.1—2006)进行, 野外采样时间为2017年和2018年6月下旬至8月上旬；为避免潮汐影响, 每日取样时间为3:30 —12:00 am。

土壤环刀带回实验室后立刻对田间持水率进行测定, 其可表示为^[25]：

(1)

式中, w为土壤样品取样质量, w_s为样品饱和质量, w_d为饱和样品在干燥沙土之上48 h后的质量, w_n为土壤样品在105℃烘干至恒重质量, w₀为环刀质量。

2.2 大尺度横向水文结构连通参数化方法

应用无标度网络模型, 基于野外监测数据建立大尺度横向水文结构连通参数化方法。网络中某一结点V_i上边的数量为该结点的中心度(Degree), 其满足幂律分布, 即：

(2)

或：

(3)

对于离散点而言, 该分布满足

(4)

式中, C和γ为两个参数, k为中心度, 因此：

(5)

式中, ζ(γ)为西塔函数^[26], 用于解析无穷级数的变复元。区域中每一个样地为一个结点, 依据结点土壤水分含量计算该结点的连通替换值：

(6)

式中, w_ni和z_ni分别为结点V_i的土壤含水率和田间持水率, I_ni为连通替换值。I_n>0的结点为连通结点, 结点V_i的属性为其边的数量(k_i), 则大尺度水文连通指数(Landscape Hydrological Connectivity Index, LHCI)数值为指定置信水平(本文采用0.8和0.9)上无标度网络结构的参数估计值：

(7)

式中, α为置信水平。

2.3 水文连通变异性分析

多重分形分析常用于不均匀分布的定量表征^[27], 其中广义分形维数可以表示为：

(8)

(9)

式中, δ为尺度系数, n为尺度δ中的样本数, u_i(δ)为质量概率, 其可以表示为：

(10)

当D(0)、D(1)和D(2)分别为容量维数、信息维数和关联维数, 信息维数和关联维数相对容量维数减少越多, 远距离的空间异质性越高；反之则异质性越低^[28-29]。

以q为变量, 可以得到奇异指数α(q)和维度分布函数f(q), 其可以分别表示为：

(11)

(12)

(13)

式中, 维度分布函数和奇异函数的拟合曲线被称为多重分型谱系, 该拟合曲线为开口向下的二次函数图像或者图像的一部分；当曲线拟合点集中在图像峰值右侧时, 说明变异性来自于较大的值, 当曲线拟合点集中在图像峰值左侧时, 说明变异性来自于较小的值^[30], 此外, 奇异函数最大值和最小值之差(Dδ)为变异系数, 值越大说明变异性越强^{[28, 31]}。本研究的分形谱计算取[-20, 20]作为q的取值区间。

3 结果与分析 3.1 横向水文连通分析

研究区域大尺度横向水文结构连通网络如图 2所示, 潮间区12个样地皆为连通结点, 平均中心度为1.83；其中V₄₈和V₄₉的中心度最高, 分别为7和6；依据网络结构该区域可分为A、B和C3个子区域, 各子区域内部连通平均中心度分别为1.40, 1.25和1.00；此外子区域A和C均向子区域B连通。盐沼区19个样地中有14个连通结点, 平均中心度为0.89；其中V₁₃和V₁₇结点的中心度最高, 分别为5和4；依据网络结构该区域可分为4个子区域；其中子区域A和子区域C向其他区域连通。生态补水区12个样地皆为连通结点, 平均中心度为1.50；其中V₃₈和V₃₉结点的中心度最高, 皆为5；依据网络结构该区域可分为3个子区域, 其中子区域C向子区域A连通, 子区域C和B则为双向连通。河滨区12个样地中连通结点为6个, 平均中心度为0.42；依据网络结构该区域可分为2个子区域, 仅有一条连通边。各个区域中, 未连通结点所分布的生境是水文学上较为脆弱的区域, 建议持续关注。

各区域大尺度横向水文结构连通参数估计的概率分布如图 3所示, 中心度参数估计如表 1所示。盐沼区和河滨区的两个参数C较为接近且数值较小, 生态补水区和潮间区的参数C较高, 说明潮间区和生态补水区分布概率随度数增长相对较为缓慢, 即结点中心度数值分布更加平均；而潮间区γ系数较高, 说明潮间区分布函数平缓区值较高, 即结点中心度出现大数值的概率较高。

当置信水平为0.8时, 潮间区、盐沼区、生态补水区和河滨区的结点中心度估计值分别为2.49、0.97、1.81和0.97；当置信水平为0.9时, 4个区域中心度估计值分别为2.27、0.83、1.54和0.83, 与各个区域实测的平均结点中心度接近但略高于实测值。综上, 基于实测值和参数估计值的结中心度分布, 四个区域大尺度水文连通由强到弱排序为潮间区>生态补水区>盐沼区≥河滨区。

本文建立的大尺度横向水文结构连通的方法以区域为研究对象, 引入无标度网络从而消除了因采样点距离不一致而可能产生的对结果的干扰；所得到的结果主要反映了各区域内部水文过程及其趋势。综合大尺度水文连通指数和结中心度分布可知, 潮间区在日潮影响下, 水文联系频繁, 因此连通性最强。盐沼区, 仅有向海一侧及紧邻潮沟的区域水文连通较强, 其他区域表层土壤受潮汐影响较小, 导致水文连通较弱。同理, 河滨区仅有靠近黄河一侧因受河水脉冲作用较为明显而水文连通较强。前人对水分获得性研究中得到的水分获得性较强的区域和本文水文连通性较强的区域吻合^[29], 这进一步说明潮汐和河流作用能够增强横向水文结构连通。

3.2 不同尺度横向水文结构连通变异性分析

已有研究结果^[14-16]显示, 小尺度水文连通指数在潮间区、盐沼区、生态补水区和河滨区分别为0.49, 1.22, 0.61和1.16。中尺度水文连通指数在4个区域分别为0.87、0.55、1.53和0.48, 其中生态补水区显著高于其它3个区(P < 0.01), 潮间区显著高于盐沼区(P=0.04)和河滨区(P=0.02)。

表 2给出了空间变异性分析结果, 由表可知, 潮间区水文连通指数在小尺度、中尺度和大尺度上, D2相对D0分别降低了0.09%, 2.21%和2.33%, 该区域小尺度的水文连通长距离的空间变异性较小, 而中尺度和大尺度水文连通的长距离空间变异性则相对较高；该区域, 变异系数在小尺度、中尺度和大尺度依次增高, 在大尺度上该区域的水文连通变异性最强。盐沼区3个尺度上D2相对D0分别降低了5.00%, 9.85%和6.55%, 说明该区域中尺度的水文连通长距离的空间异质性明显, 而在另外两个尺度上则相对较弱；该区域大尺度的变异系数大于其余两个尺度, 该区域的水文连通在大尺度上的变异最为明显。河滨区3个尺度上D2相对D0分别降低了13.99%, 26.24%和6.56%, 该区域中尺度的水文连通长距离的空间异质性较为明显, 而大尺度空间异质性较弱；该区域, 小尺度、中尺度和大尺度的变异系数依次增大, 说明在大尺度上该区域的水文连通变异性最强。生态补水区在小尺度上的D0、D1和D2的值相等, 中尺度和大尺度上D2相对D0分别降低了3.44%和1.29%, 说明该区域3个尺度上水文连通长距离的异质性均相对较小；该区域中尺度的变异系数最小, 大尺度的变异系数则比其他两个尺度高出一个数量级。

对比不同区域可知, 潮间区和生态补水区在3个尺度的上, D2相对D0的降幅都较小, 而且在3个尺度上生态补水区和潮间区的变异系数也相对较低, 说明盐沼区和河滨区的水文连通的空间变异性更加显著。

3.3 变异性来源分析

基于多重分形谱的黄河三角洲滨海湿地水文连通尺度效应分析如图 4所示。对于潮间区而言, 水文连通指数的变异性在小尺度主要来自于较小的水文连通指数值, 而在中尺度和大尺度, 倾向于受到较大指数值的影响, 这是由于该区域小尺度水文连通指数很少出现较大的数值, 而在中尺度和大尺度水文连通值分布较为平均。对于盐沼区而言, 3个尺度上的水文连通变异性几乎都来自于较小的值, 这是由于盐沼区水文连通弱, 3个尺度的指数值皆集中在较小的范围内, 仅有中尺度和大尺度临近潮沟区域偶尔出现较大数值。对于河滨区而言, 河滨区仅在临近黄河河道一侧, 大尺度水文连通指数值较高, 而其余两个尺度值都偏低, 因此, 在小尺度和中尺度水文连通的变异性几乎都来自于较小的指数值, 而在大尺度, 较大的指数对空间异质性有一定的影响却不明显。对生态补水区而言, 3个尺度水文连通指数的值都较高, 较大的值之间差异也相对更加明显, 因此变异性主要来自于较大的水文连通值。此外, 小尺度和中尺度方法描述了一个样地与周边较小范围的空间区域土壤水分连通状况, 在大的潮汐或补水条件下不同区域各个样地倾向于保持一致；而大尺度方法则侧重于描述区域内的水文过程, 这是导致尺度差异的根本原因。

对比不同区域可知, 盐沼区和河滨区水文连通的变异来源主要为数值较小的水文连通指数, 生态补水区则主要为较大的指数值；而潮间区的变异来源则比较平均。盐沼区和河滨区的水文连通较弱, 而生态补水区和潮间区水文连通较强, 因此各个区域中, 水文连通较强的区域变异来自于连通指数较高的值, 而水文连通较弱的区域则与之相反。

4 讨论

不同尺度之间的差异来源主要有两个：(1)水文结构连通的参数化方法和(2)水文过程差异。从水文过程的角度来看, 潮间区受到日潮影响较为明显, 而盐沼区的水淹时间较短, 水淹深度较浅；河滨区主要受河水脉冲影响, 土壤含水量与河道的距离负相关；生态补水区在人为管理下, 维持着相对稳定的干湿条件。

已有研究已经在小尺度和中尺度建立了较为成熟的, 适用于黄河三角洲滨海湿地的横向水文结构连通参数化方法^{[4, 14-16]}：首先在与本文采样点相同的位置设置10 m×10 m的中尺度样地；然后将每个中尺度样地分为100个1 m×1 m的单元, 即小尺度样地, 记录每个单元土壤含水量及整个样地的田间持水率。在样地中, 某一单元C_ij的土壤质量含水率(w_ij)和田间持水率(z_ij), 然后计算每个小尺度样地的连通替换值：

(2)

与大尺度计算方法相同, 连通替换值大于0的单元为连通节点(V_ij)^[32-34], V_ij在水平投影面上以8-邻域模式与相邻的连通节点连接。从节点V_ij到其某一相邻连通节点V_rc的点连通指数分量：

(14)

而节点V_ij的小尺度连通指数可以表达为：

(15)

中尺度水文连通指数(TOFCI)可以表示为

(16)

式中, EHD为两个相邻连通节点的欧式-水文距离, 其可以表达为：

(17)

式中A为样地面积(m²), m为样地中FL的数量, FL_i(m)为第i个FL的值, 其可以表达为：

(18)

从方法上看, 小尺度的参数化方法仅考虑了小范围内, 一个样地及与其相邻的8个1 m²样地间土壤含水量差异, 描述了空间区域表层土壤通过水分运动而相互连通的状况, 再通过对不同区域所有小样点指数的统计学分析比较各区域水文连通差异。这一方法默认小区域内地形的同质性, 其仅对土壤含水量和田间持水量的变化十分敏感。稳定的潮水活动和规律的人为补水使得潮间区和生态补水区的土壤含水量普遍高于田间持水率, 因此该方法统计学分析所得的结果必然是这两个区域水文连通较高；同理, 盐沼区整体土壤含水量较低, 因此小尺度水文连通较弱。

在中尺度上, TOFCI指数在10 m×10 m样地综合考虑了土壤含水量分布以及地形因素, 结合图论理论, 量化了目标区域水文连通的方向及大小；之后同样根据多个样地计算数值的统计学分析评估不同区域的水文连通强弱, 其优点在于不仅能够给出各个研究样地的水文连通值, 还可以直观的表达样地内部各个单元之间的连通方向^[33]。在黄河三角洲, 生态补水区除土壤水分含量较高以外, 靠近水面一侧和远离水面一侧土壤含水量异质性明显, 加之该区域地形较为平坦, 三个因素叠加使得TOFCI指数值都相对较高且方差较小；在潮间区, 虽然潮汐活动保证了较高的土壤含水量, 但在一个样地中土壤含水量缺乏空间异质性, 加之潮沟的存在使得地形异质性较高, 因此该区域TOFCI指数相对生态补水区较低。对盐沼区和河滨区土壤水分含量较低, 大部分单元为非连通单元, 因此其整体中尺度TOFCI指数较小, 水文连通较弱。

在大尺度上, 区域大小、样点分布、采样空间分辨率等多个因子会导致不同区域研究结果缺乏可比性, 这也是目前大尺度横向水文连通多借助水文模型来进行研究的原因。本文引入了无标度网络模型消除上述因子对研究结果的影响, 结合样地的分布和土壤含水量构建了水文结构连通网络, 最后通过对网络结构的参数估计便可量化各个区域的水文连通并对其进行比较。这一方法直接以区域为研究对象, 通过基于样点分布和连通边的权重的空间聚类, 能够给出连通子集以及需要重点关注的连通较弱的区域。

黄河三角洲植物以草本为主, 根系分布较浅, 表层土壤的水文连通与植物个体生长、群落分布乃至植被景观结构的关系十分密切^[33-34], 本文在3个尺度上构建的水文连通算法也着眼于表层土壤, 因此能够应用于后续对水文连通与植物群落关系的研究。此外, 研究区域中河滨区的横向水文结构连通受季节和天气影响较为明显, 而潮间区和盐沼区则主要受潮汐驱动, 是否会有明显的季节异质性需要更多的监测数据。

5 结论

黄河三角洲滨海湿地不同区域横向水文结构连通强弱在3个尺度上存在细微的差异, 一个较为统一的结论是：生态补水区和潮间区的水文连通强于而盐沼区和河滨区。此外, 在水文连通较强的区域, 其远距离的空间变异性主要来自于水文连通指数较高的样点, 反之则来自于指数较小的样点。

不同尺度的差异主要来源于水文过程本身的尺度效应及各个尺度上水文连通量化方法的区别。从水文过程的角度来看, 潮间区受到日潮影响较为明显, 而盐沼区的水淹时间较短, 水淹深度较浅；河滨区主要受河水脉冲影响, 土壤含水量与河道的距离负相关；生态补水区在人为管理下, 维持着相对稳定的干湿条件。从方法上看, 小尺度的参数化方法仅考虑了小范围内, 空间区域表层土壤通过水分运动而相互连通的状况, 其对各个采样点的土壤含水量和田间持水率异质性十分敏感。中尺度水文连通参数化方法综合考虑了土壤含水量分布以及地形因素, 量化了目标区域水文连通的方向及大小；之后同样根据多个样地计算数值的统计学分析评估不同区域的水文连通强弱。大尺度水文连通参数化方法基于直接以区域为研究对象, 参数化的结果直接受到水文过程的影响而在一定程度上忽略了土壤含水量和田间持水率的空间异质性。