积雪特性参数分析及雪深模型建立

付 强，彭 莉，汪恩良，李天霄，刘 东，崔 嵩

积雪特性参数分析及雪深模型建立

付 强，彭 莉，汪恩良，李天霄，刘 东，崔 嵩

（东北农业大学水利与土木工程学院，哈尔滨 150030）

通过2015～2016年冬季积雪观测试验，测量自然条件下雪深、雪层温度、雪层密度、雪层液态含水率等积雪特性参数，分析积雪特性参数影响因子和发育变化规律，构建径向基神经网络雪深模型。结果表明，空气、地表温度及水汽压是影响雪参数发育三大主要气象因子。在时长为124 d覆雪期中，雪深最大平均融雪速率达3.20 mm·d-1，出现在融雪后期。雪层密度随雪深增加波浪式上升，其中稳定期全层密度0.1～0.4 g·cm-3；雪层液态含水率0.5%～3.5%，分层液态含水率中层较大，积雪表层和底部波动较小；全层积雪孔隙率平均值为72.3%，积底层孔隙率廓线下降最快，极差最大。经两次优化后，雪深模型最大绝对误差为0.614 1 cm，最大相对误差为5.85%，模拟精度控制在6 mm以内。

雪深；雪密度；液态含水率；雪孔隙率；雪深模型

降雪约占地球表面全部降水量5%[1]，随区域地表水质持续下降，积雪供给人类生存用水比重越来越高，尤其是农业用水。积雪高反照率、低导热率、大热容量等特性可有效阻隔雪面上低温入侵[2-3]，大幅减少大气和地表间热量交换，减少土壤内部热量散失，避免土壤过度降温，利于农作物发育，保护植被不被冻伤，影响土壤冻结深度、速率、冻胀量和热质迁移等地表自然过程[4-5]。积雪深度作为表征积雪物理特性基础参数[6]，在计算雪水当量[7]、模拟融雪径流[8]、估计积雪导热系数反射、透射和吸收率，以及建筑物雪荷载计算等方面均为重要参变量。

目前融雪模型研究方向主要集中于大尺度水文模型应用[9]，由于雪深变化规律呈非线性趋势，一般线性拟合难以呈现其特征。许剑辉等通过空间插值方法对雪深作多元非线性回归分析[10]，白淑英等选用长时间序列作大尺度雪深时空规律性研究[11]，许剑辉等借助遥感参数反演结合水文融雪模型研究区域雪情[12-14]，Godio和Tiwari等基于能量平衡模型以冻土-雪被连续体为整体对象作水热运移分析[15-16]，以小尺度精细试验结合统计学分析及相关算法建立雪深模型的区域性研究较少。

松嫩平原冬季覆雪时间长达5个月，掌握松嫩平原积雪特性规律，对保障区域粮食和水资源安全具有重要意义[17]。因此，本文通过积雪特性观测试验，分析雪层密度、雪层液态含水率、雪层孔隙率等积雪参数变化规律，采用径向基神经网络构建雪深模型，模拟雪深动态变化，预测雪深和雪水当量，为预报雪情和制定春汛防治措施提供理论依据。

1 材料与方法

1.1 区域概况

松嫩平原属全国三大积雪区（青藏高原地区、新疆天山地区、东北地区）之一，是典型季节性积雪区。其耕地面积占中国耕地面积1/9，粮食商品量占全国1/6，是我国最大商品粮生产基地。夏季高温潮湿，冬季寒冷干燥，多年平均降水量为529 mm，多年平均积雪覆盖期为110 d，最大冻深达180 cm，可为积雪试验研究提供条件[18]。该区域内积雪时间通常为每年11月初至次年3月中旬。经人工取样分析可知，该地区0～30 cm土层为黑色壤土，30 cm以下为黑色黏土。

1.2 试验方案

积雪特性观测在黑龙江省哈尔滨市东北农业大学水利与土木工程学院综合试验场开展，自每年初雪开始，至次年积雪完全融化结束。选取3组10×10 m2无植被残茬覆盖开阔平坦均质壤土裸地为试验组，处于自然积雪覆盖状态。无任何积雪覆盖10×10 m2裸地作为对照组。

据2013～2015年观测经验，本试验观测时间设为2015年11月8日到2016年3月11日每日8:00，采取3点取样法，从等边三角形顶点位置开始取样，随雪深增加，积雪剖面每5 cm划分为1层，每层取5点测量求平均值。使用Snow Fork雪特性分析仪测量各层雪密度及液态含水率，探针式温度计测定雪层温度，量雪尺和自动气象站雪深传感器测雪深。试验场内架设有锦州阳光PC-3型自动气象站，可采集：空气温度、空气湿度、大气压强、风向风速、雪深、总辐射、散辐射、直接辐射、长波辐射、净辐射等气象数据，步长设为30 min（见图1）。

在2015～2016年观测数据中，选取2015年12月～2016年2月积雪较稳定时期分析雪密度、液态含水率和雪孔隙率参数特性；结合历年观测资料，构建径向基神经网络雪深模型时，选用2013～2016年数据为训练样本，2016～2017年数据为预测样本。

1.3 方法

1.3.1 积雪孔隙率计算

积雪孔隙率测量一般可通过积雪其他物理参数计算。理想状态下只需积雪密度和冰密度，但自然条件下积雪实际由冰、水汽和液态水3种形态物质组成。水汽质量可忽略不计。融雪期随温度升高，积雪中液态水体积不断增大，加上冰晶作用，直接影响孔隙率。因此，计算时必须考虑积雪中液态水体积[19]。

式中，V为雪体积，Vi为雪中固态冰体积，Vw为雪中液态水体积，Va为雪中空隙体积。m为雪质量，mi为雪中固态冰质量，mw为雪中液态水质量。

图1 试验区域地理位置及场地布置Fig.1 Geographic position and experimental site layout of the measured fields

故有，积雪孔隙率

得推导公式：

式中，冰密度取ρi=0.9168 g·cm-3；雪密度ρ和体积含水率Wvol均由Snow Fork雪特性分析仪测得。

1.3.2 气象因子主成分提取

由于气象指标较多且具相关性，采取主成分分析法对气象因子作降维处理，将多个气象指标重新组合为包括空气温度、净辐射、水汽压、相对湿度、日总辐射量及每日测量地表温度数据等综合指标[20]。日总辐射量以日辐射值总和作为代表值，其他指标均以日平均值作为代表值，对以上气象数据作归一化处理。主成分分析实质是计算指标组成方差-协方差矩阵特征向量及其对应特征值。经运算，六大气象指标主成分特征值及贡献率见表1。

采用最小二乘回归算法确定原始变量与主成分间线性关系，将气象因子成分载荷除以对应所属主成分特征值平方根，取累计贡献率在85%以上，前3个主成分为代表，气象因子系数矩阵见表2。

表1 主成分特征值及贡献率Table1 Principal component eigenvaluesand contribution rates

表2 主成分得分系数矩阵Table2 Principal component scorecoefficient matrix

由表1、表2可知，第一主成分主要反映空气温度、地表温度及水汽压特征，累计贡献率为40.412%；第二主成分以净辐射和日总辐射量特征为主，二者累计贡献率达67.390%；第三主成分则突出显示相对湿度影响，三个主成分累积贡献率达85.672%，集中代表气象因子主要特征。

2 积雪特性参数分析

2.1 雪深变化分析

图2 试验期雪深和降雨变化曲线Fig.2 Curve of snow depth and rainfall during the accumulation-thaw period

绘制自然状态下裸地雪深变化曲线如图1所示。根据寒区季节性积雪特点，将研究周期划分为积累期，稳定期，融雪期。

由初雪到12月22日，雪深逐渐攀升至最高点168 mm，以新降雪过程为主，此阶段为积雪期。大中降雪集中在12月上中旬，雪层较松软，雪粒孔隙率大，随时间推移，雪粒抱团粘连，雪层中新雪和细粒雪迅速沉降，孔隙率降低，雪剖面高度随之下降。

由于雪深明显上升且降雪量偏大，雪层受上覆积雪压力及自身重力作用密实化过程不断增强，达最大雪深后第二阶段，稳定期。太阳辐射微弱，气温极低，土壤已基本冻结，雪盖仅在表层发生单向融化，外界输入能量较少，积雪自身密实化过程缓慢，雪深以锯齿状上下浮动为主，此阶段平均积雪融化速率为0.45 mm·d-1，雪情状态较稳定。

融雪期自2月4日至积雪完全融化，期间虽偶有阵雪，气温已开始逐步回升至零上，雪盖上界面受太阳和大气热辐射持续增强，下界面土壤逐步解冻，感热、潜热通量陆续回升，积雪剖面底部深霜层逐渐消融，雪层开始发生双向融化，融雪水不断下渗，积雪内部不同雪层之间水量和能量交换活跃，平均融雪速率高达3.20 mm·d-1，少量融雪水出流。

2.2 雪层密度变化分析

随新雪累积消融和变质作用，各时期不同雪层参数时空分布变化显著。选取2015年12月～2016年2月积雪较稳定阶段观测数据作雪密度及含水率等特性参数分析，除去特异值计算三组裸地日平均值，试验周期内雪层密度分布见图3。

图3 试验期积雪各层密度分布Fig.3 Variations in the density of each snow cover layer

由图3可知，积雪各层密度均呈两侧波动中间平稳变化趋势，最小雪密度0.0585 g·cm-3出现在不稳定积雪期上层，最大雪密度0.2471 g·cm-3出现在融雪期下层，平均雪密度为0.05～0.2 g·cm-3。不稳定积雪期气温尚未完全降至零下，温度日差较大，新降雪疏松多孔，土壤尚未冻结，水汽蒸凝和夜间重冻结过程剧烈，雪层冻结消融交替周期短，随表层低密度雪持续补充及自身重力挤压密实化作用，中层和底层密度增长速率偏快，尤其是底层。

融雪期气温回暖，昼夜温差大，辐射量差异显著，不断促进融雪水运移，雪水大量出流，雪层空隙被雪水填满，呈多孔状态，密度分布呈振荡性特点，在不稳定积雪期和融雪期雪密度上下波动明显，与中心点振荡间距偏大。而平稳积雪期降雪量逐步减少，气温均在零下并平稳降低，太阳辐射对上层影响减弱，底层雪-土界面已冻结，热通量交换几乎为零，仅中部受上下双层挤压效应和融雪水运移影响，密度上升较活跃。降雪层空隙中融雪水和固态冰晶间达动态平衡，相变过程稳定。雪深达峰值后以小幅度波动变化为主。积雪进入缓慢密实化阶段，相应密度呈缓慢上升趋势。

从12月下旬开始积雪底层逐渐出现深霜颗粒，雪层剖面自上至下可分辨新雪、细粒雪、中粒雪、粗粒雪及深霜，积雪密度随粒径增大而增大。随雪深及负温度梯度增大，粗颗粒松散型积雪占比上升，到次年1月末，全雪层中近1/3厚度底部积雪转化为霜层，整个融雪期深霜层厚度可达积雪层厚度40%，至2月上旬深霜层占比不断攀升，附加降雨促进温度梯度转化，大量水汽由高温层向低温层迁移，水汽压较小晶体表面发生升华，而水汽压较大晶体表面则发生凝华，晶体结构变化。两场降雨后，剩余积雪层除雪盖表面为不连续新雪或疏松雪层，80%变质为深霜颗粒，已完全冰晶化。

2.3 雪层液态含水率变化规律分析

融雪水在雪层中运移直观表现形式为液态含水率变化，由于存在雪粒晶体表面张力及晶体空隙毛管吸附力，液态水暂时吸附于固态雪层中，与冰晶颗粒不断发生相变，同时融雪水受自身重力影响发生垂向浸润下渗，地面坑洼有坡度处发生侧向渗漏，不断运移使积雪各层之间液态含水量产生差异，随积雪深度增加，日均液态含水率显著增大，极值达3.627%，出现在底层。当融雪水超过积雪持水能力时，才可能发生出流，对应图中实际未出流。融雪期一般不发生出流，为黑龙江地区积雪发育状态一般规律。

图4 试验期积雪各层液态含水率分布Fig.4 Variations in the liquid water contents of each snow cover layer

2.4 积雪孔隙率变化分析

由图7、8可知，整个覆雪周期内全层积雪孔隙率呈下降趋势，融雪期震荡幅度显著增大。全层积雪孔隙率平均值为72.3%，由积累前期80%到积累后期降至70%，进入融化期积雪孔隙率在60%～70% ，最低达49%，雪层空隙中融雪水流动使孔隙率高低值交错明显，经计算全周期内孔隙率日变化速率为-2.8%·d-1。3层中，底层孔隙率平均在62%，融雪期更是低至55%。孔隙率减小，一方面源于气温剧烈升高积雪融化加速，冰晶融冻变质，冰晶间空隙被压缩；另一方面由于大量雪融水充斥孔隙，减小积雪孔隙率。可验证孔隙率是积雪持水能力关键指标之一。

中国天山西部山区冬季积雪类比试验，积雪孔隙率变化范围均在60%～80%[21]，哈尔滨地区积雪孔隙率交替增减下降变化规律与此试验结果一致。相比于日本北部、芬兰属温带海洋性气候，纬度低，气温相对较高，中国新疆、天山西部等内陆干旱地区属温带大陆性气候，受温带季风影响。哈尔滨地区纬度偏高，冬季寒冷，雪量较天山地区偏小，融雪期辐射强烈升温快，孔隙率变化跨度更大，融雪后期下降程度最大。

计算图5中分层积雪孔隙率并作雪层特征参量分布统计如表3。积雪各分层孔隙率变化界面处能量交换活跃，水热运移剧烈，最先开始融化出流，孔隙率变化幅度最大。

图5 积雪分层孔隙率变化趋势Fig.5 Variations in the temperatureof each snow cover layer

表3 积雪分层孔隙率特征参数值Table 3 Character parameters of snow porosity in each layer

3 雪深模型建立

积雪深度是实现雪水当量和融雪径流等指标计算和模拟基础，在掌握积雪特性参数和相关气象参数详细变化规律基础上，利用径向基神经网络算法模拟雪深[22]。

3.1 径向基神经网络原理

径向基函数RBF（Radial Basis Function）能逼近任意非线性函数，处理系统难以解析的规律性，泛化能力和学习收敛能力良好，广泛应用于多变量插值[23]。径向基网络RBFNN包括输入层、隐层、输出层，拓扑结构如图6所示[24]。

由图6可知，信号由输入层节点传递到隐层，位于隐层节点径向基函数对输入向量作非线性变换，输出层是线性，通过线性加权为作用于输入层RBF激活模式提供响应。

记RBFNN隐层神经元（隐层节点）个数为K,则网络输出形式如下：

图6 RBFNN拓扑结构Fig.6 Topology structureof RBFNN

其中：w0为偏差（bias），wk(k=1,2,...,K)为隐层和输出层之间连接权值。

3.2 雪深模型建立及优化

利用雪深数据构建径向基神经网络。选用2013～2016年雪深360组数据为训练样本作网络训练。通过自相关分析和偏相关分析，分析数据系列分布特征。

由图7可知，时间序列18～70 d数据处于置信区间内，最终数据趋于零（积雪完全消融），可认定该时间序列自相关函数拖尾。由图8可知，本数据序列为2阶截尾，偏相关系数基本处于上下界内。选择2016～2017年雪深数据为试验样本作模拟检验。

图7 自相关分析Fig.7 Auto correlation analysis chart

图8 偏相关分析Fig.8 Partial correlation analysis chart

调用径向基函数net=newrb(P,T,goal，spread,MN,DF)。其中，P为输入矩阵，T为目标输出矩阵，均方误差目标goal取0.01，扩展系数spread默认为1，MN为隐层神经元最大个数，取为样本值124，DF默认为25。经多次计算确定参数值为goal=0.001.spread=2，MN=124，DF=1，程序运行到116次，网络误差为0.000961351＜0.001，符合要求。

在雪深自模拟基础上，加入气象因子对模型第1次优化，3个主成分作为RBFNN输入因子。以2013～2016年3个融雪期气象数据对网络不断训练，程序运行至86次时，网络误差为0.000947081＜0.001，符合要求。在加入气象因子基础上，加入积雪特征参数二次优化模型。对积雪期实测逐日分层雪特性参数积雪温度、雪密度和雪液态含水率按照各自分层深度求积分，以积分值为代表值归一化处理后作为输入因子添加在原有网络中，建立二次优化模型模拟雪深。以2013～2016年3个融雪期为训练样本不断训练网络，程序运行至72次时，网络误差为0.000924261＜0.001，符合要求。

选择2016～2017年雪深数据为样本分析模型检验。

表4 3次模型模拟结果误差对比分析Table 4 Error analysis of three model simulation results

3.3 模型误差评价

比较3次预测结果模拟值，间隔10 d数据结果对比见表4。

由表4可知，绝对及相对误差均随优化程度增加而逐级降低。3个模型绝对误差平均值分别为：1.488 8、0.628 9和0.165 2 cm，最大值分别为3.649 6、1.961 8和0.614 1 cm。绝对误差平均值分别为：12.46%、4.06%和1.15%，最大值分别为47.63%，11.02%和5.85%。经过二次优化，模型最终精度误差范围可稳定控制在6.141 mm以内，误差占比在1.15%左右。说明本模型适用于雪深预测。

4 结论

通过积雪累积消融过程水热运移微观分层研究，初步确定积雪特性参数雪密度、液态含水率、孔隙率等变化规律和水热运移特点。结论如下：

a.东北寒区季节性积雪区在积雪消融过程中密度偏大，可由0.05 g·cm-3升至0.3 g·cm-3；液态含水率上升，平均在0～4%。孔隙率呈波浪式下降特征，廓线最低峰值达49%，较海洋性积雪区和干旱性积雪区偏小，表现出受季风影响的大陆性积雪区特点。

b.受上界面能量输入影响最大积雪表层最先发生能量转换，剖面传递过程使中层和底层表现一定滞后性，底层受上层覆雪和下垫面地热双重影响在整个积雪周期中能量转换最为活跃，变化幅度最大。

c.雪深模型模拟结果说明，在东北季节性积雪区，以径向基神经网络算法为基础建立雪深模型可行，二次优化雪深模型绝对误差控制在6 mm以内，最大相对误差为5.85%，在试验采集期，要不断细化和扩展雪深相关参数观测，设置不同下垫面条件等。未来可扩大研究区域、增加时间序列长度并结合遥感反演，多个角度优化网络并丰富模型边界条件，提高雪深模型稳定性。

[1] 赵求东,刘志辉,秦荣茂,等.融雪模型研究进展[J].新疆农业科学,2007,44(6):734-739.

[2]Domine F,Barrere M,Sarrazin D,et al.Automatic monitor⁃ing of the effective thermal conductivity of snow in a low-Arctic shrub tundra[J].Cryosphere Discussions,2015,9(3):1633-1665.

[3] Thackeray C W,Fletcher C G.Snow albedo feedback:Cur⁃rent knowledge,importance,outstanding issues and future di⁃rections[J].Progress in Physical Geography,2016,40(3):392-408.

[4] Endrizzi S,Gruber S,Dall A M,et al.GEOtop 2.0:simu⁃lating the combined energy and water balance at and below the land surface accounting for soil freezing,snow cover and terrain effects[J].Geoscience Model Development,2014,7(6):2831-2857.

[5] Yi Y,Kimball J S,Rawlins M A,et al.The role of snow cover and soil freeze/thaw cycles affecting boreal-arctic soil carbon dynamics[J].Biogeosciences Discussions,2015,12(14):11113-11157.

[6] Nosaka M,Sasaki H,Murata A,et al.Bias correction of snow depth by using regional frequency analysis in the nonhydrostatic regional climate model around Japan[J].Scientific Online Letters on the Atmosphere Sola,2016,12(16):5-9.

[7] Rittger K,Bair E H,Kahl A,et al.Spatial estimates of snow water equivalent from reconstruction[J].Advances in Water Resources,2016,94(3):45-63.

[8] Stojkovic M,Jacimovic N.A simple numerical method for snowmelt simulation based on the equation of heat energy[J].Water Science and Technology,2016,73(7):1550-1559.

[9] Manzan C A,Muir D,Kirk J,et al.Temporal variation in the deposition of polycyclic aromatic compounds in snow in the Athabasca Oil Sands area of Alberta[J].Environmental Monitoring and Assessment,2016,188(9):1-12.

[10] 许剑辉,舒红,李杨,等.基于MODIS雪盖数据的北疆雪深多元非线性回归克里金插值[J].国土资源遥感,2015,27(3):84-91.

[11] 白淑英,史建桥,沈渭寿,等.近30年西藏雪深时空变化及其对气候变化的响应[J].国土资源遥感,2014,26(1):144-151.

[12]许剑辉,舒红.基于Triple-Collocation的地面观测与卫星遥感数据融合的雪深反演[J].武汉大学学报:信息科学版,2015,40(4):469-473.

[13]唐志光,李弘毅,王建,等.基于多源数据的青藏高原雪深重建[J].地球信息科学学报,2016,18(7):941-950.

[14] 王子龙,胡石涛,付强,等.积雪参数遥感反演研究进展[J].东北农业大学学报,2016,47(9):100-106.

[15] Godio A.Multi Population Genetic Algorithm to estimate snow properties from GPR data[J].Journal of Applied Geo⁃physics,2016,131(1):33-44.

[16] Tiwari S,Kar S C,Bhatla R.Examination of snowmelt over Western Himalayas using remote sensing data[J].Theoretical and Applied Climatology,2016,125(1):227-239.

[17]李天霄,付强,孟凡香.黑龙江省农业水土资源时空匹配格局研究[J].东北农业大学学报,2017,48(6):51-58.

[18]付强,蒋睿奇,王子龙,等.不同积雪覆盖条件下冻融土壤水分运动规律研究[J].农业机械学报,2015,46(10):152-159.

[19] 洪雯.利用Snow Fork测量积雪孔隙率及其应用[C].第28届中国气象学会年会,2011.

[20] Freedman D,Pisani R,Purves R.Statistics:Fourth edition[M].New York:W.W.Norton&Company,2007.

Analysis on snow characteristic parameters and construction of snow depth modelin

FU Qiang,PENG Li,WANG Enliang,LI Tianxiao,LIU Dong,CUI Song
(School of Water Conservancy and Civil Engineering,Northeast Agricultural University,Harbin 150030,China)

By snow observation experiment in 2015-2016,snow parameters of depth,temperature,density and liquid water content were collected under the natural conditions,the influence factors and the development law of snow parameters were analyzed,and the snow depth model was constructed by radial basis neural network.The results showed that air temperature,surface temperature and water vapor pressure were the three main meteorological factors affected snow parameters.In the 124 d during the whole snow cover period,the maximum snow melting rate of the snow depth reached 3.20 mm·d-1and appeared in the late snowmelt period.With the increase of deep snow density wave risen,the snow density varied between 0.1 to 0.4 g·cm-3during the stable period,and the average liquid water content of snow increased from 0.5%to 3.5%.The bottom snow layer exhibited larger parameter variations compared to those in the surface and middle layers.The average porosity of snow was 72.3%in the entire snowpack,ranged from 4%to 19%during the accumulation period and from 7%to 25%during the snowmelt period,in which the porosity of the bottom layer displayed the fastest decline.The maximum absolute error of the snow depth model after two times optimization was 0.6141 cm,the maximum relative error was 5.85%,and the simulation precision could be controlled within 6 mm.

snow depth;snow density;liquid water content;snow porosity;snow depth model

2017-12-03

国家自然科学基金（51679039）

付强（1973-），男，教授，博士生导师，研究方向为农业水土资源系统分析。E-mail：fuqiang 0629@126.

S161.6

A

1005-9369（2017）12-0036-10

时间2017-12-18 13:40:07 [URL]http://kns.cnki.net/kcms/detail/23.1391.S.20171218.1339.004.html

付强,彭莉,汪恩良,等.积雪特性参数分析及雪深模型建立[J].东北农业大学学报,2017,48(12):36-45.

Fu Qiang,Peng Li,Wang Enliang,et al.Analysis on snow characteristic parameters and construction of snow depth modeling[J].Journal of Northeast Agricultural University,2017,48(12):36-45.(in Chinese with English abstract)