滇池宝象河流域氮磷流失空间格局解析

严长安， 杜展鹏， 姚 波， 雷 英， 王 凯， 陈静敏， 蒋洪强

1.昆明市生态环境科学研究院， 昆明市环境系统模拟与风险评估重点实验室， 云南 昆明 650032

2.生态环境部环境规划院， 国家环境保护环境规划与政策模拟重点实验室， 北京 100012

氮、磷是维系陆生、水生生态系统结构和功能的两大关键性营养元素，增加活性氮、磷输入负荷有利于提高生态系统的初级生产力，但过量氮磷输入将导致区域水环境质量下降、水生态系统退化和水生态系统服务降低等问题[1]. 因此，有效控制氮磷流失是改善水环境的重要措施之一，而高分辨率的水污染源排放清单是研究氮磷流失的重要基础. 高分辨率的排放清单研究主要集中于大气环境中[2-7]，空间分辨率主要从1°×1°到1 km×1 km不等，受限于网格化数据的获取，高分辨率的水污染源排放清单主要体现在流域层面[8]. 周文婷等[9]以2008年的污染调查结果为依据，研究了茂名市的高州水库水源地的污染排放情况；邱斌等[10]以整个海河流域为研究单元，采用较低的分辨率对流域内氮磷等污染物的排放进行核算和分析，发现流域内平原区的污染程度要远高于山区. 近年来很多研究[11-13]主要集中在污染源排放清单，发现高分辨率的水污染源排放清单解析不仅有利于精准治污，有效控制氮磷污染物的流失，而且有利于精细管理，有效改善水生态环境[9]. 但针对水环境，高分辨率的污染源排放清单只是研究水环境的基础，如何有效控制污染物流失量才是水质持续改善的关键因素，因而基于高分辨率排放清单的污染物流失量的研究对于有效控制污染物流失和改善水环境具有十分重要的理论价值和实践指导意义.

滇池是云贵高原面积最大的淡水湖泊，一直以来接纳了流域内生产生活中排放的大量氮磷污染物，水质常年处于中度或重度富营养化状态. 随着社会经济的快速发展，入滇污染源排放逐渐加剧，不仅会导致滇池水质的进一步恶化，而且将成为流域甚至昆明经济可持续发展的重要短板，因此现阶段滇池已经成为全国重点治理的湖库之一. 宝象河是滇池流域最主要的入湖河流之一，其入滇水质的优劣直接影响到滇池水生态环境[14-18].

因此，该研究基于第二次全国污染源普查数据及核算方法，采用2018年为基准年，对宝象河流域的污染源进行核算，建立宝象河流域高分辨率的氮磷排放清单，通过构建宝象河流域LODEST模型，科学估算流域非点源污染入河系数和宝象河流域氮磷流失量，并对其空间格局进行解析，以期为宝象河流域的水环境精准控污和精细化管理提供有效的决策支撑，同时也为滇池流域入湖污染负荷控制与削减研究提供重要的科学依据.

1 材料与方法

1.1 研究区概况

宝象河属金沙江水系，是昆明古六河之一，位于昆明市的主城区——官渡区，该区位于昆明主城东南部，地处滇池北岸、牛栏江下游，在102°38′E～103°03′E、24°54′N～25°17′N之间. 该流域面积344 km2，地形高程在 1 884.3～2 728.5 m之间，多年平均降雨量为937.1 mm，降雨时空分布不均匀，约85%以上的降雨集中在汛期. 宝象河流域(见图1)属于典型的农业-城市混合型流域，以林地、耕地和建设用地为主. 流域内不同土地利用类型相对集中，宝象河上游区域人类活动较少，河道环境以及径流特征都接近天然状态；宝象河中游流经大板桥街道和阿拉街道，以农业用地为主，主河道周围集中分布着耕地和少量村镇建筑；宝象河下游河道进入昆明主城区，河道周围为城市建设用地，区域内河道均被人工硬化，失去了原始河道的特征，最终于宝丰湿地入湖口注入滇池.

图1 宝象河流域区位图

1.2 数据来源

该研究所需数据的主要来源： ①水文水质数据，宝象河流域干海子水文站多年逐日流量，来源于《长江流域西南诸河水文年鉴(2008—2018年)》；水质数据(TN、TP)来源于课题组试验观测，自2018年1月开始开展了宝象河干海子水文站的水质监测，时间跨度为2018年1月1日—12月31日，采样频率为每周一次，全年共监测了51次. ②污染源数据，城镇生活点源数据来源于昆明市第六水质净化厂.

1.3 子流域划分

为了建立宝象河流域高分辨率氮磷排放清单，将宝象河流域划分为8个子流域，其步骤如下： ①按照最陡坡度(5%)原则和最小集水面积阈值要求，充分考虑气象、土壤、土地利用的空间异质性，利用SWAT 2012软件处理DEM(比例尺为1∶50 000)生成流域河网水系和节点，识别8个子流域；②通过现场勘探落实流域河网水系和节点，同时考虑城市排水片区的影响，人工修正子流域边界，最终确定了8个子流域划分方案(见图2)，其平均空间分辨率达6 km×6 km.

注： 1～8为宝象河子流域编号.

1.4 污染源排放量核算方法

a) 城镇生活点源. 城镇生活点源污染物排放量根据《第二次全国污染源普查生活源产排污系数手册》[18]中的公式〔见式(1)〕来核算. 由于昆明主城区存在城镇生活污水收集系统不完善、雨污分流和清污分流不彻底的问题，城镇生活点源不能完全收集进入污水处理厂. 因此，城镇生活点源的排放量主要包括了污水处理厂的尾水排放量和未收集城镇生活点源排放量. 根据《集中式污染治理设施产排污系数手册》中的公式，将式(1)进一步优化〔见式(2)〕.

PE=Pup-PR

(1)

PE=CwpQwp×10-2

(2)

b) 未收集城镇生活点源. 由于宝象河流域所在区域现阶段无法实现城镇生活点源的完全收集，存在未收集城镇生活点源，因此，可通过城镇污水处理厂实际污水收集率来计算未收集城镇生活点源排放量〔见式(3)〕，可通过《生活污染源产排污系数手册》中的公式〔见式(4)～(6)〕来核算城镇生活点源产生量.

Pnup=Rsc×(1-Pup)

(3)

Pup=kRqc×365×10-5

(4)

q=2.74×QyRy

(5)

(6)

式中：Pnup为未收集点源排放量，ta；Rsc为城镇污水处理厂实际污水收集率，%，根据《昆明市环境保护与生态建设“十三五”规划》，目前昆明主城区旱季污水收集处理率为92.3%；k为城镇综合生活污水折污系数，取值参考《生活污染源产排污系数手册》；R为城镇常住人口数，104人；q为人均日生活用水量，L；c为城镇综合生活污水平均浓度，mgL；Qy为城镇综合生活用水量，104m3；Ry为用水人口；为城镇污水处理厂进水口浓度平均值，为入河(海)排污口排污浓度平均值，mgL.

c) 工业点源. 从第二次全国污染源普查数据中可直接获得工业点源污染物排放量，利用ArcGIS 10.2软件可以将不同行政区的工业点源划分至宝象河流域，得到宝象河流域的工业点源排放量.

d) 规模化畜禽养殖点源. 规模化畜禽养殖点源污染物排放量采用《第二次全国污染源普查农业源普查核算方法》中的公式〔见式(7)〕来核算.

P=Am×(1-Aupi)

(7)

式中：P为规模化畜禽养殖污染物产生量，ta；Am为第m种动物产生的污染物，ta；Aupi为第i种污染物的处理率，%.

e) 种植业非点源. 种植业非点源主要是指化肥流失导致的污染，主要影响因素是化肥施用强度、区域降雨量、地形、土壤等. 根据研究区统计年鉴中化肥施用量情况，结合官渡区及盘龙区的农田面积以及宝象河流域内总农田面积，可得到流域范围内种植业非点源排放量.

1.5 非点源入河系数估算方法

非点源入河系数是估算污染物流失量的关键，目前入河系数的测算方法主要有实地监测[19]、统计模型[20]和机理模型[21]等，受限于流域内详细数据及参数的获取，该研究基于传统方法，提出了非点源入河系数的估算方法. 根据河道断面污染物径流通量与流域污染物排放量估算非点源污染物入河系数(NICi)〔见式(8)〕.

NICi=(RLi-WPPi-IPi)NPDi

(8)

式中：NICi为第i种污染物非点源入河系数；RLi为第i种污染物径流通量，ta；WPPi为污水处理厂第i种污染物排放量，ta；IPi为第i种污染物未收集点源，ta；NPDi为第i种污染物非点源产生量，ta.

污染物断面径流通量的计算需要断面流量与污染物浓度数据. 一般河流径流通量估算方法可以分为统计估算法和模型估算法. LOADEST模型被广泛应用在河道断面径流通量计算中，并且有很好的适用性，已经在宝象河流域得到了较好的应用[22-23]. 因此，该研究运用LOADEST模型估算流域氮磷径流通量.

LOADSET模型是美国地质调查局基于FORTRAN编程语言开发的河流污染物负荷估算模型. 该模型运用离散的水质数据和连续的流量数据估算河流污染物径流通量〔见式(9)〕，并从11种估算河流污染物径流通量的回归方程中优选出拟合度最高的方程，用来模拟不同时间尺度(日、月、年)下河流断面的污染物径流通量.

(9)

2 结果与讨论

2.1 氮磷排放清单

基于第二次全国污染源普查数据及核算方法，建立了宝象河流域高分辨率氮磷排放清单. 由表1和图3可见，2018年宝象河流域TN和TP排放量分别为 1 456.92和191.16 t，排放强度分别为4.24和0.56 t(km2·a)，略高于滇池流域〔分别为3.60、0.38 t(km2·a)〕[8]，一方面是因为近10年滇池流域的氮磷流失量增加，另一方面是由于宝象河流域人类活动强度增大，导致宝象河成为滇池污染最为严重的入湖河流之一[22].

从表1也可以看出：TN主要来源于种植业和城镇生活，占比分别为52.24%和47.53%；TP主要来源于种植业，占比为90.66%，城镇生活点源次之，占比为9.21%. 该研究的流域污染源排放结构与已有的研究[8]结论相似，但主导污染源已变为城镇生活源，体现了污水处理厂对城镇生活源的削减贡献，尤其是TP. 因此，宝象河流域不仅需要继续加强对城镇生活点源的控制，而且应重视农业源的控污问题.

表1 2018年宝象河流域TN和TP排放清单

不同子流域的TN和TP排放量对流域水环境变化的贡献率存在显著性差异(见图3)，依据TN和TP的等标污染负荷对8个子流域污染类型进行划分，即把TN等标排放量与TP等标排放量的比值作为划分标准. 由图3可见，62.5%子流域(即1号、3号、6号、7号、8号子流域)的TN等标污染排放负荷高于TP，表明TN为流域主导污染物类型. 其中，7号和8号子流域的TN和TP等标污染负荷比均超过5号子流域，原因是在宝象河流域西南片区，城镇人口密集，相应的城镇生活污染贡献率较高.

图3 2018年宝象河各子流域TN和TP排放量及其等标污染负荷比

该研究对TN和TP排放量占比较大的3种污染源进行了空间分析(见图4). TN和TP具有一致的空间排放格局，城镇生活点源主要集中在宝象河下游子流域，主要是缘于污水处理厂尾水的排放；未收集城镇生源点源呈现出西南城区明显高于中部城区的特征，主要是因为西南城区是宝象河流域乃至官渡区人口最集中、经济最发达的中心区域，人口密度远高于中部的城郊混合区；种植业非点源主要分布在东北部农田为主的子流域，主要来源于化肥施用.

图4 宝象河流域TN和TP主要污染源排放量空间格局

2.2 非点源氮磷入河系数

基于河流污染物径流通量的估算原理，通过输入宝象河干海子断面的水质监测资料和逐日流量数据，LOADEST模型将根据AIC准则和SPPC准则进行参数率定和优化选择，确定了TN和TP径流通量的最优回归方程，并对回归方程的有效性进行检验[24-26]. TN和TP径流通量回归方程分别如式(9)(10)所示.

ln LoadTN=aTN0+aTN1lnQ+a2lnQ2+aTN3dTNtime

(AIC=-1.198, SPPC=26.685)

(9)

ln LoadTP=aTP0+aTP1lnQ+aTP2dTPtime

(AIC=0.084, SPPC=-5.033)

(10)

式中：LoadTN、LoadTP分别为污染物TN和TP径流通量，kgd；Q为河流断面流量，m3d；dtime为估算模型根据研究对象(TN、TP)水文监测数据的频率复核出的时间累积；a0、a1、a2、a3均为污染物(TN、TP)径流通量最优回归方程的参数；AIC为Aksike信息准则；SPPC为Schwarz后验概率准则.

由表2可知，TN径流通量回归方程的相关系数(R2)为 0.940 2，TP的R2为 0.785 1，证明污染物径流通量回归方程的总体拟合程度较好；显著性检验P值均小于0.05，证明优选出的各污染物径流通量回归方程均具有统计学意义；概率曲线相关系数(PPCC)均大于0.90，证明污染物径流通量回归方程的残差服从正态分布；残差序列相关系数(SCR)取值范围为 0.095 0～0.159 0，证明残差不存在序列相关性. 综上，LOADEST模型优选出的各污染物径流通量回归方程较为有效，适用于对宝象河TN和TP径流通量的估算.

表2 氮磷径流通量回归方程检验结果

通过估算，2018年宝象河流域干海子断面的TN和TP径流通量分别为270.49和11.19 ta，与李娜等[22]研究的宝象河氮磷多年平均径流通量的估算值(239.2 ta)较为相近. 从逐日的时间尺度来看，宝象河干海子断面TN和TP的平均日径流通量分别为741.07和30.662 kgd，最小值分别为292.55和11.817 kgd，最大值分别为 2 169.20 和120.08 kgd. 如图5所示，全年TN和TP径流通量具有显著的季节性特征，径流通量主要集中在6—9月，夏季较高，夏季农业氮磷肥施用量大，大量的氮磷被雨水冲刷进入附近的水体，而春冬季较低. 随着逐月降雨量的增多，TN和TP径流通量逐渐增大，8月达到最大值，之后直到12月呈逐月减少的趋势. 因此，降雨是流域氮磷污染物流失的主要驱动因素[20]，应重点关注汛期流域氮磷流失的控制.

图5 2018年宝象河流域TN和TP逐日负荷

根据河道污染物径流通量与区域污染物排放数据，采用非点源污染入河系数计算公式，得到宝象河流域TN和TP的非点源入河系数(见表3). 结果表明：宝象河流域TN的非点源入河系数为0.297，95%置信区间为0.281～0.312，TP的入河系数为0.048，95%置信区间为0.041～0.053. 与高伟等[27]基于灰水足迹理念所测算得到的长江流域非点源氮磷淋溶流失系数(分别为0.143和0.021)相比，该研究基于径流通量得到的非点源氮磷入河系数偏高，可能的原因是：①该研究核算非点源产生量偏小，可能与未考虑大气沉降中的氮磷输入量有关，宝象河流域来源于大气沉降的氮约占流域非点源污染负荷的38%[28]；②该研究测算的未收集生活点源偏低，作为城中村较多的官渡区，城镇生活污水收集率要低于昆明市主城区. 尽管该研究估算的非点源氮磷入河系数偏高于同类研究，但仍然处于非点源氮磷的淋溶流失范围内[29]，并不影响流域氮磷流失量的空间格局解析. 因此，该研究提出的从径流通量角度核算流域氮磷流失量的方法能够较好地揭示流域氮磷流失空间分布规律，其估算结果具有一定的可靠性.

表3 宝象河流域TN和TP的非点源入河系数

2.3 氮磷流失量及空间格局解析

2018年宝象河流域的TN和TP流失量分别为432.28和18.57 t. 苏斌[30]于2019年通过试验观测估算了2017年宝象河入湖TN径流通量(687.12 t)，发现该研究的氮磷流失量估算值偏低，可能主要与2018年降水量显著减少有关. 由于滇池流域入湖河流中的氮磷主要来自于城市面源和农业面源[30-32]，氮磷通过降水-冲刷作用进入河流，2018年降雨量的减少致使宝象河氮磷径流通量偏小. 非点源氮磷流失量除受降雨因子影响外[20]，还受到流域地形因子、土地利用类型和土壤类型的影响[29,31,33]，降水和地形因子的不均匀性最终导致了宝象河流域氮磷流失的空间异质性.

由图6可知，TN流失的空间格局主要集中在1号和8号子流域，流失量分别达到111.91和147.25 t，两个子流域的流失量占比整个流域流失量的60%以上. 究其原因主要是：1号子流域位于宝象河流域的上游区域，是流域内种植业非点源污染的主要区域，且受到降雨影响最大[29]；8号子流域位于入湖口位置，是流域内人口最为密集的地区之一，城镇生活点源污染较为严重，是流域内TN流失最为严重的区域. 此外，4号和5号子流域因种植业与人口两个原因叠加，致使其成为流域内TN流失的次级重点区域. TP流失的空间格局主要集中在上游区域的1号、2号、4号子流域，流失量分别达到3.06、1.27和1.76 t，3个子流域的流失量占比高于70%，成为全流域TP污染最为严重的区域. 这主要与农村生活和农业生产息息相关，因而进一步证实了TP的主要污染来源为农村农业面源，如何有效控制农村农业污染源已成为当前亟待解决的关键问题.

图6 2018年宝象河流域的TN和TP流失量空间分布特征

从图6也可以看出，流域氮磷流失量空间分布特征具有一定的相关性，1号子流域的氮磷流失量是流域内污染最为严重的区域，TN的流失量占比为26%，TP的流失量占比为36%. 已有研究[28]显示，降雨因子和地形因子是导致流域氮磷流失量较大结果的最重要的因素. 因此，有效防控1号子流域的氮磷污染，是控制整个流域氮磷污染的最为重要的途径之一. 此外，氮磷流失量空间差异最大的是流域入湖口所在的8号子流域，其TN的流失量较大，而TP的流失量却较小，导致这样结果的原因主要是8号子流域承接了污水处理厂的尾水，其TN的负荷较大，而TP的负荷较小，如何提升污水处理厂TN的出水标准是解决流域TN污染的关键.

该研究通过与高伟等[8,34]研究成果的对比，进一步证实了近10年宝象河氮磷的主要污染源依然在宝象河上游和入湖口附近的流域[35]. 然而，高伟等[8]所计算的滇池流域的高分辨率排放清单，由于是基于DEM数据，所提取的子流域存在一定的偏差性. 该研究以此为基础，对流域出口进行了实地勘验，并重新确定了流域出口，使得子流域划分更加贴近实际情况，对氮磷流失的空间异质性和污染源类型的划分更为精确，因而将更有利于提出切合实际的污染物有效控制的研究建设方案.

由于流域是一个复杂的复合系统，对污染物的核算存在诸多的影响因素，无法对其进行一个完全准确的核算. 在排放因子核算方面，虽然该研究已经尽量增加了现场观测和实际测量次数，但核算中仍存在一定的偏差性. 同时，该研究仅针对各污染源的空间分异进行了分析，未能对年内和年际变化进行有效分析. 在数据核算方面，该研究基于第二次全国污染普查数据，对流域氮磷污染进行核算，但是由于未对大气沉降[36-37]和水土流失[38-40]等情况进行充分考虑，致使研究结果可能存在估值偏小的情况.

3 结论

a) 该研究基于第二次全国污染源普查数据，建立了2018年宝象河流域高分辨率的氮磷排放清单，分别为 1 456.92 和191.16 t. 通过测算，流域内种植业非点源是贡献最大的污染源，其次是城镇生活点源和未收集点源. 通过解析，氮磷排放量呈现显著的空间异质性，流域内单位面积的氮磷排放强度总体呈现外高内低的分布. 依据TN和TP的等标污染负荷子流域污染类型进行划分，发现TN为流域主导污染物类型，对于该区域的水环境安全具有较大风险.

b) 基于宝象河流域氮磷径流通量分别为270.49和11.19 ta，TN和TP径流通量年内变化特征与降水量有着较高相关性，具有显著的季节性差异，主要集中在6—9月，呈现夏季较高、春冬季较低的特征. 此外，通过测算分析，宝象河流域TN和TP的非点源入河系数分别为0.297和0.048.

c) 该研究证实了流域氮磷流失量主要受降雨因子、地形因子、土地利用类型等影响，2018年宝象河流域的TN和TP流失量分别为432.28和18.57 t，降雨和地形的不均匀性造成了氮磷流失量呈现显著的空间异质性，1号子流域是整个流域内污染最为严重的区域，而TN流失量集中分布在宝象河上游的1号子流域和宝象河下游的8号子流域，TP流失量集中分布在流域上游农业生产生活最为活跃的子流域.