入渗性能是土壤的主要物理性质, 对土壤调节径流、涵养水源能力有重要作用, 也是评价土壤侵蚀的重要指标, 研究土壤入渗能力对防治水土流失、改造及恢复植被有重要意义^[1]。研究发现土壤容重、孔隙状况、颗粒组成、水稳性团聚体及有机质含量等都与土壤入渗能力密切相关^[2-6]。土壤理化性质是导致不同研究区域土壤渗透性能存在差异的主要因素之一, 但不同研究区主导因子不同, 如在青藏高原地区^[5]容重、有机质、< 0.1 mm的微粒是影响土壤入渗的主要因素, 而在岩溶坡地^[2]则为毛管孔隙度和黏粒含量。同样, 不同土地利用类型的土壤入渗性能差异也与其土壤理化性质的差异具有很大相关性。一般来说林地丰富的根系活动导致其孔隙状况优于草地、容重显著低于草地从而具有良好的土壤渗透性^{[2, 7-8]}, 但风沙地^[3]由于土质疏松、结构松散, 荒草地入渗速率较大, 林地对土壤的改良作用强于草地, 使得土壤质地变好、非毛管孔隙度减少、蓄水能力增强, 入渗速率较草地小, 丁康等^[9]发现长武塬边坡地的侧柏林地由于土壤黏粒含量低不利于形成团聚体且土壤密度较草地大, 入渗性能较差；灌木林地地表植被较林地丰富, 土质疏松, 入渗性能优于林地^[10]；与阔叶林相比, 针叶林有机质含量高、孔隙较多, 入渗能力较强^[11]；梯田旱地土壤黏粒含量较高, 有机质和孔隙度较低, 入渗能力较灌木林地差^[12]。此外学者对退耕年限^[13-14]及演替进程^[15]影响下土壤渗透性能的研究发现：随种植年限的增加土壤理化性质得到改善, 入渗性能逐渐提高；随演替进展土壤饱和导水率显著提高, 植被演替到顶级群落土壤饱和导水率最高。可见不同研究区甚至同一研究区不同研究对象的土壤入渗特征及影响因素存在差异, 但以上研究主要集中在比较不同植被类型、退耕年限、土壤类型的土壤渗透特征, 对同一植被不同坡位下的研究缺乏。不同地形部位土壤微气候的变化会导致土壤条件产生差异, 从而影响植被的恢复进程与生长状况^[16]。坡位是重要的地形因子, 对林地有机质、氮等养分及水热状况的分布有着重要影响^[17], 植被的多样性、生物量等存在差异, 从而导致不同坡位土壤容重、机械组成等性质发生变化, 进而影响土壤的入渗特性。黄土高寒区位于青藏高原与黄土高原交界地带, 生态环境脆弱, 是植被重建与生态退耕的重点区域。该地区太阳辐射强烈, 不同坡位水热条件差异明显, 导致植被生长也存在明显差异, 土壤水分是限制该区域植物生长和恢复的主要因子。为减轻水土流失、保护和改善生态环境, 该地区2000年开始大规模实行退耕还林工程, 水土流失得到极大改善, 但在实施过程中不考虑土壤水分条件, 在坡面由上至下多配置相同植被, 可能会导致植被生长不良, 无法充分发挥对环境的改善作用。研究典型植被类型在不同坡位下的土壤入渗特征及其主要影响因子可以更好地了解当前植被分布的合理性, 为验证典型植被拦蓄降雨、涵养水分能力提供数据支撑, 为后续植被恢复与调整做出指导。

1 材料与方法 1.1 研究区概况

试验地位于青海省大通县桥头镇安门滩小流域, 地处黄土高原西部向青藏高原的过渡带, 海拔高程2448—2562 m。属高原大陆性气候, 降水少蒸发大, 多年平均降水量506.5 mm, 且主要集中在6—9月份, 多年平均蒸发量1762.8 mm, 最大蒸发量出现在4—6月, 属于典型的黄土高寒区小流域。主要树种为2000年退耕后种植的青海云杉(Picea crassifolia)、华北落叶松(Larix principis-rupprechii)、祁连圆柏(Sabina przewalskii)、中国沙棘(Hippophae rhamnoides)、柠条(Caragana korshinskii)等, 主要草本植物有密生苔草(Carex crebra)、垂穗披碱草(Elymus nutans)、鹅绒委陵菜(Potentilla anserna)等。土壤类型为黄土母质上发育的山地棕褐土和栗钙土。

1.2 研究方法 1.2.1 样地布设与土样采集

研究区位于青海省东部安门滩小流域, 2000年进行退耕还林, 行间距为3 m×2 m。2018年7月至8月, 在流域内选取青海云杉林、华北落叶松林、祁连圆柏林3种典型植被类型设置样地进行调查, 并选择荒草地作为对照(表 1)。按照样地海拔高度进行坡位分类, 样地实际海拔与坡顶海拔之比为小于1/3的定义为下坡位, 1/3—2/3的为中坡位, 大于2/3则为下坡位。针对每一种典型植被类型在不同坡位分别布设2个10 m×10 m的样地进行每木检尺, 确保同一坡位各样地保持在同一海拔高度上。

每个样地沿任一对角线的三等分点处挖取2个60 cm深的土壤剖面, 用环刀在垂直方向上每隔20 cm分层采集原状土样, 每层3个重复, 共360个样品。同时, 每个土壤剖面分层用“四分法”取扰动土500 g左右装入布袋带回室内自然风干, 共120袋。

1.2.2 土壤理化性质测定

室内用定水头法测定^[8]土壤渗透速率；环刀法和烘干法测定土壤容重、土壤含水量、总孔隙度、毛管孔隙度和非毛管孔隙度；利用Mater Sizer 2000激光颗粒分析仪测定土壤机械组成, 并按国际制分为黏粒(< 0.002 mm)、粉粒(0.002—0.02 mm)、砂粒(0.02—2 mm)；湿筛法测定＞3 mm、3—2 mm、2—1 mm、1—0.5 mm、0.5—0.25 mm、0.25—0.1 mm、< 0.1 mm各粒级的水稳性团聚体含量；重铬酸钾稀释加热法测定有机质含量。

1.3 数据分析

选取Kostiakov方程、Horton方程、通用经验方程3种常用的模型, 在SPSS 18.0中对不同植被类型的土壤入渗过程进行拟合, 模型表达及参数意义见前人研究^{[4, 12, 18]}。在SPSS中进行ANOVA方差分析、Duncan多重比较、Person相关分析, 并利用多元回归、通径分析提取主导因子。在Origin 8.5中作图。

2 结果与讨论 2.1 不同植被类型土壤入渗特征

由图 1可知, 青海云杉林地、华北落叶松林地、祁连圆柏林地、荒草地60 cm土层内稳渗速率均值分别为3.59、2.04、3.01、1.45 mm/min, 存在显著差异性, 林地均显著高于荒草地(P < 0.01)。说明相较荒草地, 林地具有显著改善土壤结构, 增加入渗的能力。通过比较3种与荒草地的基本理化性质发现(表 2), 林地容重小于荒草地, 孔隙状况优于荒草地, 这是因为林地发达的根系、丰富的地表植被和枯枝落叶增加了土壤孔隙, 提高了有机质含量, 促进了团聚体的形成, 对土壤的改善作用较荒草地显著, 从而提高了土壤的渗透性能^[18]。这与姚淑霞等^[19]研究结果相同, 而徐学选等^[9]在黄土高原研究发现退耕年限较短草地表层的饱和导水率大于林地。可能是由于研究区内人类活动剧烈, 频繁放牧, 导致荒草地被践踏、结构变紧密, 入渗能力降低。

流域内不同坡位光照水分的再分配不均, 制约植被的生长, 也影响土壤理化性质的变化。由图 2可知, 同一坡位3种林地渗透速率不存在显著差异性(P＞0.05), 但仍存在上坡位和下坡位青海云杉>祁连圆柏>华北落叶松, 中坡位青海云杉＞华北落叶松＞祁连圆柏。这是由于上坡位、下坡位林地总孔隙度、毛管孔隙度表现为：青海云杉＞祁连圆柏＞华北落叶松, 容重表现为：华北落叶松＞青海云杉＞祁连圆柏；中坡位3种林地的土壤孔隙度不存在显著差异(P＞0.05), 但华北落叶松林地总孔隙度、毛管孔隙度较小, 容重显著大于青海云杉林地、祁连圆柏林地, 这是由于华北落叶松为落叶树种, 林下枯落物丰富, 但地表植被缺乏, 根系分布较少, 导致土壤中孔隙较少, 导水性能较差, 加之林地内人类游憩活动频繁, 践踏导致土壤容重较大。

3种林地在不同坡位间的土壤渗透速率不存在显著差异(P>0.05), 但仍表现出青海云杉林地、祁连圆柏林地初渗速率和稳渗速率及华北落叶松林地稳渗速率上坡位＜中坡位＜下坡位。这是由于研究区下坡位土壤水分条件较好, 植被生长旺盛, 有机质与水稳性团聚体质量分数较高, 土壤结构良好, 入渗能力较强(表 2)。吴钦孝等^[20]在黄土丘陵区发现由于降雨侵蚀产生的泥沙细颗粒沿坡面逐渐沉积, 堵塞土壤孔隙, 导致森林及荒坡草灌小流域沿坡面向下土壤的渗透速率减小；纳磊等^[21-22]在晋西黄土区研究发现由于不同坡位土壤水分条件的差异, 导致下坡位植被生长较上坡位茂盛, 同一植被类型均表现为下坡位土壤入渗速率较上坡位高。后者与本研究的结果一致, 可能是由于研究区大规模进行退耕还林, 人工整地使坡度减缓, 泥沙颗粒沿坡面的沉积较弱, 且研究区海拔较高、辐射强烈使得下坡位水分条件较中、上坡位好, 植被生长茂盛, 根系穿插导致土壤孔隙状况良好, 所以下坡位入渗速率较大。

2.2 不同植被类型土壤入渗过程模拟

图 3表明土壤入渗速率随时间推移逐渐减小, 最终趋于稳定。除荒草地外, 青海云杉林地、华北落叶松林地、祁连圆柏林地在入渗前期入渗速率有明显的降低过程。除青海云杉林地中坡位和华北落叶松林地上坡位的土壤入渗速率在40 min左右达到稳定外, 其他土壤的入渗速率均在15 min左右即达到稳定。各林地下坡位土壤的稳渗速率明显高于上坡位和中坡位。将不同林地上、中、下坡位的土壤入渗速率(y)随时间(t)的变化过程选取不同模型进行拟合, 结果如表 3所示：青海云杉林地中Horton模型拟合R²最大介于0.862—0.918, 通用经验模型(R²介于0.614—0.839)拟合效果最差, 3种模型对入渗性能最好的下坡位入渗过程的模拟拟合程度最好。华北落叶松林地中, Horton模型(0.827—0.910)拟合程度最好, 与Kostiakov模型(0.853—0.899)拟程度的差异不大, 通用经验模型(0.702—0.785)拟合效果最差。对于祁连圆柏林地, 3种模型拟合程度为Horton模型(0.895—0.930)＞Kostiakov模型(0.883—0.984)＞通用经验模型(0.756—0.982)。对荒草地入渗过程模拟效果是Kostiakov模型(0.938)＞Horton模型(0.938)＞通用经验模型(0.830)。

综上所述, Horton模型对不同植被类型的土壤入渗过程模拟的拟合精度均在0.8以上；Kostiakov模型除青海云杉林地上坡位拟合精度稍差(0.727), 其他部位均大于0.85, 但在t趋向0和趋向正无穷时, 稳渗速率趋向无穷大和0, 不符合实际情况；通用经验模型拟合精度较差, 决定系数在0.614—0.982之间。这与刘目兴等^[1]在三峡库区研究结果一致。对比不同林地模型的模拟效果发现, 对祁连圆柏林地的拟合程度最高。刘洁等^[23]在元谋干热河谷研究发现通用经验方程最适于描述该流域的土壤入渗过程, 闫东峰等^[18]发现Kostiakov模型、Philip模型和Horton模型均适合丹江口水库地区, 李志等^[24]在武功山研究则发现通用经验模型更适合草甸土壤。说明不同研究区域、不同研究对象模型的适应性不同。

2.3 土壤水分入渗的影响因素 2.3.1 土壤理化性质对入渗的影响

本研究选取与土壤结构密切相关的16个因子(表 4)与土壤初渗速率、稳渗速率进行相关分析可以发现：

稳渗速率、初渗速率均与土壤容重呈极显著负相关关系(P < 0.01), 与总孔隙度、非毛管孔隙度、毛管孔隙度呈极显著正相关关系(P < 0.01)。容重是土壤的基本物理性质, 它的大小间接反映了土壤的紧实程度、透气透水能力^[25]。容重越大土壤越紧实, 则孔隙度随之减少, 土壤透气透水能力下降, 入渗速率减小^[12]。土壤孔隙度、容重与土壤入渗能力对应较好, 土壤容重小、孔隙度大的土壤入渗能力较强, 由表 1可知坡面向下林地土壤总孔隙、毛管孔隙度增大, 容重减小, 土壤渗透性能增强。

稳渗速率、初渗速率与＞3 mm、＞1 mm、＞2 mm、＞0.5 mm、＞0.25 mm、＞0.1 mm水稳性团聚体呈极显著正相关关系(P < 0.01), 与 < 0.1 mm团聚体呈极显著负相关关系(P < 0.01), 其中与＞0.25 mm的团聚体相关性最高, 表明入渗速率受＞0.25 mm团聚体的影响最大。土壤团聚体是反映土壤结构优劣的主要指标之一, 良好的团聚结构对保持土壤水分、肥力有着重要作用进而促进植物生长根系发育。土壤团聚结构越好, 其孔隙越多、持水性、通透性越好^[9]。二者还与有机质呈极显著正相关关系(P < 0.01), 相关系数较高(0.665、0.663)。土壤有机质对促进团粒结构的形成、维持土壤结构稳定有着重要作用, 有机质含量越多, 土壤团粒结构就越多, 从而具有良好的透气透水能力, 促进土壤入渗^[6]。

稳渗速率、初渗速率与土壤机械组成、土壤含水量均没有显著相关性(P＞0.05)。与高婵婵等人研究结果相同^[8]。

2.3.2 影响土壤入渗性能主导因素筛选

在SPSS中对以上16个与土壤渗性性能密切相关的因子进行主成分分析, 结果如表 5所示。表 5表明, 16个因子可概括为4个独立的主成分, 累计贡献率82.959%, 可用来解释对土入渗性能的影响。第1个主成分因子主要由水稳性团聚体决定, 贡献率为45.219%, ＞2 mm、＞1 mm、＞0.5 mm、＞0.25 mm的水稳性团聚体有较高载荷, 该类土壤有机质、水稳性团聚体含量高, 同时有较高的非毛管孔隙度和总孔隙度；第2个主成分因子主要由土壤机械组成决定, 方差贡献率为18.599%, 粘粒、粉粒含量高, 砂粒含量极低(-0.966), 三者载荷均较高；第3个主成分因子主要由土壤孔隙状况决定, 总孔隙度、毛管孔隙度、容重在第3主成分上载荷较高；第4个主成分因子主要由土壤含水量决定, 载荷量为0.940。

选择由主成分分析中载荷较高的11个因子共60组数据作为自变量分别和因变量稳渗速率(Y₁)、初渗速率(Y₂)进行多元逐步回归, 筛选出主要影响因子, 经检验二者符合正态分布。结果显示：

式中, X₁：容重, X₂：＞0.25 mm团聚体, X₃：总孔隙度, X₄：毛管孔隙度

可知, 影响稳渗速率的主要因子为容重、＞0.25 mm团聚体、总孔隙度, 影响初渗速率的主要因子为＞0.25 mm团聚体、总孔隙度、毛管孔隙度。为了判断已筛选的因子分别对稳渗速率和初渗速率的影响程度, 进一步进行通径分析。

由表 6可知, 初渗速率与相关因子的直接通径系数绝对值大小为＞0.25 mm团聚体质量分数 < 毛管孔隙度 < 总孔隙度, 总孔隙度的直接通径系数最大, 为0.939, 对初渗速率产生的直接正效应最强, 说明总孔隙度越大, 初渗速率越大；同时通过毛管孔隙度、＞0.25 mm团聚体间接对初渗速率产生的负效应最弱(-0.178)。＞0.25 mm团聚体对初渗速率的直接正效应最弱(0.404), 其通过总孔隙度和毛管孔隙度对初渗速率的间接正效应较强(0.361)。毛管孔隙度对初渗速率有着极显著的负向作用, 可见土壤毛管孔隙度越多, 初渗速率有减小趋势；但毛管孔隙度与初渗速率呈极显著正相关关系(0.544), 这是毛管孔隙度对初渗速率直接效应(-0.446)受到总孔隙度间接正效应(0.846)综合影响, 导致毛管孔隙度对初渗速率的直接效应不能通过二者的相关系数大小显示出来。通过计算决策系数可知, 毛管孔隙度决策系数绝对值最大, 且是负值, 说明是初渗速率的主要限制因子。(Durbin Watson2.213, 剩余通径系数=0.468)

由表 7可知, ＞0.25 mm团聚体质量分数与总孔隙度对稳渗速率表现出极显著的正向作用, 说明二者对稳渗速率的增大有明显的促进作用。容重对稳渗速率具有极显著的负向作用, 表明随着容重增大稳渗速率有减小趋势。＞0.25 mm团聚体的直接通径系数最大, 为0.460, 对稳渗速率产生的直接正效应最强, 同时通过容重、总孔隙度间接对稳渗速率产生的正效应最弱(0.315)；总孔隙度对稳渗速率的直接正效应最弱(0.271), 但其通过容重、＞0.25 mm团聚体对稳渗速率的间接正效应最强(0.518)。容重对稳渗速率的直接负效应较弱(-0.323), 通过＞0.25 mm团聚体、总孔隙度间接对稳渗速率产生的负效应较强。由决策系数可知, ＞0.25 mm团聚体单独作用排序第1, 容重排序第2, 总孔隙度排序第3, 所以＞0.25 mm水稳性团聚体是影响稳渗速率的主要决策因子。沿坡面向下, 容重减小降低了土壤紧实度, 提高总孔隙度和＞0.25 mm水稳性团聚体含量, 进而提高了土壤的稳渗速率^[15]；总孔隙度、＞0.25 mm水稳性团聚体含量沿坡面向下逐渐增大, 土壤容重减小孔隙性提高, 进而促进土壤入渗。(Durbin Watson 1.767, 剩余通径系数=0.426)

通径分析发现土壤物理因子与初渗速率、稳渗速率的剩余通径系数较大, 均在0.4以上, 说明存在重要因素未考虑, 如作物根系分布、根长密度、土壤生物等。研究表明根系作为植物的重要功能器官, 它的构型及生长状况反映了植被对环境的适应及对土壤的固持, 改良了土壤的紧实度, 强化了土壤的抗侵蚀能力, 是水土流失防治的关键^[26]。根系一方面在土壤中的穿插、分割等作用增加了土壤中的大孔隙, 改善土壤渗透性能；另一方面根系的代谢还会产生有机化合物胶结土壤颗粒, 促进团聚体的形成, 进而影响土壤的渗透性能^[27]。王国彬等^[28]研究发现土壤有效孔隙度和根长密度是影响表层土壤稳定入渗的主要因子；王意锟等^[29]发现土壤的渗透指标均随根长密度的减少而变小；闫东锋等^[18]研究发现土壤的入渗速率与根长密度、根体积密度及根系生物量均存在显著或极显著的相关关系, 相关系数均在0.85以上, 且数量多、密度大的0.5—1 mm、1—2 mm径级的根系改善土壤入渗性能的作用最强；徐少君等^[30]发现根系的根长密度、根表面积密度越大, 土壤抵抗水蚀的能力越强, 且径级越小的根系与土壤抗水蚀增强系数的相关性越大；王鑫皓等^[31]研究发现倾斜根、垂直根影响土壤水分的纵向运动, 数量越多土壤稳渗速率越大, 水平根系能阻止水分下渗, 数量越多稳渗速率越小。除根系在土壤中的缠绕串联改善土壤理化性质外, 土壤动物在土中的筑穴、挖掘、促进有机物分解等活动也进一步了改善土壤理化性质, 增强土壤通透性, 进而影响土壤中水分入渗^[29]。王意锟等^[32]发现在毛竹林中土壤动物的数量会随着土壤密度的增大而减少, 土壤动物的优势类群与土壤的有机质含量、pH、全N含量具有正相关关系。土壤入渗特征与多数大中型节肢动物、小型节肢动物的密度存在正相关关系, 这些动物的减少不利于土壤水分入渗^[29]。毛竹林的长期覆盖会限制土壤动物活动, 恶化土壤性质, 减弱土壤入渗, 最终导致林地退化^[29]。微生物是土壤动物的食物来源, 它们的生物量及活性会因动物取食而提高^[33]。土壤微生物作为土壤生态系统重要的组分之一, 参与众多土壤过程, 土壤酶活性作为生物催化剂, 具有重要作用。研究发现土壤渗透特性指标与细菌、放线菌含量和蛋白酶、脲酶活性极显著或显著相关, 初渗速率、平均渗透速率与多酚氧化酶活性显著相关, 说明土壤生物因子能够改善土壤环境, 加速土壤风化和养分循环, 最终提高土壤渗透性能^[34]。因试验时间关系未对以上因子进行详细调查, 将在后续研究中进行完善。

3 结论

(1) 不同植被类型的入渗速率存在显著差异, 初渗速率和稳渗速率均表现为青海云杉>祁连圆柏>华北落叶松>荒草地；同一坡位三种植被类型渗透速率不存在显著差异性, 但仍表现出上坡位和下坡位青海云杉>祁连圆柏>华北落叶松, 中坡位青海云杉＞华北落叶松＞祁连圆柏；同一植被类型在不同坡位的土壤入渗不存在显著差异, 青海云杉林地、祁连圆柏林地土壤渗透性能及华北落叶松林地稳渗速率均表现为上坡位最差、下坡位最好, 华北落叶松林地初渗速率在中坡位最大、下坡位最小.

(2) 对不同坡位的3种林地及荒草地的入渗过程进行模拟, 拟合精度最高的是Horton模型, 均在0.8以上；Kostiakov模型除在青海云杉林地上坡位拟合精度稍差(0.727), 其他部位均大于0.85；通用经验模型拟合精度总体较差, 决定系数在0.614—0.982之间, 但对青海云杉林地和祁连圆柏林地的中、下坡位及荒草地的拟合精度较高。

(3) 土壤入渗性能与孔隙状况、水稳性团聚体质量分数、有机质含量均有极显著正相关性, 与容重呈极显著负相关性。主成分分析结果表明影响因子可归为4个主成分, 分别为水稳性团聚体因子、土壤机械组成因子、土壤孔隙状况因子、土壤水分因子, 累计贡献率为82.959%。多元分析可知影响初渗速率的主要因子为＞0.25 mm团聚体、总孔隙度、毛管孔隙度, 影响稳渗速率的主要因子为容重、＞0.25 mm团聚体、总孔隙度。通径分析可知, 总孔隙度是影响初渗速率最主要的决策因素, 而毛管孔隙度是主要的限制因子；＞0.25 mm团聚体质量分数是影响稳渗速率的主要决策因素。二者的剩余通径系数较大(0.468、0.426), 说明渗透速率的影响因子还有部分未考虑。