杨 凡,马 巍,陈 欣,王云飞,王军亮

(1.中国水利水电科学研究院 水生态环境研究所,北京 100038; 2.云南省水利水电勘测设计研究院,昆明 650021)

0 引 言

pH值作为水体最重要的理化参数之一,在水环境质量评价和环境保护中都具有重要的指示性作用。湖泊pH值常年偏高不仅会影响水生动植物的栖息生长环境,同时也有助于湖泊蓝藻水华的形成,使得水质进一步恶化。正常情况下,河流、湖泊、水库等天然水体的pH值多在6～9范围内波动。近年来滇池等高原湖库频繁出现pH值异常升高(pH值>9)的现象,查明pH值异常升高的原因并解决其异常升高问题正成为滇池环境保护工作中不可缺少的环节。

目前,国内学者对于河流、湖泊、水库pH值升高影响因素的研究方法多为简单的时空特征分析、皮尔逊相关分析或单变量线性回归[1-4],并未直接有力地指出导致水体pH值异常升高的主要因素。相关分析通常只能用来判断变量之间共变程度的大小和方向,分析过程中不区分自变量和因变量,其结果并不代表变量之间的因果关系[5],也不能说明变量之间影响的强弱,且简单相关分析并不适用于非正态数据,不同河流的水体流入湖泊以后水质数据不一定符合正态分布,因此,在进行湖泊水质数据的相关分析以前需要进行正态性检验。为了综合分析滇池pH值异常升高的影响因素,并比较相关分析和多元线性回归在湖泊多因素综合分析的结果上是否存在差异,本文以相关分析以及多元线性回归结合的方法来对滇池pH值异常升高的原因进行系统性研究[6],识别其主要影响因素,分析其驱动机理,可为高原湖泊及水库pH值异常升高成因识别与科学施治提供正确方向。

1 研究区域概况

滇池位于云南省昆明市西南面,隶属金沙江流域一级支流普渡河水系,是云贵高原面积最大的天然淡水湖,为我国四大断陷构造湖泊之一。滇池分为内湖(草海)和外湖(外海),其中草海位于滇池北部,湖面面积约11.2 km2;外海是滇池的主体部分,湖面面积约298 km2,占滇池总面积的96.7%。

受牛栏江—草海应急补水工程、新老运粮河入湖河口前置库水体净化生态工程、环湖截污工程等影响,滇池草海湖区水质相对较好[7-8],而外海日益成为滇池水质污染相对较重、最需要治理和关注的部分。滇池外海常规水质监测数据评价结果显示(其水质监测点分布见图1),2016—2019年滇池外海年均综合水质类别依次为Ⅴ类、劣Ⅴ类、Ⅳ类、劣Ⅴ类,年际波动变化特征突出,水质类别控制性指标主要有化学需氧量(COD)、总氮(TN)、总磷(TP)、pH值等。

图1 滇池常规水质监测站点布置示意图Fig.1 Layout of conventional water quality monitoring stations in Dianchi Lake

2 数据来源与研究方法

2.1 数据来源

本文采用的水质数据主要来自昆明市环境监测中心,监测时间为2016年1月至2021年12月,常规水质监测点包括晖湾中、罗家营、观音山东、观音山中、观音山西、海口西、滇池南、白鱼口共8个国控监测点,监测频次为1次/月。水质指标主要包括pH值、TP、TN、COD、溶解氧(DO)、氨氮(NH3-N)、叶绿素a 7个常见理化因子。由于2020—2021年期间获取的pH值均为整数,无法满足本次研究的精度要求,故本文只使用2016—2019年期间的水质监测数据。

2.2 研究方法

2.2.1 斯皮尔曼相关分析

斯皮尔曼相关系数,全称斯皮尔曼秩相关系数,其计算内容为等级数之差,用于分析变量之间单调关系的强度。与皮尔逊相关分析相比,该方法的适用范围更加广泛,对数据本身没有严格的要求,当出现以下情况时,使用斯皮尔曼相关分析更加恰当:①变量之间为非线性关系,或变量本身不是正态分布;②有一个或多个变量为序数类型;③变量有明显的异常值。斯皮尔曼相关系数r具体公式为

(1)

式中:Ri和Si分别为2个变量第i个值的秩;R、S分别为2个变量的平均值。r的取值范围为[-1,1],当r>0时,变量之间呈正相关关系,r<0时,变量之间呈负相关关系,r的绝对值越大说明变量之间的相关性越强。

斯皮尔曼相关系数的显著性检验方法为t检验,其具体公式为

(2)

式中n为样本数据量,n-2为自由度。

2.2.2 多元线性回归

相关分析通常用来判断变量之间关联程度的大小和方向,分析过程中不区分自变量和因变量,其结果并不代表变量之间的因果关系,也不能说明变量之间影响的强弱。而多元线性回归能够根据模型回归的结果分析各个自变量对因变量影响的显著性。本文采用多元线性回归方法对滇池水质数据进行分析,以便更深入地探究滇池pH值异常升高的主要因素。具体公式为

y=β0+β1x1+β2x2+…+βn-1xn-1+ε。

(3)

式中:y为因变量;x1,x2,…,xn-1为自变量;β0,β1,β2,…,βn-1为模型当中的n个未知参数;ε为均值为0且不可观测的随机误差项。

本文的相关分析与多元线性回归分别采用Origin 2021与IBM SPSS Statistics 26统计软件完成。

3 结果与分析

3.1 水质数据特征分析

滇池外海2016—2019年各监测点位pH值均在7.5以上,且8个采样点pH值的中位数和均值均处于8.5以上,有24.59%的pH值超过了《地表水环境质量标准》(GB 3838—2002)规定的限值(6～9),处于异常升高的状态。滇池外海各站点pH值异常升高程度详见图2。

图2 滇池外海各监测点位pH值统计情况Fig.2 Statistics of pH values at monitoring points in the outer waters of Dianchi Lake

据2016—2019年期间滇池外海各监测点的pH值实测数据资料统计分析,出现pH值>9的情况多发生在夏(6—8月份)、秋(9—11月份)、冬(12、1、2月份)3个季节(见图3)。为了研究外海pH值异常升高的原因,考虑到不同季节自然环境条件的差异,将夏、秋、冬3个季节的水质数据进行单独分析处理。

图3 2016—2019年滇池外海pH值异常升高现象年内分布情况Fig.3 Monthly distribution of anomalously high pH values in the outer waters of Dianchi Lake from 2016 to 2019

2016—2019年期间外海夏、秋、冬季pH值正态性检验结果如表1及图4所示,3个季节pH值假设检验显著性均<0.05,原数据不属于正态分布,故均采用斯皮尔曼相关分析来计算pH值与其他水质数据之间的相关性。

表1 滇池外海夏、秋、冬季pH值正态性假设检验摘要Table 1 Summary of hypothesis testing for pH normality in summer, autumn and winter in the outer waters of Dianchi Lake

图4 滇池外海夏、秋、冬季pH值正态性检验结果Fig.4 Results of pH normality test in summer,autumn and winter in the outer waters of Dianchi Lake

为了得到更精准、更稳定的回归结果,需要对序列的平稳性进行检验。选用扩展的迪克富勒(Augmented Dickey-Fuller,ADF)检验验证所选理化因子是否具备单位根,若不存在单位根,则说明序列稳定,验证结果见表2。根据检验结果可知TN和TP的单位根存在的可能性较高,不通过显著性为10%的检验;在经过一阶差分以后,每一个序列的单位根存在的可能性都极低,均通过了显著性为1%的显著性检验。因此,在进行回归以前应对模型的自变量和因变量进行一阶差分。

表2 2016—2019年相关水质参数单位根检验结果Table 2 Unit root test results for water quality parameters from 2016 to 2019

3.2 pH值异常升高影响因素分析

将外海2016—2019年期间pH值异常升高的季节水质数据进行斯皮尔曼相关分析后得到各水质指标之间的相关系数(见图5)。从相关分析结果可以发现,夏季pH值与DO呈非显著正相关,秋季pH值与叶绿素a呈非显著正相关,夏、秋两个季节均无水质指标与pH值呈显著性正相关;冬季pH值与DO和叶绿素a存在显著的正相关关系,相关系数分别为0.23、0.42。

图5 夏、秋、冬季水质指标相关分析结果Fig.5 Correlation analysis of water quality indicators in summer, autumn and winter

将各季节所选参数一阶差分以后进行多元线性回归,其中pH值为因变量,其余理化参数为自变量。

根据表3可知,夏季pH值回归模型调整后R2为0.348,F检验显著性为0.00,该模型对夏季pH值具有一定的解释能力,模型具有显著的统计学意义;秋季pH值回归模型调整后R2为0.215,F检验显著性为0.00,该模型对秋季pH值的解释能力较差,但模型具有显著的统计学意义,回归模型的因变量与自变量之间存在显著的线性关系;pH值冬季回归模型调整后R2为0.493,F检验显著性为0.00,该模型对滇池pH值异常升高的现象具有一定的解释能力,且模型具有显著的统计学意义。

表3 夏、秋、冬季回归模型摘要Table 3 Summary of regression models for summer,autumn and winter

由表4可知,所有系数的容差均>0.1,且方程膨胀因子(Variance Inflation Factor,VIF)<5,说明模型中不存在严重的共线性问题。夏季回归模型中TN和DO系数的显著性分别为0.000、0.002,它们对于回归模型均有较大的影响。DO是模型中系数>0且具有显著性的唯一的自变量,其未标准化系数和标准化系数分别为0.090和0.290,这意味着当DO浓度增长1个单位时pH值增长0.09,且消除了量纲和等级的差异后,模型中DO对pH值升高的贡献是所有变量中最大的[9]。

表4 夏、秋、冬季回归模型系数Table 4 Regression model coefficients for summer,autumn and winter

秋季回归模型中仅有TN和叶绿素a系数通过了显著性为0.05的检验,它们在回归模型中的显著性分别为0.000和0.024,但TN系数为负数,因此叶绿素a才是使得秋季外海pH值异常升高的主要因素。秋季回归模型中叶绿素a的未标准化系数和标准化系数分别为2.126和0.231,说明当秋季叶绿素a浓度上升一个单位时pH值增长2.126,且消除了量纲和等级的差异后,叶绿素a对pH值升高的贡献是所有变量中最大的。

冬季回归模型中氨氮和叶绿素a系数的显著性分别为0.001、0.000,它们对于回归模型均有较大的影响。叶绿素a是模型中系数大于0且具有显著性的唯一的自变量,其未标准化系数和标准化系数分别为4.046和0.444,这意味着当叶绿素a增长1个单位时pH值增长4.046,且消除了量纲和等级的差异后,叶绿素a对pH值升高的贡献远超其他变量[9]。

图6和图7为各季节回归模型的标准化残差P-P图及残差散点图。从结果来看,秋季模型回归结果较好,夏季和冬季的残差存在一定的偏态,但并不严重,对回归结果的稳定性不会造成太大的影响,3个模型整体效果都较好。

图6 夏、秋、冬季回归模型标准化残差P-P图Fig.6 P-P plots of standardized residuals of regression models for summer, autumn and winter

图7 夏、秋、冬季回归模型散点图Fig.7 Scatterplot of regression model for summer,autumn and winter

根据相关分析和多元线性回归的结果可知,滇池pH值异常升高主要的影响因素为叶绿素a和DO。但仅有冬季的DO和叶绿素a与pH值存在显著性正相关关系,也就是说仅从相关分析的结果来看,叶绿素a和DO只会在冬季对滇池pH值的升高有显著的贡献,这显然与观测到的数据不符合。而多元线性回归模型能够在消除理化因子等级和量纲的差距以后识别它们对因变量(pH值)的影响大小,因此,多元线性回归模型更适合用来分析湖泊中某一指标随多变量的变化情况。

从多元线性回归的结果来看,不同理化因子在不同季节对滇池水质pH值升高现象的贡献有所不同,说明各类理化因子对pH值的影响还与气象等外部因素有关。多因素的影响使得滇池pH值异常升高现象的主要推动因素存在一定季节性的差异,这样的差异导致了冬季叶绿素a对pH值升高的贡献明显大于其他季节。但这种差异只能说明在冬季叶绿素a对pH值升高的作用比其他理化因子更加明显,而夏季和秋季则没那么明显,并不表示叶绿素a在冬季比其他季节更能使pH值升高。

4 讨 论

河流、湖泊等天然水体的pH值主要受到水中游离CO2以及碳酸盐含量的控制。水体中的碳酸平衡方程为:

CO2(aq)+H2O=H2CO3(aq);KH=10-1.5;(4)

(5)

(6)

式中KH、K1、K2为化学平衡常数。

4.1 pH值与叶绿素a的响应关系

藻类作为水生态环境中的初级生产力,在生态系统中有着重要的作用,在一定程度上氮磷等营养盐对生态系统的良性循环有不可忽视的作用。而叶绿素a是所有藻类均含有的组成成分,其含量高低从侧面上反映了水体中藻类的数量。

藻类细胞(C106H181O45N16P)的光合作用与水体pH值之间存在着密不可分的关系。在藻类生长较为旺盛的水域,需要消耗更多的CO2来进行光合作用,此时藻类吸收CO2的速度远超大气中CO2溶于水体中的速度。当藻类光合作用需求不能被水环境内游离的CO2所满足时,只能用HCO3-替代CO2来进行光合作用,其具体反应式为:

(7)

(8)

将式(4)、式(5)、式(6)变形后得到:

lg(H2CO3)=lg(CO2(aq))=

lgKH+lgP(CO2)=-4.9;

(9)

lg(HCO3-)=lgK1+lg(H2CO3)+pH=

-11.3+pH;

(10)

lg(CO32-)=lgK1+lgK2+lg(H2CO3)+2pH=

-21.6+2pH。

(11)

图8 不同pH值时水体中碳酸的存在情况Fig.8 Presence of carbonic acid in the water column of different pH values

图9为滇池外海2016—2019年期间pH值与叶绿素a均值的分布情况。不难发现图9(a)中pH值较高的区域正好与图9(b)中叶绿素a浓度较高的区域相吻合,与回归模型结果较为符合。

图9 滇池外海pH值及叶绿素a浓度分布情况Fig.9 Distribution of pH value and chlorophyll-a concentration in the outer waters of Dianchi Lake

4.2 pH值与DO的响应关系

在富营养化水体中,DO浓度受到水生生物光合作用和呼吸作用的影响而呈现一定规律的周期性变化,最高时可达饱和溶解氧的3倍。如式(7)与式(8)所示,当藻类进行光合作用时会释放出氧气,使得藻类生长旺盛水域的DO随水生生物光合作用强度增加而上升,且这一过程也伴随着pH值升高。藻类呼吸作用会消耗水中的氧气,使得DO降低,同时释放CO2,抑制水中pH值升高。

DO在一定程度上也影响了水的氧化还原性以及氧化还原反应发生的方向。水中的溶解氧与H+反应产生水,使得水体pH值升高。

DO和pH值之间的响应关系还与多藻型浅水湖泊中生物的死亡、腐解、排泄等生命代谢活动有关。有研究表明,富营养化浅型水体中的pH值和DO存在周期性的协同变化,并且这种变化是各种生命活动以及气象条件等多种内外因素相互作用的必然结果[10-13]。

4.3 pH值与氮、磷的响应关系

从相关分析和回归分析的结果来看,滇池中的氮、磷元素并不会直接使得滇池pH值升高,但这两种元素都是湖泊藻类生长必不可缺的营养元素,滇池中充足的氮磷营养盐使得水体中的藻类迅速生长并且消耗掉水中的碳酸释放OH-,从而间接使得滇池pH值升高,甚至出现异常偏高的现象。

Redfield定律认为,维持藻类健康生长所需的最佳氮磷比为16∶1,当氮磷比高于这一比值时认为磷成为限制性元素,反之则氮为限制性元素[14-15]。滇池外海的氮磷比长期在16∶1上下波动(见图10),给藻类的生长生存空间提供了良好的水质环境。

图10 滇池外海2016—2019年氮磷比变化情况Fig.10 Change of nitrogen-to-phosphorous ratio in the outer waters of Dianchi Lake from 2016 to 2019

研究表明[16],滇池藻类对磷的利用导致磷从沉积物向上覆盖水体的迁移更大,当蓝藻的可溶活性磷浓度因蓝藻增殖的营养需求而下降时,磷从沉积物扩散到上浮水体中后形成稳定的强化流,沉积物释放的可溶性活性磷有70.2%会被蓝藻吸收。因此,外源排放并不是滇池氮磷污染物的唯一来源,沉积物中的磷释放同样可能导致藻华爆发并致使水质恶化[17-18]。

5 结 论

(1)滇池外海pH值异常升高现象主要出现在夏、秋、冬三季,其中夏季pH值与DO呈非显著正相关,秋季pH值与叶绿素a呈非显著正相关,夏、秋两个季节均无水质指标与pH值呈显著性正相关;冬季pH值与DO和叶绿素a均显著正相关,其相关系数分别为0.23、0.42。

(2)夏季回归模型中DO对滇池pH值升高有显著性贡献,且消除了量纲和等级的差异后,模型中DO对pH值升高的贡献是所有变量中最大的;秋季回归模型中叶绿素a对滇池pH值升高有显著性贡献,且对pH值升高的贡献是所有变量中最大的;叶绿素a是冬季模型中系数>0且具有显著性的唯一的自变量,消除了量纲和等级的差异后,叶绿素a对pH值变化的影响远超其他变量。

(3)相关分析和多元线性回归的结果对pH值异常升高主要影响因素的指向不一致,但相关分析的识别结果较为模糊,相关系数无法体现因果关系且只能体现两组数据之间的共变性,而多元线性回归结果能较为精准地对目标变量进行识别,因此,在分析水体pH值异常升高的影响因素时应以多元线性回归结果为主。

(4)综合外海3个季节模型的回归结果来看,叶绿素a对滇池外海pH值异常偏高的直接贡献最大,在夏、秋、冬3个季节的标准化系数分别为0.260、0.231、0.444。因此,浮游植物异常增殖导致的叶绿素a浓度偏高是滇池外海pH值异常的主要原因,如要限制滇池pH值升高必须要治理滇池湖内的氮、磷等营养盐含量,且在减少外源污染物排放的同时,应强化对内源污染物的控制和治理。

(5)由于数据的局限性,本文仅从水质方面来分析滇池外海pH值异常升高的原因。根据现有数据分析结果可知,各类理化因子对pH值的影响还与气象等外部因素有关,但仅从现有的分析无法判断这部分因素对pH值升高影响的大小。不过其余因素也只能间接对pH值产生影响,因此,从总体上看,气象等外部因素并不影响“浮游植物异常增殖导致的叶绿素a浓度偏高是滇池外海pH值异常的主要原因”这一结论。