柴达木盆地湖泊天然生态水文特性辨识

2023-02-04陈一迪诸葛亦斯石岳峰张馨予

水资源保护 2023年1期

陈一迪，诸葛亦斯，石岳峰，杜强，张馨予，聂睿

(1.三峡大学水利与环境学院，湖北宜昌 443002；2.中国水利水电科学研究院水生态环境研究所，北京 100038；3.黄河生态环境科学研究所，河南郑州 450000)

位于我国西北部的柴达木盆地是青藏高原的重要组成部分，是我国湖泊分布最密集的地区之一，也是全球盐湖最集中的地区[1-2]，属于全球气候变化敏感、生态环境脆弱的地区，区域内的气象水文变化直接影响青藏高原的水文循环系统。近年来由于受气候变化和人类活动等诸多因素的影响，盆地内的湖泊濒临绝境、甚至消亡[3]。位于柴达木盆地中部的东、西台吉乃尔湖为那棱格勒河的尾闾湖泊，湖区内富含硼、钾、锂、镁等资源[4]。盐矿企业自2003年陆续入驻湖区进行开发，为保护企业生产安全，分别在东、西台吉乃尔湖的西侧和东侧筑起拦水坝，使上游来水在两湖之间蓄积形成鸭湖滞洪区，出现了“水上雅丹”的景观。生物多样性随着鸭湖面积的不断增加而增加，趋于成为一个稳定的生态系统。由于鸭湖水位上升过快会造成洪水，因此多余的湖水沿西台吉乃尔湖东侧的苦水沟流向西台吉乃尔湖西北侧的一里坪，致使常年处于干涸状态的干盐湖开始出现湖表卤水。同时，上游来水减少导致东、西台吉乃尔湖持续萎缩，逐渐向干盐湖转化[5-7]。盐湖资源的开发逐渐改变了天然生态系统，而资源开发与生态保护之间的平衡关系也是目前研究的热点。在人为活动逐渐强烈的情况下，东、西台吉乃尔湖区形成了人为次生湖泊生态系统，人为生态系统能否替代天然生态系统是目前亟待研究的问题，为更好地解决该问题，需要准确识别湖泊的天然状态。

目前对于湖泊天然状态的研究集中在基于长序列观测数据下对湖泊面积、水位等动态变化特征及影响因素的分析。柴达木盆地气候恶劣、地形复杂、基础设施相对落后，观测站点稀疏，缺少长序列的实测数据，是典型的缺资料地区[8]。补充缺测数据常用的手段是构建流域模型[9-11]。目前，也有学者针对缺资料地区高原尾闾湖泊进行了一些研究，如Wang等[12]对缺资料地区尾闾湖泊流域基于遥感流量估算方法模拟湖泊水量；Zhang等[13]以尕海为研究对象，基于遥感影像及卫星测高数据构建了湖泊面积—水位曲线；Wang等[14]以艾比湖流域为研究对象，分析河流演变对湖泊水位变化的响应。这些研究对了解缺资料地区湖泊动态变化研究方法具有重要的参考价值，但是缺资料地区湖泊天然状态的辨识研究成果较少。

鉴于此，本文以柴达木盆地的东、西台吉乃尔湖区为研究区域，利用未开发时期的1990—2002年全年Landsat遥感影像数据(1996年、1997年影像数据缺失)，分析湖泊面积的动态变化及水量平衡，辨识缺资料地区的湖泊在不同保证率下的天然状态，以期为解决柴达木盆地内盐湖资源开发布局与生态空间矛盾提供参考。

1 研究区概况

青海省海西州东、西台吉乃尔湖为柴达木盆地内最大河流——那棱格勒河的尾闾湖泊(图1)，位于东经93°49′～94°01′、北纬37°24′～37°36′，平均海拔2 683 m，两湖近似三角形，呈西北-南东向分布。湖区及周边人口极少，地质构造复杂，基础设施相对落后，是典型的缺资料地区。根据距离湖区最近的小灶火气象站统计资料，该地区年均气温4.4℃，平均风速3.6 m/s，平均降水量29.8 mm，蒸发量2 600 mm，为典型的高原大陆性干旱气候。

(a) 1990年

(b) 2020年图1 东、西台吉乃尔湖不同时期遥感影像Fig. 1 Remotesensing images of East-Taijiner Lake and West-Taijiner Lake in different periods

盐矿企业为保证生产安全，在东、西台吉乃尔湖的西侧和东侧自建拦水坝，截断了东台吉乃尔河与两湖的直接水力联系，导致上游来水在两湖之间蓄水形成鸭湖。由于季节性洪水及采矿活动使湖泊丧失了储水调节功能，迅速向干盐湖转化。区域内水文特征发生长期不可逆转的不利转变，同时鸭湖面积逐年增加，水位逐渐升高，对企业生产设施带来防洪隐患，严重影响东、西台吉乃尔湖的周边生态环境。

2 研究数据与研究方法

2.1 数据来源

选择Landsat卫星系列影像为主要数据来源，为提高湖泊水体数据提取的准确性，优先选择云量不超过10%的遥感影像。由于1990年以前遥感数据严重不足，选取东、西台吉乃尔湖区域1990—2002年(除1996年和1997年)全年遥感影像，计算每年的湖泊平均面积。

选取距离东、西台吉乃尔湖最近的小灶火气象站实测数据进行气候要素变化分析，包括降水量、蒸发量、风速及气温的逐日数据，资料系列为1990—2002年。

那棱格勒河水文站于1958年设立，1964年撤销，此后没有监测资料，水文资料严重不足，因此选取附近的奈金河纳赤台站、格尔木河格尔木站、香日德河香日德站作为参证站在实测数据的基础上进行水文数据的插补延长。

2.2 数据处理与分析

利用ENVI 5.3软件对遥感影像依次进行假彩色合成处理、几何校正、辐射定标和大气校正[15]。选用改进的归一化差异水体指数(modified normalized difference water index， MNDWI)对水体信息进行提取，该方法被证实可以更好地应用于干旱区水体信息提取[16-17]，其表达式为

(1)

式中：IMND为MNDWI值；ρgr、ρMIR分别为绿光波段和中红外波段的反射率。

为分析水量平衡变化原因，引入累积距平值Si[20]检验实测数据系列平均值发生突变的情况，表达式为

(2)

当实测值大于长系列平均值时斜率为正，反之为负，持续的正负斜率用来鉴别系列平均值的中间突变。

2.3 水量平衡计算

东、西台吉乃尔湖属于干旱区的半封闭式湖泊，湖泊补给量主要来自降水、地表水入湖和地下水补给，排泄量主要是湖面的蒸发损失和地下水渗流[18-19]。由于研究区资料有限，将地下水的补给与渗流量合并为地下水交换量，其量值为库容差与其他水平衡项的差，正项表示补给量多于渗流量，反之则为渗流量多于补给量。水量平衡方程为

V=P+R+W-E

(3)

式中：V为湖泊水量变化量；P为湖面降水量；R为入湖径流量；W为地下水交换量；E为湖面蒸发量。

由于东、西台吉乃尔湖缺乏湖面降水量的气象观测资料，选用小灶火气象站的降水量观测值和遥感影像提取的湖面面积计算湖面降水量：

P=P小A

(4)

式中：P小为小灶火气象站的降水量观测值，m；A为湖面面积，m3。

由于湖面未布设蒸发皿或蒸发池，因此利用小灶火气象站的小型蒸发皿(20 m2)观测的蒸发量来计算湖面实际蒸发量，满足世界气象组织观测湖泊水面蒸发量的要求[21]，则东、西台吉乃尔湖的湖面蒸发量E计算公式为

E=E小F

(5)

其中

F=1-0.011(1-r1)-0.379B′e1.001

式中：E小为小灶火气象站的小型蒸发皿观测蒸发量，mm；F为卤水蒸发折算系数，反映水体含盐量对水面蒸发的影响；r1为空气相对湿度；B′e为含盐量，用波美度表示，1波美度相当于1 L水含10 g氯化钠。

2.4 湖泊天然状态识别方法

通过对湖泊水量的计算得出湖泊上游来水补给的变化和地表水地下水的交换量，结合湖水的补给路径和湖区地形，可以分析东、西台吉乃尔湖在天然状态下的湖泊面积、湖容和两湖的空间布局关系。将1990—2002年提取出的湖泊面积分为19组，计算各组出现的次数及频率，将各组频率相加，其累积频率就是湖泊面积保证率，即湖泊面积在多年期间可以得到满足的程度，可以绘制湖泊面积保证率曲线。选取不同的保证率可以确定天然状态时的湖泊面积，根据水量平衡方程和长序列的气象水文数据计算得出湖容，结合补水路径和地形进一步确定两湖的空间布局关系。东、西台吉乃尔湖的盐湖资源开发可以根据该保证率确定开发程度，从而保证天然湖泊生态系统的稳定和湖泊水量的平衡。

3 结果与分析

3.1 气候变化特征

选取小灶火气象站1990—2002年降水量、气温、蒸发和风速4个气象因子绘制逐月变化曲线，如图2所示。从图2可见，降水量集中在4—9月，最大降水量出现在7—9月。降水年际变化倾向率为每10年11 mm，呈下降趋势，1997年降水量达到最大值37 mm，在2000年产生突变，以2000年为分界线，1990—2000年年均降水量为27.33 mm，2000—2002年年均降水量为11.67 mm。蒸发量为降水量的80倍，集中在4—9月，6—8月达到峰值。以26 mm/a的速率呈上升趋势，1990—1997年蒸发量大部分低于均值，1992年蒸发量最小，为2 319.5 mm，2000年后蒸发量逐渐增加，2002年达到最大值2 892.1 mm。气温在每年3、4月左右达0℃以上，11月下降到0℃以下，7—8月气温达到峰值。2000年7月气温最高，为20.18℃，1995年1月气温最低，为-13.36℃，1990—1997年气温变化幅度不明显，1997年以后气温以每10 a 0.6℃的速度上升。风速的峰值集中在5—7月，以西-西北风和西北风为主，共占研究时段内的51.64%。1990—2002年风速年均下降0.16 m/s，1991年年均风速最大，为39.72 m/s，1998年风速最小，为32.65 m/s。

(a) 降水量和蒸发量逐月变化

(b) 气温和风速逐月变化

(d) 气温和风速年际变化图2 小灶火气象站气象因素变化Fig.2 Variation ofmeteorological factors at Xiaozaohuo Meteorological Station

3.2 水量平衡计算结果

由于缺少实测湖盆数据，通过构建湖盆数字高程模型得到面积-库容曲线，根据遥感影像提取出的湖泊面积确定水量变化量V。具体步骤为：通过DEM信息插值计算遥感影像提取出水体边界线高程值，对不同时间序列的湖泊边界图层进行叠加，生成湖泊等水位图，再通过ArcGIS软件构建湖泊的不规则三角网(triangulated irregular network，TIN)模型，计算不同湖泊面积对应的库容[23]。东、西台吉乃尔湖的水量平衡方程内各项变化情况见图3。

(a) 蒸发量

(b) 降水量

(d) 湖容变量

(e) 地下水交换量图3 水量平衡各项变化情况Fig.3 Changesin items of water balance

水量平衡的计算结果表明，东台吉乃尔湖的地下地表水交换量平均值为0.05亿m3，补给量最大值为0.53亿m3，渗流量最大值为0.78亿m3；西台吉乃尔湖的地下地表水交换量平均值为0.01亿m3，补给量最大值为0.34亿m3，最大渗流量为0.32亿m3，地下水交换量呈现不确定性。两湖的库容差呈现正常的波动状态，2000年以后变化增大，地下水交换更加剧烈，是由于2000年后两湖的面积都发生了突变，水平衡各项需要进行调整从而确定新的平衡点。

3.3 湖泊面积动态变化

图4和图5分别为东、西台吉乃尔湖面积年内和年际变化。由图4图5可见，两湖总面积在140 km2以上，1990—2002年两个湖泊总面积变化具有波动性特征，不同时段的年内波动振幅与均值有显著差异。1990—2002年，东台吉乃尔湖在1990年2月面积最小，2000年后面积呈扩大趋势，最大值为261.70 km2；西台吉乃尔湖呈萎缩趋势，面积最大值出现在1992年7月，2001年1—3月湖泊干涸，除此之外面积最小值出现在2000年2月。东台吉乃尔湖面积变化可以分为3个阶段：1990—1992年湖泊面积以33 km2/a的速率扩张；1993—1999年湖泊面积逐渐萎缩，每年缩小6 km2；2000—2002年湖泊面积变化较剧烈，湖泊面积年均扩张2.19 km2。西台吉乃尔湖在1992年急剧扩张，此后面积由118.41 km2逐渐下降到25.55 km2，即使1994年和2001年有短暂的扩张，也没有改变湖泊总体萎缩趋势。两湖总面积以1.3 km2/a的速率萎缩。

(a) 湖泊面积逐月变化

(b) 湖泊面积月均值变化图4 东、西台吉乃尔湖面积年内变化Fig.4 Annual changes of area of East-Taijiner Lake and West-Taijiner Lake

(a) 东台吉乃尔湖

(b) 西台吉乃尔湖图5 东、西台吉乃尔湖面积年际变化Fig.5 Interannual changes of area of East-Taijiner Lake and West-Taijiner Lake

3.4 湖泊边界动态变化

那棱格勒河发源于昆仑山北部，位于柴达木盆地西南，由西南流向东北，多年平均流量32.92 m3/s，出山口后有60%～80%的水量渗漏，少量河水到下游与台吉乃尔河汇流，最终在湖盆海拔最低处形成东、西台吉乃尔湖[24]。选择1990年、1992年、1994年和2000年中面积最大月份的湖泊边界，分析东、西台吉乃尔湖的空间变化情况。东、西台吉乃尔湖是发育在柴达木盆地中部的构造断陷湖泊，湖区及周边地质构造复杂，北部为雅丹地貌群，南部为河流区，东、西侧为盐滩-荒漠区，湖盆边缘平均海拔2 700 m，主要是砂质泥岩和砂质黏土，湖盆内部平均海拔2 680 m，以黏土和盐类沉积为主[25]。

图6为湖泊边界变化情况，由图6可见，东台吉乃尔湖西南侧地势平坦，湖岸线变化明显，其东部为近50 km2的沙质干盐滩，厚15～20 m，坡度较陡，湖岸线变化不明显；西台吉乃尔湖西北高东南低，因此湖岸线变化主要发生在西北侧，上游来水在南侧向湖泊内补给，湖岸线变化不明显。当两湖水量充足时均呈现为三角形，当水量减少时，东台吉乃尔湖依然以近似三角形的形状向内缩小，由于东南侧海拔较高，因此面积萎缩时湖岸线变化较明显，而西台吉乃尔湖则逐渐变为椭圆形。由此可见，湖泊水量的变化不仅影响湖泊面积，还影响湖泊的形态分布。

图6 湖泊边界变化情况Fig.6 Change of lake boundary

4 讨论

4.1 水量平衡

1990—2002年东、西台吉乃尔湖的月水量平衡计算结果(图7)表明，两个湖泊库容差变化不一，没有持续萎缩，处于正常波动状态。2000年后两湖的水量变化剧烈，采用累积距平法进行分析，结果表明，东台吉乃尔湖的湖面蒸发量、降水量、径流量和库容的平均值分别为0.208亿m3、0.003亿m3、0.159亿m3和1.496亿m3。其中，蒸发量、径流量和库容的变化趋势几乎一致，1990年开始呈下降趋势，在2000年前后达到最小值，此后逐渐上升，而降水量的变化趋势完全相反，说明东台吉乃尔湖库容发生变化的主要原因是蒸发量和径流量。西台吉乃尔湖的湖面蒸发量、降水量、径流量和库容的平均值分别为0.099亿m3、0.002亿m3、0.085亿m3和0.360亿m3。库容变化与降水量和蒸发量的变化趋势相似，1992年开始急剧上升，1999年前后达到最大值，此后逐渐下降，而径流量较稳定，呈先下降后上升的趋势，说明西台吉乃尔湖的水量变化主要取决于降水量和蒸发量。由于2003年以前湖区处于天然状态，并且径流量、降水量、蒸发量的突变几乎在同一时期发生，表明引起湖泊水量的变化主要原因是气候变化。

(b) 西台吉乃尔湖图7 水平衡各项累积变化Fig.7 Cumulative changes in items of water balance

4.2 湖泊天然状态识别

由于东、西台吉乃尔湖区域年均降水量仅为30 mm左右，蒸发量是降水量的80倍，所以降水量的影响可以忽略不计。湖泊的径流量年内变化出现两次峰值：第一次峰值在4月，由于气温上升使冰雪融水量增加，第二次在7月，与汛期降水量导致地表产流量增加有关。相关研究表明，由于冰雪融水参加径流过程会导致年内径流过程有两次峰值，这也是高寒地区水文过程的特征之一[26-27]。同时柴达木盆地近年增温速度明显高于青藏高原整体水平[27]，那棱格勒河流域位于昆仑山脉附近，周围有冰川分布，气温的变化对流域积雪的冻融过程贡献较大，因此两湖泊面积从3月开始随着气温的增加逐渐扩大，至7月达到峰值。

相关性分析结果表明，两湖泊年内面积变化与气温和径流量均呈正相关关系，并且相关性系数均通过了置信度为0.01的双侧检验。选择气温和径流量作为自变量，湖泊面积作为因变量，分别对东、西台吉乃尔湖湖泊面积年内变化建立多元线性回归方程：

A1=142.243+48.35R+0.246T

(6)

A2=64.049+28.695R+0.736T

(7)

式中：A1、A2分别为东、西台吉乃尔湖湖泊面积，km2；T为气温，℃。

拟合面积的峰值出现在8月，主要是两湖径流量在8月份达到最大值所致。根据回归系数t检验结果，径流量对年内变化的贡献率较大；经方差分析，两湖的线性回归模型总体达到极显著水平(p<0.01)；Durbin-Waton检验值表明模型拟合效果较好，预测值的方差独立。

根据年蒸发量等气象因素与湖泊面积年际变化相关性分析结果，选择相关性较高的因子对两湖年际面积变化建立多元线性回归方程：

A1=717.179-0.146E-2.209T-5.920w
+41.220R

(8)

A2=542.801-0.108E-9.606T-4.561w
+6.651R

(9)

式中w为风速，m/s。

经分析，方程拟合结果较符合实际，西台吉乃尔湖2000年的拟合值大于实际值是由于模型中未考虑降水因素，而2000年年均降水量骤减至8.7 mm，对湖泊面积产生了一定的影响。经方差分析，模型总体达到极显著水平(p<0.01)；Durbin-Waton检验值表明模型拟合效果较好，预测值的方差独立。

图8为东、西台吉乃尔湖面积累积频率分布。可见，东台吉乃尔湖面积为180～220 km2，西台吉乃尔湖面积为130～140 km2。由于地形的原因，天然状态下两湖的位置保持在湖区的海拔最低处，当人为筑起拦水坝阻隔上游来水时，会在人为规定的空间范围内海拔最低点处蓄水形成新湖泊，即鸭湖，其平均海拔为2 685 m。因此，天然状态时的湖泊面积和形态随水量的变化而改变。表1为选取25%、50%、75%和90%的湖泊面积保证率情况下辨识的东、西台吉乃尔湖湖泊面积、湖容及空间形态的结果。

表1 不同湖泊面积保证率下的湖泊天然状态Table 1 Lake States under Different Guarantee Rates

图8 东、西台吉乃尔湖面积累积频率分布Fig.8 Accumulative frequency distribution of area of East-Taijiner Lake and West-Taijiner Lake

5 结论

a.在湖泊面积变化方面，东台吉乃尔湖面积以2.19 km2/a的速率呈增加趋势，面积为180～220 km2，西台吉乃尔湖面积以3.50 km2/a的速率缩小，面积为130～140 km2。两湖面积在7月出现峰值，4月为次峰值；选取相关性较高的气象水文因子分别对东、西台吉乃尔湖面积的年际和年内变化建立了线性回归方程，拟合度较高。

b.在湖泊空间布局方面，东、西台吉乃尔湖呈西北-南东向分布，湖泊水量充足时两湖均近似为三角形，水量减少时，东台吉乃尔湖以原有的湖泊形态向内缩小，而西台吉乃尔湖逐渐萎缩成椭圆形。由于湖区周围地形和上游来水补给路径的原因，两湖边界变化较明显，主要在东台吉乃尔湖的西南侧和西台吉乃尔湖的西北侧发生变化。

c.研究时间范围内湖泊水量维持在一个有正有负的相对平衡的状态，结合由长序列的湖泊面积变化绘制出的湖泊面积保证率曲线，以25%、50%、75%和90%的保证率为例，识别天然状态的湖泊面积、湖容及空间形态。该保证率曲线可为盐湖资源开发过程中保护天然湖泊生态系统提供一定的参考。