利用弹性网对PLS佐证分析城市水生态足迹及驱动因素——以贵阳市为例

许国钰, 任晓冬, 杨振华, 肖 杰, 刘立冰, 高 洋

水资源作为一种战略性资源，关系着国家安全与人民的生产生活。随着经济发展、工业化和城市化不断推进，水资源短缺，水污染，生态恶化等问题，已经成为制约城市发展与生态协调的一个重大问题。开展针对水生态足迹的研究已成为实现水资源可持续开发利用的一个重要环节。生态经济学家William[1]提出的生态足迹概念为可持续发展提供了理论依据。Mathis等[2-3]对生态足迹计算原理和方法进行进一步的拓展，并广泛应用于不同空间尺度的定量评价体系，并由此衍生出水足迹，水生态足迹等概念。黄林楠等[4]建立水生态足迹和水资源承载力模型，并对整个水资源账户进行细化。潘真真等[5]基于水生态足迹并将概念引申为水生态占用，并提出生态系统供给及净化服务功能的水生态占用研究。在此之后，国内学者对水生态足迹的研究，从单一的核算水足迹及其历年变化特征；发展到胡永红等[6]运用ARIMA模型的大样本建模对水生态赤字进行预测。杨振华等[7-9]分别采用GRNN模型、SVR(支持向量回归机)模型对时间短、变异性大的小样本容量对城市水生态足迹进行预测，揭示水生态足迹的时间序列规律。虽然取得了良好的预测精度，但是缺乏对影响水生态足迹驱动因素的分析。本文以地质条件特殊而导致工程性缺水的城市贵阳为例，计算其2002—2016年水生态足迹时间序列，引入STIRPAT模型对影响其水生态足迹的驱动因素进行分析。虽然传统的线性回归方差最小，但是由于数据拥有大量预测变量，而且存在严重的相关性，导致矩阵病态。之后采用偏最小二乘法(PLS)及弹性网对STIRPAT进行修正，PLS回归兼具主成分分析和多元回归的优点，克服了变量间存在强烈的交互相关而导致多元共线性[11]。弹性网是一种是对PLS的算法优化，它是对回归系数采用联合惩罚函数的有偏估计。并且在筛选变量方面更加精准。本文在丰富了水生态足迹理论的实践及应用的同时，运用PLS及其弹性网对影响城市水生态足迹驱动因素进行分析，且进一步探讨如何减缓贵阳市水生态足迹，对城市水资源的可持续利用提供一定的决策参考。

1 研究区概况

贵州省贵阳市地处云贵高原东斜坡地带，东经106°71′—107°17′，北纬26°11′—27°22′。属于明显的亚热带高原性季风气候，温和湿润。贵阳市地处长江流域的乌江水系和珠江流域的红水河水系，境内虽有98条河流。但是由于贵阳市地处我国喀斯特地貌的中心地带，水循环受“二元”结构影响导致地表蓄水能力差[12]。独特的地质条件导致用水率极低，加上山高坡陡峡谷深，很难建造水库，河水大量流失，属于典型的工程性缺水[12]，对耗水严重的工农业及相关的产业的约束性强。喀斯特地貌地表起伏度大，土层稀薄，多裂缝隙及落水洞发育，易与地下水形成联动污染[12]。贵阳市已探明矿种有52种，主要有煤(9.00×108t)、铁(2.40×107t)、硅、钒土、磷、硫、汞等矿产资源。其中铝土矿的保有储量4.30×108t，磷矿4.60×108t，是全国3大磷矿基地之一。2016年3大产业结构比重为15.61∶39.47∶44.9，第一、第二产业占国内生产总值的比重过高，农业和大规模工业行业耗水量过大。

2 研究方法

2.1 基本理论

2.1.1 水生态足迹 生态足迹理论是用于比较自然资源消耗量与给定区域内人口与生态承载力的研究，水生态足迹是由生态足迹理论发展而来。水生态足迹是指维持人口对资源消费和废弃物消纳所必需的水域和水资源用地面积[13]。城市水生态足迹又直指城市人口消费所必须的水域和水资源用地面积，由两部分组成：城市的水量(工业、农业、生活)生态足迹；城市的水污染(氨氮、COD)生态足迹,包括维持城市废弃物消纳所必须的水域和水资源用地面积。但实际在应用层面，城市用水量生态足迹和水质生态足迹两部分叠加即可求出城市水生态足迹[5,10,14]。其计算模型为：

EFWf= WFWF+EFWQ+rw×

(1)

式中：EFwf——水生态足迹(hm2)； EFWF——水量生态足迹(hm2)； EFWQ——水质生态足迹(hm2)；rw——全球水资源均衡因子；Pw——全球平均水资源生产能力(m3/hm2)；Uwr——城市用水资源总量Uwr=Uwri+Uwra+Uwrc；Uwri——工业用水；Uwra——农业用水；Uwrc——生活用水；UCOD，UNH3——COD和氨氮排放量；PCOD，PNH3——水域对污染物COD和NH3的吸纳能力(t/hm2)。

2.1.2 水生态承载力核算 水生态承载力是反映水生态系统的供给能力，是指支持水生态系统在保持健康下维持人类活动的阈值。城市的水生态承载力可通过供需两方面的比较来计算，用于测度城市水生态可持续发展状况，其计算模型为：

(2)

式中：ECw——城市水生态承载力(hm2)；φ——区域水资源产量因子；rw——全球水资源均衡因子；Q——城市水资源总量(m3)；Pw——全球平均水资源生产能力(m3/hm2)。

2.1.3 水资源可持续利用指数 水资源可持续利用指数是指水生承载力与水生态足迹的比值，实际上是反应水资源的生态压力指数。水资源可持续利用指数计算公式为：

EQ=ECw/EFw

(3)

式中：EQ——水资源可持续利用指数，EQ>1时说明水资源可持续利用程度高，EQ<1说明水资源缺口大，水资源可持续利用程度低。

2.2 STIRPAT模型(环境压力模型)

STIRPAT模型[14-15]是通过对人口，财产，技术3个影响因素对环境压力之间的关系进行评估。其具体形式为：

(4)

式中：Ii——环境压力、本文采用贵阳市水生态足迹来表征；Pi——人口总量；Ai(affluence)本文引入的贵阳市人均GDP；Ti——技术水平；η——模型的系数；μ,θ,φ——公式中几个因素的人文驱动力的指数；εi——模型误差项。由于STIRPAT模型是一个具有多个自变量的非线性模型，在公式两边取对数处理后：

lnEFi=η+μlnPi+θlnAi+φlnTi+lnεi

(5)

式中：EFi——因变量；Pi，Ai，Ti——自变量；η——常数项；εi——误差项，把公式(5)进行多元回归，根据弹性系数的概念，模型中3个自变量，如发生1%的变化，就会分别引起因变量各自发生μ%，θ%和φ%的变化。本模型是对人文驱动因素与环境压力的一种非线性关系假设。

本文将公式(5)中技术指标Ti分解为结构化指标Ni(本文指的是第一二产业占总产值比重，因为贵阳市工农业用水比重超过总用水量的60%)和城市化指标Mi(城市化率)。即lnN和lnMi两项，为了验证各驱动因素与水生态足迹间是否存在倒U形的EKC曲线[15-16]，将模型调整为公式(6)：

lnEFi=η+μlnPi+〔θ1lnAi+θ2(lnAi)2〕+〔κ1lnNi+κ2(lnNi)2〕+〔ω1lnMi+θ2(ω2lnMi)2〕lnεi

(6)

再对公式(6)求一阶偏导数，就得出富裕度(即人均GDP)，对水生态足迹的弹性系数为θ1+2θ2(lnAi)，结构化(即第一二产业对占总产值比重)对水生态足迹的弹性序数为κ1+2κ2(lnNi)，城市化对水生态足迹的弹性系ω数为1+2ω2(lnMi)。

2.3 偏最小二乘法(PLS)

PLS允许自变量存在多重相关性，并且在样本点个数少于变量个数的条件下进行回归建模[11]；自变量X经标准化处理后的数据矩阵记为E0=(E01,…,E0P)n·q，设p个自变{x1,…,xp}，因变量Y经过标准化处理后的数据矩阵记为F0=(F01,…,F0q)n·q，并设有q个因变量{yi,…,yq}。为了研究与X和Y的统计关系，取n个样本点，由此构成了自变量与因变量的数据表X=(x1,…,xp)n*p和Y=(y1,…,yq)n·q。用偏最小二乘法回归分别在X与Y中提取出t1作为x1,…,xp的线性组合和u1作为y1,…,yq的线性组合。在第一个成分t1和u1被提取后，偏最小二乘法回归分别实施X对t1的回归以及Y对t1的回归。算法将会在方程达到了满意精度后终止。反之，将利用X被t1解释后的残余信息，以及Y被t1解释后的残余信息进行第二轮的成分提取，直到能达到一个较为满意的精度为止。该方法同时采用自变量投影重要指数VIP判断驱动因子的重要性。一般认为，VIP大于1.00，表示该自变量是因变量的重要确定驱动因公式为：

(7)

2.4 岭回归，LASSO，弹性网

3 数据来源及模型参数

本文数据来源于2001—2017年《贵阳市统计年鉴》《贵阳市水资源公报》《贵阳市国民经济和社会发展统计公报》。城市化率，学术界有人口比重指标法，农村城镇化指标，现代城市化指标法，调整序数法，城镇土地指标法。本文采用的是人口比重指标即城镇人口占总人口的比例。本文关于水生态足迹模型计算参数，其中计算参数参照黄林楠等[4]的研究成果，全球均衡因子在40 a内只有微小变化可采用静态值。本文取2000年的值作为计算参数，水质生态足迹中的计算参数根据《地表水质量环境标准(GB3838-2002)》，表1为公式(1)—(2)中计算参数。

表1 水生态足迹计算参数 t/hm2

4 结果与分析

4.1 贵州省贵阳市水生态足迹时间序列

根据图1所示结果，2002—2016年贵阳市总的水量生态足迹呈波动变化。2002—2008年，水量生态足迹，大幅上升后逐步下降，变化幅度达7.75×105hm2。2009—2016年贵阳市水量生态足迹整体呈缓慢上升趋势，年均升幅不足6.10×104hm2；分别从工业、农业、生活用水生态足迹来看，工业用水生态足迹自2009年起大致呈下降趋势，但仍然是影响贵阳市总的水量生态足迹的最大因素。由于灌溉用水效率提升，农业用地萎缩，导致农业用水生态足迹呈逐年下降趋势。同时，由于十二五以来，贵阳市推进快速城市化战略，导致城市人口激增，生活用水生态足迹从2002年的3.39×104hm2暴涨至2016年的7.13×105hm2，年均增幅超过17.175%。

图1 贵阳市2002-2016年水量生态足迹变化特征

由图2所示，2002—2016年贵阳市水质生态足迹总体呈先下降后又波动式上涨。其原因是贵阳市2005年起实施新的排放标准和工业企业改进排污技术。导致水质生态足迹中COD与氨氮排放量下降。2011—2016年，污水中的COD排放量基本稳定，但是由于城市人口迅猛增长，生活污水中氨氮排放量增加，导致氨氮生态足迹上涨，造成水质生态足迹有强烈反弹上升的趋势。

由图3所示，2002—2016年贵阳市水生态足迹在1.28×107～1.60×107hm2浮动，呈先上升后微降再上升的趋势。水生态承载力的变化幅度也十分明显，2006，2011，2012，2013年这4个年份是贵阳市水生态承载力的低点，主要是由于贵州省出现大旱等极端天气，导致水资源总量波动，进而影响水生态承载力波动。

水资源可持续利用指数又是由水生态足迹和水生态承载力共同决定的。2002—2016年贵阳市水资源利用指数都小于1.0，其变化范围0.315～0.624。表明贵阳市水生态系统常年处在超负荷状态，水资源持续利用面临很大压力。总体而言，贵阳市对水资源的消费量远远超过生态系统的承载力范围，总体处于不可持续的状态。随着贵阳市经济发展，城市化的快速推进，水生态足迹将进一步增加，贵阳市的可持续发展形势更为严峻。

图2 贵阳市2002-2015年水质生态足迹变化特征

图3 贵阳市2002-2016年水生态足迹、水生

4.2 基于STIRPAT模型的LS回归分析

贵阳市水生态足迹变化的驱动因素主要包括总人口，富裕度，结构化和城市化水平。因此本文选取这4个指标，根据公式(6)构建贵阳市水生态足迹变化的STIRPAT模型。表2模型拟合度介于24%～66.56%。模型1是以lnEF为因变量、lnA和(lnA)2为自变量的STIRPAT方程，其调整R2=0.24，DW统计量为1.378，sig.<0.05，(lnA)2对lnEF存在显著影响。在模型2中添加自变量lnP，lnP对lnEF不存在显着性影响。模型3添加了自变量lnN，在模型3中，lnA，lnN，(sig.<0.05)对lnEF解释提升，调整R2>0.24。模型4添加自变量(lnA)2(sig.<0.01)和(lnN)2(sig.<0.05)，(lnA)2对模型有显着性影响。模型5，6，7是在模型3，4的基础上考虑lnN和(lnN)2，lnM和lnM2，但这几个自变量都对lnEF影响不显著。模型8是把所有影响lnEF的驱动因素添到回归模型中，调整R2=0.665 6。由表3得出，对lnEF有影响的自变量主要是lnA，(lnA)2,和lnN，(lnN)2。在8个模型当中lnA的弹性系数在-0.546～-0.06之间，(lnA)2的弹性系数在0.186～0.316之间。表明lnA对lnEF存在负影响，(lnA)2系数为正，说明贵阳水生态足迹不存在环境EKC假说。lnN和(lnN)2的弹性系数从分别是1.145 1～6.325和-1.128～4.086。因为(lnN)2的弹性系数出现负值，说明第一二产业占总产值比重在某种程度上存在环境EKC假说。普遍认为自变量的方差膨胀因子VIF>10意味着变量之间具有严重的共线性。尤其是模型7.8共线性更为严重，完全不适合用最小二乘法LS进行无偏估计。表明用上述方法进行水生态足迹驱动因素分析结果的准确性受到质疑。

表2 贵阳市水生态足迹变化的STIRPAT模型分析

注：*sig<0.05；**sig.<0.01；***sig.<0.001。

表3 自变量投影重要指数(VIF)值

4.3 偏最小二乘法(PLS)分析结果

根据PLS回归普遍用于解决自变量的多重共线性问题，对上述引入所有影响水生态足迹的驱动因素(lnP，lnA、(lnA)2，lnN,(lnN)2，lnM)的模型8进行修正，其结果由表4所示，当t=1，即提取一个主成分时，调整R2=0.05，说明提取一个主成分时能达到对因变量5%的解释；当t=2时，调整R2=0.28，表明提取2个主成分能对因变量28%的解释。同理当t=5为最佳成分数，建立基于PLS的STIRPAT模型，调整R2=0.714，达到最大值。提取5个主成分对因变量产生71.4%的解释，Press的最小值达到5.306 2，此时lnP的弹性系数为0.918，当lnP发生1%的变化时，将会引起lnEF发生0.918的变化，此时lnA，(lnA)2，lnN，(lnN)2,lnM，(lnM)2发生1%的变化，分别会引起lnEF发生-0.492，0.282，1.223，-0.003 7，0.778和-1.271的变化。但是PLS模型属于未标准化回归方程，其系数大小反映5个自变量对因变量的影响，因此需要进一步引入VIP值分析各因变量对自变量影响。

VIP值结果如图4所示，lnP,(lnA)2>1表明总人口，人均GDP及二次项是贵阳市水生态足迹变化的重要影响因素，其变化按VIP重要性排序为：人均GDP的二次项>人口。但是R2介于0.05～0.714表示各模型拟合效果不够理想，此时主成分并不能解释自变量和因变量之间的关系，基于PLS对STIRPAT模型建立的水生态足迹与其驱动因素之间的线性回归模型不合理，与文献[21-22]研究结果不一致。

表4 基于PLS的STIRPAT模型

图4 贵阳市水生态足迹变化的重要影响因素VIP值

4.4 基于弹性网对PLS的优化

偏最小二乘法适用于观测数少于预测变量数以及预测变量之间存在共线性的情况，但是计算过程需先提取主成分，在一定程度上会同时损失预测变量和因变量的信息，提取因子的个数越少，损失的信息越多。弹性网(岭回归和LASSO回归为其特例)是对回归系数采用联合惩罚函数的有偏估计，既能达到变量选择的目的，又能保证群组效应，结合了岭回归和LASSO的优点[17-19]。本研究纳入2002—2016年贵阳市水生态足迹相关数据，利用STIRPAT模型探究影响贵阳市水生态足迹4个驱动因素，研究水生态足迹和人口，富裕度，结构化，现代化之间的相关性，并探讨是否存在倒U形环境库兹涅茨曲线关系。本研究的模型中共纳入15条记录，8个自变量，包含信息相对较少，所以使用弹性网进行数据拟合，所有统计采用R软件中glmnet包完成。

表5结果表明，岭回归的惩罚系数初始值较高，随着lamda的增加，系数的变化成波动状，并非逐渐减小。同时按照倒U形结构假设，二次项的系数为正，一次项系数为负，但在岭回归结果中lnN的系数为负，这与模型规律不符。LASSO和弹性网中在相同的alpha变化情况下，各回归系数的压缩趋势类似，均能逐渐收缩系数，达到变量选择的目的，但LASSO回归的回归收缩速率较快。此外，从对应的数据表中可见，LASSO回归中系数中出现从0上升后回到0的波动，同时lnA和(lnA)2，lnN和(lnN)2，以及lnM和(lnM)2构成3个不同群组，但三者的回归系数均在alpha细微的增加情况快速收缩到0，模型结果无法合理解释。弹性网的结果显示回归系数收缩平滑，并未显示出群组效应，在惩罚系数较大时仍有部分指标系数不为0。其中人口因素lnP在不同的模型的不同惩罚系数中均持续不为0，表示其重要作用。最后对在alpha=0.5的条件下进行弹性网的交叉验证，根据glmnet包的建议，选择lambda在一个标准误差的情况下为最优条件，选择lambda=0.008是为最优的模型结果。此时lnP,(lnA)2，lnN，(lnN)2，lnM，(lnM)2在STIRPAT模型当中是有意义的，此时lnP，(lnA)2，lnN，(lnN)2，lnN，lnM，(lnM)2的弹性系数分别为0.484，0.428，-0.097，-0.187，0.036。由岭回归、LASSO、弹性网3种回归的系数变化图可知(图略)，当lnP发生1%变化时，lnEF发生0.484%，0.428%，-0.097%，-0.187%，0.036%的变化。由于lnN,(lnN)2的系数一个出现正值，一个出现负值表明贵阳市的水生态足迹与第一二产业占总产值比重间存在环境EKC假说。lnM,(lnM)2的系数同样是一个正值，一个负值，这同样表明贵阳市水生态足迹与城市化水平之间存在环境EKC曲线。

5 结 论

(1) 贵阳市的水生态足迹是远大于水生态承载力的，这与该市地处喀斯特发育严重地区有关。其特点是土层稀薄脆弱，岩石渗透严重，地表蓄水能力弱，工程性缺水明显，导致整个市水资源利用效率低。农业水生态足迹受到自然条件的巨大影响。工业主导的水生态足迹存在大量的COD以及氨氮排放的，同时由于地表与地下属典型的喀斯特地貌，土层渗透系数较高，地表水与地下水形成联动污染，导致水体污染严重。而且在喀斯特地型状况下，岩体破碎，平缓沉积地形缺乏，虽然整个城市划分6个区，但是建成区域异常密集且坡度大，城市人口密集导致供水紧张。由于贵阳市正处于快速城市化和社会经济的高速发展期，虽然在2002—2008年生活用水水生态足迹有小幅波动，但是整体处于飞速上涨期。贵阳市水生态足迹的潜在因素表现在人口的增加、产业结构不够合理和城市化水平的过快增长。

(2) 根据STIRPAT建模分析贵阳市水生态足迹变化的重要驱动因素，其LS(最小二乘法)回归结果如下：人均GDP及二次项，第一二产业占总产值比重是贵阳市生态足迹变化最重要驱动因素，但是对水生态足迹有重大影响的总人口、城市化指标被排除在外。VIF值的计算结果表明各个驱动因子之间存在多重共线性。之后采用PLS(偏最小二乘法)对模型进行修正，当t=5时，R2=0.714，得出总人口，人均GDP及其二次项是影响贵阳市水生态足迹的总要因素。按VIP重要性排序为：人均GDP的二次项>人口，一个城市的水生态足迹仅仅只有这2个影响因素是不符合常理的。于是本文引入弹性网回归对其模型进行佐证，其回归系数采用联合惩罚函数的有偏估计，达到了变量选择目的，又能保证群组效应。lambda=0.008是为最优的模型结果，lnP,(lnA)2,lnN,(lnN)2,lnM，(lnM)2在STIRPAT模型当中是有意义的。即贵阳市总人口，人均GPD，第一二产业占总产业比值，城市化率都对贵阳市水生态足迹产生影响。而且第一二产业占总产业比值与贵阳市的水生态足迹间存在环境EKC假说，城市化与贵阳市的水生态足迹间也存在环境EKC假说。

(3) 由于第一第二产值占总产值比重与城市化率和贵阳市水生态足迹存在环境EKC假说。虽然受整个论文思路所限，本文没有进一步讨论环境EKC假说其临界点到底在哪里。但是这2个驱动因素出现了环境EKC假说，就表明整个贵阳市在推进快速城镇化的发展过程中需要调整步伐。进一步优化产业结构，优化水资源配置，按序调控。在优先保证生活用水前提下，合理规划和保障社会经济发展的用水，保证一定的生态用水。把工业、农业的发展规模、区域产业结构定位与水资源的配置结合起来，引导与水环境承载力相协调的发展布局。并推动涉水产业绿色转型与发展，积极发展生态绿色农业，大力降低农业面源污染。工业区应积极发展循环经济，开展清洁生产，中水回用，实行用水梯级循环，不断降低废水污染物排放量，淘汰废水排放量大的落后产能。城市倡导节水生活方式，另在生活污染为主的河流，加大污水收集范围，提高生活污水处理率，减少污水中氨氮排放量。并实施河道垃圾清理、清淤疏浚、生态修复等综合措施，全面降低污染物入河负荷，逐步改善水质。