1961—2016年松嫩平原月尺度下水分盈亏态势时空分异特征

路 中，雷国平，张露洋

(东北大学土地管理研究所，辽宁沈阳 110169)

松嫩平原处于全球环境变化速率最大的东亚季风区，既是气候变化敏感区［1］，又是气候变暖的强增温区域［2-3］。干旱与半干旱是该地区的显著特点，区域内大部分地区处于水分亏缺状态［4］。特别是随着全球平均气温的升高［5-6］，松嫩平原水分亏缺呈加剧态势，春旱与春夏连旱的频率和强度增加，这直接影响了区域内的粮食生产和经济发展［7］。因此，在全球气候变暖、干旱发生频率和强度逐渐加剧的大背景下，厘清该地区作物生长周期内水分盈亏的演变规律，对于该地区的粮食生产安全和经济可持续发展具有重要意义。

目前关于区域内水分盈亏状况的研究方法大致可以分为2类:一类是运用气象干旱指数研究区域水分盈亏状况，包括综合气象干旱指数(CI)［8］、标准化降雨蒸散指数(SPEI)［9-10］、帕默尔干旱指数(PDSI)［11-12］、Z指 数［13］以 及 干 旱 侦 测 指 数(RDI)［14-15］等;另一类是从水土资源平衡的角度研究区域水分盈亏情况，依据流域水量平衡模型［16-17］和傅抱璞(Budyko)模型［18-19］等方法求得实际蒸散量，并根据作物需水量计算水分盈亏量和水分盈亏率，进而对研究区域内的水分盈亏情况进行研究。总体来讲，第1类方法可以有效估算长时间范围内研究区域的水分盈亏情况，且可以分析作物生长周期内不同时间段的水分盈亏情况，但该方法仅依据气候条件进行区域水分盈亏情况的研究，未考虑土地资源对于水资源的需求情况，因此评价结果的精确性受到一定影响;第2种方法充分考虑了水资源与土地资源之间的相互关系，但学者多以年为时间尺度计算年平均实际蒸散量，这样可以有效排除土壤含水量变化带来的影响，但对于生长周期内实际蒸散量的求解成为难点。

基于此，选择季节性旱涝表现明显的松嫩平原作为研究区域，利用土地利用数据、GLDAS土壤含水量数据和长时间序列气象数据，结合区域农田水量平衡模型求解月尺度上的实际蒸散量，在此基础上构建水分盈亏模型，进行月尺度上的松嫩平原水土资源平衡效应研究，定量研究松嫩平原作物生长周期内不同月份水分盈亏态势，以期为松嫩平原水土资源合理利用、农业生产环境协调布局以及实现区域可持续发展提供科学依据和决策支持。

1 研究区概况

松嫩平原位于黑龙江省西南部和吉林省西北部，西临大兴安岭东麓的丘陵台地，北部和东部以小兴安岭和长白山外围的山麓、丘陵台地为界，南部为以松辽分水岭为主的隆起带，是中国最大的平原之一，地理位置介于 42°49"～49°12"N，121°38"～128°33"E之间，主要包括吉林省和黑龙江省59个县(市、区)，总面积为21.6万km2。该区域主要地貌类型为低山丘陵、冲积平原、台地和洼地，气候类型为温带大陆性半湿润、半干旱季风气候，全年平均气温呈自北向南逐渐增高趋势，等值线近纬向分布，大部分地区的年降水量为400～600 mm，自东向西逐渐减少，6—8月降水量占全年降水量的60%～70%，易形成干旱和洪涝灾害。松嫩平原地势平缓，土质肥沃，黑土和黑钙土占60%以上，其他土壤类型包括草甸土、暗棕壤、沼泽土、盐土、碱土、风砂土、栗钙土和水稻土等。松嫩平原是中国重要的商品粮生产地区之一，耕地面积为13.55万km2。粮食作物以春小麦、玉米、高粱、谷子为主，局部地区栽种早熟的粳稻。经济作物以大豆、甜菜、亚麻为主。

2 数据来源与研究方法

2.1 数据来源

研究所用数据包括:2000和2010年2个基准年的30 m分辨率地表覆盖数据产品GlobeLand30;1961—2016年5—9月研究区附近48个站点的气象数据，主要包括降水、气温(最高气温、平均气温和最低气温)、相对湿度、降雨天数、风速和日照时数等;国家地球系统科学数据共享平台(http:∥www.geodata.cn)提供的中国1∶400万土壤类型图(1980s与2000s);全球陆面同化系统提供2种版本的数据集，即GLDAS-1和GLDAS-2，2种版本的数据集可提供气压、气温、风速、蒸散发、降水、降雪和土壤湿度等数据，时间为1961—2016年4—10月;基于监督分类和目视解译相结合的方式进行松嫩平原不同时间节点的耕地信息提取，得到30 m分辨率的耕地利用分布图，时间节点为1980、1995、2005和2015年。

2.2 研究方法

2.2.1 作物需水量估算

作物需水量是作物在土壤水分适宜的条件下进行正常生长发育且能达到高产潜力值时所需要的水量，是耕地水分盈亏计算过程中必不可少的因素［20］。为了估算大区域范围内月尺度上作物需水量的空间分布，采用世界粮农组织(FAO)公布的作物需水量公式计算，具体公式为

式(1)中，TE，t为作物生长期内t月的需水量，mm;Kt为作物生长期内t月的作物系数;TE0，t为作物生长期内t月的潜在蒸散量。Kt表示作物需水量与其潜在蒸散量之间的比值，其可靠性决定了作物需水量计算结果的准确性。为了有效估算研究区域内的作物需水量，参考黄志刚等［1，21］的研究成果，确定研究区域内的主要作物系数(表1)。

表1 松嫩平原主要农作物生育时期内各月份作物系数Table 1 The crop coefficient of each month during the growth period of main crops in Songnen Plain

对于潜在蒸散量，采用目前应用比较广泛的Penman-Monteith模型［22］进行计算，具体公式为

式(2)中，Δ为饱和水汽压曲线的斜率，kPa·(℃)-1;γ为干湿计常数;Rn为地面净辐射，MJ·m-2·d-1;G为土壤热通量;u2为2 m高处风速，m·s-1;es为饱和水汽压，kPa;ea为实际水汽压，kPa;T为平均气温，℃。

2.2.2 作物实际蒸散量估算

作物实际蒸散量指在实际降雨量和土壤含水量条件下，农田土壤、土壤水分流失和植物蒸腾的总耗水量。选用区域农田水量平衡方程求解月实际蒸散量，具体公式为

式(3)中，Pt为t月时的降雨量，mm;TES，t为作物生长周期内t月的实际蒸散量，mm;Qt为t月时的径流量，mm;Sd，t和Su，t分别为t月时的降雨入渗量和毛管上升水量，mm;St+1和St分别为t+1和t月初时的土壤含水量，mm;ε为误差项，mm。因笔者主要研究耕地的水分盈亏情况，土壤含水量取1 m深度的数据，且Sd，t和Su，t数据值较小，可相互抵消，所以将其计入误差项，计算中忽略不计。因此，该研究使用的区域农田水量平衡方程可以表述为

运用SCS-CN模型计算地表径流量［23］，为了有效排除融雪径流对于地表径流的影响，将提取的降雪数据叠加到5月降水数据进行计算，其公式为

式(5)中，Ia为降雨初损量，mm;S为可能最大滞留量，mm。其中，Ia=λ×S，λ表示区域参数，其取值范围为0.1～0.3。依据研究区域的实际情况和前人的研究成果［24-25］，该研究中λ取值0.2。S的计算公式为

式(6)中，NC为径流曲线值，其值由土壤前期的湿润程度、坡度、土壤属性和土地利用状况综合决定。参考美国《国家工程手册》和研究区域的实际环境特征来确定中等湿润(AWCⅡ)时的径流曲线值NC，Ⅱ。对于干旱(AWCⅠ)和湿润(AWCⅢ)条件下的径流曲线值NC，Ⅰ和NC，Ⅲ，具体计算公式［26］为

对于土壤干旱情况，依据次降雨前5 d的降雨量进行确定。如果前5 d的降雨量小于35.6 mm，则土壤处于干旱条件;如果前5 d降雨量大于53.3 mm，则土壤处于湿润条件;其余处于中等湿润条件。

为了确定1961—2016年逐年各月研究区域土壤含水量的变化情况，利用全球陆面同化系统提供的GLDAS-1和GLDAS-2数据集，得到各年生长周期内各月月初的土壤(深度1 m)含水量数据。GLDAS数据精度已得到了认证［27-28］，可满足研究需求。

2.2.3 水分盈亏模型

依据月作物需水量和月实际蒸散量构建水分盈亏模型，以此来计算研究区域内不同月份作物的水分盈亏状况，具体公式为

式(9)中，Wt为耕地在生长周期内t月时的水分盈亏状况，mm。

2.2.4 水分盈亏指数模型

为分析松嫩平原作物干旱严重程度，构建水分盈亏指数来表征研究区域内不同月份作物的水分盈亏等级，进一步分析松嫩平原月尺度上水土资源平衡效应。具体公式为

式(10)中，WD，t为作物在t月时的水分盈亏指数。参考杨艳昭等［17］的研究，将水分盈亏指数划分为正常(≤0.15)、轻度缺水(＞0.15～0.30)、中度缺水(＞0.30～0.45)、重度缺水(＞0.45～0.60)和严重缺水(＞0.60)5个等级。

2.2.5 趋势分析

采用一元线性回归分析法研究1961—2016年松嫩平原逐年时间序列数据的变化趋势，其公式为

式(11)中，Ts为1961—2016年各时间序列数据年际变化的一元线性回归方程的斜率，即趋势变化率;n为年数，该研究取值55;Ti为研究区内第i年的时间序列数据值。Ts＞0表明时间序列数据总体呈增加趋势，反之为减少趋势。

2.2.6 小波变换原理

Morlet小波分析是指利用固定窗口来揭示事件发生的时间和频率，它能够有效揭示时频发生的细节信息，适合分析非稳定信息［29-30］。对于给定小波函数ψ(t)，水分亏缺指数时间序列f(t)∈L2(R)的连续小波变换定义为

式(12)中，Wf(a，b)为小波实部系数;a为尺度因子，表示频域特征;b为时间因子，表示时域特征;t为复共轭函数。实际水分亏缺指数时间序列是离散的，因此利用离散小波变换对时间序列数据进行分解，得到不同时间尺度的小波系数。这些小波系数是时间序列数据在不同时间和空间位置上的投影，能用来描述水分亏缺指数时间序列数据的内在结构、性质和变化特性。

3 结果与分析

3.1 潜在蒸散量的时空变化特征分析

依据世界粮农组织提供的Penman-Monteith模型，并参考研究区域的实际情况，计算1961—2016年作物生长周期(5—9月)内各气象站点的潜在蒸散量，并利用ArcGIS软件分析潜在蒸散量的空间分布情况及其变化趋势，结果见图1。

图1 生长周期内松嫩平原潜在蒸散量的空间分布及变化趋势Fig.1 The spatial distribution and variation trend of potential evapotranspiration in growing period

由图1可见，1961—2016年松嫩平原作物生长周期内潜在蒸散量在432.9～605.8 mm之间，均值为528.1 mm。由于研究区地处大兴安岭余脉背风坡，降雨、气温、风速以及潜在蒸散量等气候因素受到地理位置影响，潜在蒸散量整体空间分布格局表现为自西南向东北逐渐减少。潜在蒸散量高值区域主要出现在白城市、大庆市西部和齐齐哈尔市北部，潜在蒸散量在552.9～605.8 mm之间，均值为575.4 mm;低值区域主要出现在黑河市、伊春市中东部和哈尔滨市中部，潜在蒸散量在432.8～506.1 mm之间，均值为472.3 mm。生长周期内各月份潜在蒸散量与总平均潜在蒸散量的空间分布状况基本一致，均表现为由西向东逐渐减少的趋势。5—9月潜在蒸散量分别占生长周期内总平均潜在蒸散量的 24.4%、22.9%、20.1%、17.7%和 14.9%。松嫩平原降雨量主要分布在7—8月，致使研究区域内容易出现春旱。

松嫩平原地区潜在蒸散量增减幅度存在明显的空间差异(图1～2)，其中松原市南部、长春市南部和四平市西北部均出现比较明显的上升趋势，平均增长幅度为0.210 mm·a-1;齐齐哈尔市东南部、大庆市北部、绥化市西部和哈尔滨市西部具有明显下降趋势，平均减少幅度为0.924 mm·a-1。

图2 5—9月松嫩平原潜在蒸散量的变化趋势Fig.2 The variation trend of potential evapotranspiration in each month of growing period

整体来看，研究区域潜在蒸散量的整体变化趋势为-0.509 mm·a-1，空间分布表现为由中部向南北逐渐增强。5月潜在蒸散量整体呈下降趋势，变化幅度为-0.257 mm·a-1，除黑河市少部分地区、长春市北部潜在蒸散量增加外，中部和西部地区潜在蒸散量明显下降。6月潜在蒸散量倾斜率与总平均潜在蒸散量倾斜率的空间变化相似，其正负高倾斜率中心均与总平均潜在蒸散量正负高倾斜率的分布相同，但6月潜在蒸散量的下降趋势稍弱，其变化幅度为-0.164 mm·a-1。7月不同区域表现出不同程度的变化趋势，中部和中东部地区下降趋势明显，下降幅度为-0.137 mm·a-1，南部和西北部部分地区为上升趋势，上升幅度为0.150 mm·a-1，研究区域整体变化趋势为-0.090 mm·a-1。8和9月潜在蒸散量变化趋势的空间分布情况较为相似，中部和中东部呈现下降趋势，其余表现为上升趋势，变化幅度分别为-0.030和0.020 mm·a-1。

3.2 松嫩平原水土平衡时空分布特征

经过实地调研了解到松嫩平原主要以旱地和水田为主，旱地和水田的不断增加以及气候因子的变化是松嫩平原水土平衡变化的主要原因，所以笔者主要探讨松嫩平原生长周期内耕地的水分盈亏情况。

3.2.1 松嫩平原耕地变化情况

利用ArcGIS软件随机选取松嫩平原500个样本点，提取解译的1980、1995、2005以及2015年耕地利用分类结果信息，结合Google Earth软件将分类结果与GlobeLand30提供的耕地信息进行对比，建立混淆矩阵进行精度评价。结果表明，解译耕地利用信息的平均精度为91.32%，总体Kappa系数为0.873，满足研究需求。

为了分析松嫩平原近35 a耕地利用变化情况，对不同时间节点的耕地利用类型进行统计，结果见表2。从表2可知，旱地面积占松嫩平原总面积的比例最大，2015年旱地占比为53.99%，水田占比为8.77%。从土地利用变化情况来看，1980—2015年旱地和水田面积分别增加9 760.99和15 166.43 km2。旱地面积的变化趋势表现为先增加后减少，主要是因为2005—2015年间松嫩平原大量旱地开发为水田所致。

表2 1980—2015年松嫩平原耕地利用类型变化Table 2 Changes of different land use types in Songnen Plain in different periods from 1980 to 2015

如图3所示，松嫩平原主要土地利用类型为耕地，主要分布于齐齐哈尔市东部、绥化市中部、哈尔滨市中西部、长春市以及西平市等地区，且旱地面积大于水田。从耕地的变化趋势来看，旱地和水田面积均呈现增加趋势，且水田变化趋势大于旱地。

图3 1980—2015年松嫩平原耕地利用类型分布Fig.3 The distribution map of cultivated land use in Songnen Plain from 1980 to 2015

3.2.2 不同月份水土平衡效应的空间分布特点

依据水分盈亏模型计算1961—2016年研究区各月平均水分盈亏指数的空间分布，其中作物需水量选取距离计算年份最近的耕地利用分布数据进行计算，不同干旱等级水分盈亏指数的空间分布见图4。

图4 5—9月松嫩平原水分盈亏情况的空间分布Fig.4 The spatial distribution of water deficit in each month of growing period

从图4可知，5月水分亏缺严重的区域主要分 布在松嫩平原中西部，正常和轻度缺水区域主要分布在研究区域的南北两端和东部地区;6月中度缺水以上的区域明显增加，分布范围向东部和北部扩散，但严重缺水和重度缺水的区域面积减少;受降雨量影响，7—8月正常和轻度缺水的区域面积逐渐增加，中度缺水以上的区域零散分布，这主要是由于土壤保水性和透水性不同而引起的地表径流量区域差异所致;9月研究区域内中度缺水的区域面积明显增加，主要分布于中部偏北地区，这是由于松嫩平原雨季结束，耕地水分亏缺情况增加。

从表3可知，5—9月正常和轻度缺水等级的区域面积占比表现为减少—增加—减少的趋势，而中度缺水区域面积占比的变化趋势刚好相反。这与研究区域内的降雨量变化趋势相符，即降雨量主要分布于7—8月，5—6月由于降雨量减少和季风影响，松嫩平原的旱情加剧。重度缺水和严重缺水区域的面积占比大致表现为逐渐减少，这主要是因为随着雨季的到来，土壤含水量逐渐增加所致。总体来说，生长周期内容易出现旱情的时间主要为5—6月。

表3 5—9月松嫩平原不同水分亏缺等级区域的面积占比Table 3 The water deficit area in each month of growing date %

3.2.3 基于作物亏缺指数的干旱频率空间分布

依据ArcGIS平台构建5 km×5 km的网格，提取近55 a生长周期内的干旱发生频率。用反距离加权插值法对干旱频率进行空间插值，得到松嫩平原生育期内近55 a各月份中度缺水以上的干旱发生频率空间分布(图5)。

图5 5—9月松嫩平原中度缺水以上的干旱发生频率空间分布Fig.5 The spatial distribution of dry and drought frequency in each month of growing period

5月中度缺水以上的干旱发生频率较高，高值区域主要分布于研究区中部和西部，发生频率均大于30%，发生频率大于50%的区域主要分布在白城市、松原市北部和大庆市东南部分;6月中度缺水以上的干旱发生频率分布与5月基本一致，但发生频率相比5月略有降低;7月大部分地区的干旱发生频率在20%以下，高值区主要分布在东南部的部分区域;8月中度缺水以上的干旱发生频率较低，大部分地区在10%以下，北部少量区域的干旱发生频率在＞20%～30%之间;9月干旱发生频率的空间差异较大，南部和中部区域的干旱发生频率大多低于20%，北部区域的干旱发生频率为＞20%～30%。可见松嫩平原西部易发生伏旱和秋旱，这与李取生等［26］的研究结论一致。辽宁西部是伏旱和秋旱的易发区，辽宁东南部和吉林东南部则是夏涝和秋涝的高发区。

3.2.4 水分盈亏指数变化趋势及周期变化特征

计算不同时期松嫩平原的平均水分盈亏指数，在此基础上分析得到水分盈亏指数的年际变化趋势(图6)。

图6 1961—2016年松嫩平原平均水分盈亏指数的年际变化Fig.6 The annual variation of average crop water deficit index in Songnen Plain from 1961 to 2016

1961—2016年5—7月松嫩平原水分盈亏指数呈下降趋势，其变化速率分别为-0.038·(10 a)-1、-0.035·(10 a)-1和-0.007·(10 a)-1;8—9月水分盈亏指数呈上升趋势，其变化速率分别为-0.01·(10 a)-1和-0.009·(10 a)-1。

从年际变化来看，不同月份的水分盈亏指数表现为波动变化状态。5月水分盈亏指数呈周期性波动，其峰值出现在1963、1980和2001年，从5 a滑动平均来看，水分盈亏指数波动幅度逐渐增大，表明研究区域内干旱发生频率降低，但干旱强度并未减弱。6月水分盈亏指数逐渐降低，峰值出现在1993和2001年，从5 a滑动平均来看，干旱发生频率变化趋势与5月相似。7月水分盈亏指数呈周期性变化，峰值出现在1968和2001年，从5 a滑动平均来看，7月水分盈亏指数波动幅度增大，说明研究区域水分盈亏情况出现两极分化。8—9月水分盈亏指数呈上升趋势，峰值出现在1980和2009年，从5 a滑动平均来看，水分盈亏指数平均值的波动幅度逐渐减弱，表明松嫩平原8—9月的旱情比较稳定，但发生概率变大。已有研究表明，松嫩平原2001年春夏季发生严重旱灾，2009年出现春旱，2008年7—9月频发大雨、暴雨，造成严重的涝灾［27-28］。笔者的研究结论与其基本相符。

采用Morlet小波函数研究1961—2016年松嫩平原生长周期内水分盈亏情况的周期变化特征，结果见图7。

图7 5—9月松嫩平原水分盈亏指数的Morlet小波系数等值线Fig.7 Themorlet wavelet coefficient isoline of water profit and loss index in each month of growing period in Songnen Plain

从图7可以看出，松嫩平原区域内各月份水分盈亏情况存在着明显的周期变化特征。在15 a以上时间尺度，区域内各月份水分盈亏情况出现明显的周期波动，其周期波动的特征具有很强的连续性，小波系数符号表现为正—负—正—负等周期振荡过程，表明1961—2016年水分盈亏情况表现出偏高—偏低—偏高—偏低等循环交替的周期振荡过程;从15 a左右时间尺度来看，小波系数负等值线没有全部闭合，反映在2016年以后一段时间内，5、6和9月为水分盈亏的偏低时段，7—8月处于水分盈亏的偏高时段，这为松嫩平原土地利用水资源的合理分配提供了有效依据;从10 a左右的时间尺度来看，区域内各月份周期变化的结构性和连续性表现较为明显，均表现为水分盈亏偏高期与偏低期交替出现，5、7、8和9月为水分盈亏的偏高期，6月为水分盈亏的偏低期;在5 a以内时间尺度上，全区域各月份平均水分盈亏情况的周期变化特征极其复杂，正负交替振荡情况更频繁，水分盈亏情况存在更多的偏低期与偏高期。

4 讨论

松嫩平原是我国重要的粮食生产基地，粮食生产安全直接关系到区域经济的可持续性发展，而水资源是制约粮食生产的重要决定因素。以往的研究通常以年为时间尺度，分析区域内年水土资源平衡状态，通常假设土壤含水量在年初和年末一致，即土壤含水量是以年为单位进行周期循环，从而规避土壤含水量变化的影响，而实际上受年降水强度的影响，土壤含水量年内变化明显。因此，该研究充分考虑年内土壤含水量变化对于水土平衡的影响，不仅可以更加准确地分析区域内水土资源平衡态势的变化情况，而且可以有效研究年内不同月份水土资源的平衡效应。

该研究发现，松嫩平原不同月份的耕地水土平衡情况表现出明显的差异性。其中，5—6月水分亏缺程度比较严重，而7—9月水分亏缺情况较为严重。这与已有研究结论［31］相符。此外，从水分亏缺的时间变化特征来看，5—6月水分亏缺状况呈现逐渐降低趋势，而8—9月水分亏缺状况呈现逐渐增加趋势，该结论可为研究区域的引水灌溉计划提供参考。

松嫩平原地处于湿润季风气候与大陆性干旱气候的过渡带，蒸发量较大，为气候敏感区域［32］，湿地资源和林地资源对于区域水土资源平衡具有重要影响［33］。随着建设农业生态文明理念的提出，山水林田湖草农业生态系统的概念不断被强调，从系统的概念研究土地利用变化及其影响具有重要的意义。因此，只有将区域内土地利用细化，从山水林田湖草系统整体角度出发，才能更为准确地研究区域水土资源平衡效应的变化，这也是未来研究的重点。

5 结论

(1)1961—2016年松嫩平原生长周期内的多年平均潜在蒸散量在432.8～605.8 mm之间，研究区潜在蒸散量的变化趋势表现为中部向南北逐渐增强，总体变化趋势为-0.51 mm·a-1，5—9月的变化趋势分别为-0.27、-0.16、-0.09、-0.03 和 0.02 mm·a-1。

(2)1961—2016年松嫩平原水分盈亏程度严重的月份为5—6月，中西部地区缺水程度比较严重;9月水分盈亏程度较为严重，缺水区域主要分布在松嫩平原北部偏东的少数地区。缺水发生频率的分布情况与水分盈亏情况分布较为一致，其高值月份为5—6月，其次为9月，分布区域一致。

(3)从年际变化规律来看，15 a时间尺度上各月份均呈现比较连续的周期性变化，其中5—7月的水分盈亏情况处于波动下降趋势，但水分盈亏指数波动幅度增大，表明旱情的发生强度有所增大，发生概率有所下降;而8—9月水分盈亏情况呈现波动上升趋势，但水分盈亏指数波动幅度减弱，表明8—9月出现旱情的概率变大，旱情强度略微减弱。