螺髻山北麓地下水化学特征与水质评价

吴君毅，刘洪，欧阳渊，李樋，张景华，张腾蛟，黄勇，段声义

（1.中国地质大学（北京），北京 100084；2.中国地质科学院研究生院，北京 100037；3.中国地质调查局成都地质调查中心，四川 成都 610081；4.成都理工大学地球科学学院，四川 成都 610059；5.华东冶金地质勘查局测绘总队，安徽 合肥 230088）

地下水作为水资源和水循环的重要组成部分，其水质的优劣关系到当地经济建设、人民的生产生活、生态平衡等一系列问题（Chen，1994；Adimalla，2020；张俊等，2021；党学亚等，2022；冯嘉兴等，2023）。螺髻山北麓横跨四川省凉山彝族自治州西昌市、德昌县与冕宁县，研究区经济发展滞后，以粗放型农业、畜牧业为主，生态环境脆弱（杨帆等，2018）。因此研究地下水化学特征、形成机制以及评价水质对该区保护和合理利用开发地下水资源，防止和控制地下水污染，促进农牧业绿色生产，保障人民身体健康，加强生态文明建设有着重要研究意义。

地下水化学特征研究是水体化学组成与地下水演化过程研究的基础（Abbas et al.，2021），起到反映地下水的补给、途经区域岩性、气象水文与环境特征、反映水化学演化的控制因素等作用（Li et al.，2015，2021；杨芬等，2021；Wali et al.，2021；杜金龙，2022）。能够有效揭示地下水与环境的相互作用机制以及离子交换过程，对地下水水质差别的原因有一定解释作用，可为合理开发利用水资源，改善水质提供科学依据。

国内外对于地下水水质的研究中，评价方法主要有内梅罗指数法、模糊综合评价法、主成分分析法和多元分析法等（Nemerow et al.，2009；Zhang et al.，2012；李连香等，2015；张志君，2020；周及等，2020；时雯雯等，2021；曾庆铭等，2021）。上述评价方法采用传统统计学方法虽各有其优势，但都不能很好的分析水质等级与评价指标之间复杂的非线性关系，可能存在一定的主观干预（李海涛，2020a）。近年来，随着计算机科学的发展与应用，人工神经网络（Artificial neural network，ANN）模型能够高效拟合非线性关系数据，对其做出分类或预测，ANN 模型在非线性数据拟合方面存在优势，能够避免确定性模型处理实际问题时的过度简化问题（陈能汪，2021）。对于利用人工神经网络进行水质评价最常用的网络模型是BP 神经网络，虽然该模型有收敛速度慢、结构复杂，且无法反映主要影响因子等缺点（崔永华等，2007），但其优点克服了主观赋权对于评价结果的影响，评价结果客观、合理且精度高，能在处理不同时间和空间水样数据组合时的效能优于其他方法（Kumar et al.，2020）。BP 神经网络因此被广泛应用于水质评价，并且通过不同的算法如萤火虫算法（颜建等，2020）、头脑风暴优化算法（李海涛等，2020b）、DPA 算法（徐康耀等，2015）等算法能够进一步提高评价的结果。

笔者在前人对研究区水文地质背景总结研究的基础上，研究地下水化学特征。利用Pytorch 搭建基于RMSprop 梯度下降算法的BP 神经网络模型，对研究区地下水水质进行综合客观评价。

1 研究区概况

1.1 区域自然条件

研究区位于四川省凉山彝族自治州西昌市、德昌县与冕宁县的交界处的螺髻山北麓（E 102°7′～102°30′，N 27°30′～27°46′），北接邛海盆地，南临螺髻山脉，东接马雄梁子，西靠安宁河磨盘山。研究区的总体地势为两山夹一河，山地河流整体为南北向。在气候上，具有年温差小、日温差大、干湿分明、雨量充沛、降水集中、日照充足等特点，属于亚热带高原季风气候。西侧的牦牛山南麓和安宁河河谷分别为构造侵蚀中山和中山宽谷平原，中部螺髻山北麓为构造侵蚀中山地貌，高山与盆地最大高差为500～1 500 m。螺髻山山岭海拔为3 000～4 200 m，呈近南北走向，在研究区构成了安宁河和则木河的分水岭。研究区西部为安宁河流域，安宁河从北向南贯穿西昌市，于米易县附近注入雅砻江，最终汇入金沙江，全长为320 km，流域面积为1.1 万km2，安宁河的主要支流有18 条，其中流域面积大于500 km2的支流主要是孙水河、海河、茨达河、锦川河。主流与支流多以直角交汇，成羽状水系。研究区东南部为则木河流域，则木河发源于普格县特尔果乡阿则木山，属于黑水河右岸一级支流，从特尔果乡一带由北向南流经五道箐乡、特补乡，螺髻山镇等地，最终在普格县的中梁山南端流入黑水河，其流长约为54 km，流域面积为668.2 km2（孟晶晶等，2018）。

1.2 区域地质背景

在大地构造划分上，研究区位于上扬子西缘的康滇断隆带中北段（刘洪等，2020）。研究区出露的地层（图1）主要为震旦系—第四系，震旦系以陆相火山碎屑岩类、陆源碎屑岩类和碳酸盐岩类为主，零碎分布在螺髻山北麓、马雄梁子地区。寒武系在本区出露较少，以碳酸盐岩为主。中生界（三叠系—白垩系）在本区大面积出露，以陆相砂泥质碎屑岩为主。新生界（古近系—第四系），以半固结的碎屑岩为主零星分布在安宁河谷东岸及大箐梁子等地，第四系冲洪积物（Qapl）沿安宁河、则木河河谷以及山谷沟谷分布。岩浆岩在安宁河西岸出露，主要为花岗岩类。研究区构造线复杂，主要的断裂包括安宁河断裂带（南北向）和则木河断裂带（北北西向）。由于地质构造复杂，地质条件多样，现代构造活动强烈，岩石稳定性差、岩面风化强烈、疏松破碎、碎屑物质丰富，形成了易发生地震、崩塌、滑坡、泥石流、水土流失的地质背景（张腾蛟等，2020）。

图1 研究区水文地质图Fig.1 Hydrogeological map of the study area

1.3 区域水文地质背景

根据研究区地层岩性、地下水类型分为6 类，分别是松散堆积层孔隙潜水-承压水区（A1 型）、松散堆积层孔隙潜水区（A2 型）、碎屑岩裂隙层间水区（B1 型）、碎屑岩孔隙-裂隙层间水区（B2 型）、岩浆岩裂隙水层间水区（C 型）和碳酸盐岩裂隙层间水区（D 型）。

A1 型在河谷区第四系冲洪积物区域，地下水量丰富，为孔隙潜水-承压水；A2 型在山麓区第四系冲洪积物中地下水量中等-贫乏只有孔隙潜水；B1 型地层有开建桥组、列古六组、官沟组、飞天山组、益门组、新村组、牛滚凼组等地层，地下水量中等-贫乏的裂隙层间水；B2 型地层有白果湾组、小坝组等碎屑岩地层地下水量丰富的孔隙-裂隙层间水；C 型为磨盘山的地区的晚三叠世花岗岩类地下水量贫乏的裂隙层间水；D 型分布在螺髻山区的观音崖组、灯影组、龙王庙组、西王庙组和二道水组等碳酸盐岩地层，地下水为溶洞暗河发育的裂隙层间水。安宁河河谷平原、邛海平原地下水资源丰富，为多层含水结构。一级阶地河漫滩砾卵石主要赋存潜水，主要为大气降水和河水补给。周边的山前洪积扇赋存承压水，水量较少。安宁河岸两岸山区的砂岩、粉砂岩、砾岩、泥岩、花岗岩和闪长岩等岩层中主要赋存碎隙岩类孔隙、裂隙水类型。盆地高地地区、台地漫滩地区，主要靠大气降水补给，雨季水位上升，旱季相反。基岩山区由于裂隙发育较少，无地下水位，水流量小，雨季、旱季流量差别大。

2 采样、分析测试方法

根据研究区地下水分布情况，设置15 个地下水露头作为采样点（图1）。水样的采集以及保存方法按照生活地下水标准检验方法（GB/T 14848-2017）。使用1.5 L 的聚乙烯瓶的容器采集水样，样品采集前，使用待采水样清洗3～4 次，再取样，样品装满不留气泡并密封，阳离子过滤酸化保存，阴离子原样过滤保存，阴阳共取30 件样品。样品编号、类型、水体状况、位置、水文地质分组如表1 所示。

表1 水样品概况统计表Tab.1 General situation of water samples

采样现场使用ISE 测定仪检测pH 值，重量法检测溶解性总固体。其他样品指标送予四川省地质矿产勘察开发局西昌地矿检测中心检验，检验方法如下：As、Se、Sb、Hg 含量使用AFS 方法检测；Sr、Cr、Cd、Co、Cu、Pb、Zn、Mn、Ni、Mo 含量使用MS 方法检测；K+、Na+、Ca2+、Mg2+、Al3+含量使用ICP 方法检测；Fe3+含量使用离子色谱法检测；含量使用纳氏试剂比色法检测；F-、Cl-、含量使用离子色谱法检测；偏硅酸含量使用比色法检测；含量使用滴定法检测；总硬度采用EDTA-2Na滴定法。

通过阴阳离子电荷平衡法检测地下水样品水化学分析结果可靠性，若无机离子平衡常数（NICB）小于5%，表明阴阳离子平衡数据可信。文中地下水样品的NICB 值为-4.7 %～4.9 %，平均值为-1.9%，测试分析数据可靠。

3 研究方法

3.1 水化学特征研究方法Piper

三线图能区分研究区地下水化学类型。Gibbs 图解法（Gibbs，1970）研究水岩作用对水化学成分的影响，将主要离子来源分为蒸发结晶作用、岩石风化作用、大气降水3 种类型。运用离子比例系数法不同的岩性对地下水化学成分的影响，结果能进一步反映水化学离子的来源（孙厚云等，2018）。

3.2 BP 神经网络

BP（Back Propagation）神经网络是机器学习中模拟生物神经网络进行学习的一种神经元连接模型，是一种单向传播的多层前馈神经网络，其主要特点是信号前向传播，误差反向传播（徐学良等，2017），以此往复拟合，通过一定规则输出结果。在水质评价研究中，水质指标多、差异大且关系复杂，受到多种因素的影响，以此构成一个典型非线性系统。通过训练学习对水质指标种类划分进行非线性拟合，从而凭借水质指标种类划分的内在规律进行评价（Abhijit et al.，2008）。

文中的BP 神经网络通过Python 的Pytorch 库搭建，并完成训练、学习以及评价。使用的激活函数为Sigmoid 函数，在特征相差比较复杂或是相差不大时效果比较好，并且函数整体平滑易于求导，通过对神经元加入激活函数能够增加神经网络模型的非线性，从而使模型能够更好的拟合非线性数据。其函数式如公式（1）所示。

文中运用的损失函数（Loss function）为Pytorch 中的交叉熵函数，如公式（2）所示：

式中：outputs 代表模型计算后的输出结果；targets 为样本标签；w为权重；n为张量维度；C为类别的数量。该函数在计算前将数据放入sigmoid 函数中使数据中间值更为敏感，体现出更高的不确定性，并且梯度下降时，可以避免均方误差损失函数学习速率下降的问题。在机器学习中通过损失函数对模型正向传播的输出与标签进行对比计算得到两者之间的误差值，从而确定反向传播的误差值，并能反映模型运行效果。

3.3 RMSprop 优化算法

BP 神经网络中的优化算法采用RMSprop 算法代替随机梯度下降算法（SDG），RMSprop 优化算法也称为均方根传递算法，其优点是能够加快梯度下降的速度以及有效减缓训练中损失曲线的山谷震荡以及鞍部停滞问题。在梯度下降过程开始，神经网络会从一个随机点开始，以此赋予每个属性的权重和偏置一个随机值，将该随机值计算的预测结果与标签对比，通过损失函数计算两者之间的误差即损失，再通过反向传播更新权重与偏置，直到预测结果接近标签值。这一过程就如同从山顶到山谷，山顶为最高损失，山谷为最低损失，导数为坡度。为到达山谷，每次求导都需要走下坡的道路，对于SDG 来说下山的方向是随机的，而RMSprop 算法会积累之前下坡的方向来决定下一个迭代下坡方向，从而优化SDG 算法。其计算过程为包括以下3 步：①计算每个参数在当前位置的梯度；公式（3）中wi为权重，bi为偏置，L(x)为损失函数。②计算更新量，通过对当前梯度计算权重均方根，以及偏置均方根；公式（4）中S为权重均方根，α为常数。③更新参数，公式（5）中 η为学习率，计算下一迭代中的权重与偏置，β 为防止S dwi为0 的极小常数。

从RMSprop 计算过程中，该算法计算更新量公式比一般梯度下降算法增加了一个常数 α来控制历史信息获取量，在设定全局学习率后，全局学习率在每次迭代中都会随衰减系数控制的历史梯度平方和而改变，从而使迭代方向在参数空间中更加平稳且快速。

以上3 个函数的关系如图2 所示，输入数据先通过Sigmoid 函数计算，再进行线性连接，其中乘以权重（w1，w2，w3）加上偏置（b1，b2，b3）输出第一隐层（H1），再通过同样的操作输出第二隐层（H2），通过线性层得到输出。将输出与标签利用损失函数比较得到误差，通过误差计算优化算法，计算结果对权重与偏重进行更新。

图2 模型函数关系图Fig.2 Model function relationship diagram

4 结果与讨论

4.1 研究区水化学特征

根据研究区地下水采集样品分析结果，地下水化学指标概况见表2。分析结果显示，pH 值平均为7.68，偏中性；TDS 值为6.7～230.4 mg/L，均值为85.91 mg/L；水离子阳离子浓度排序为Ca2+＞Mg2+＞Na+＞K+，阴离子浓度排序为，阳离子中Ca2+为优势离子，阴离子中为优势离子；As 含量为1.77～2.58 µg/L；Cr 含量为0.1～5.58 µg/L，均值为1.92 µg/L；Pb 含量为0.76 µg/L，最大值为10.26 µg/L；Mn 含量为0.11～6.33 µg/L，标准偏差为1.95；Ni 含量为0.05～9.54 µg/L，标准偏差为3.02，波动较大（表2）。

4.2 研究区水化学特征

Piper 三线图（图3）显示，研究区水化学类型主要是Mg2+·Ca2+-类型，阳离子主要分布在 Mg2+-Ca2+线上，分布于三角图左下区域，阴离子主要分布在线上，主要集中于三角图的左下角。

图3 研究区地下水水化学Piper 图Fig.3 Piper diagram of hydrochemical of groundwater in the study area

将研究区的水样数据绘制于Gibbs 图（Gibbs，1970），研究区水样基本落在Gibbs 的回旋镖内（图4），TDS 值约为100 mg/L，Na+/（Na++Ca2+）值小于0.7，值小于0.2，表明其受到人类活动的影响较少；水样主要分布于左侧以及中部偏下的位置，表明水化学离子组成总体受到岩石风化作用控制，大气降水也对其有一定程度的控制，但是没有岩石风化控制显著，而蒸发结晶作用微弱（王慧玮等，2021）。其阳离子在B2 区、D 区受到岩石风化作用控制；B2区以碎屑岩类泥岩砂岩为主，其化学离子主要来自于硅酸盐矿物；D 区以碳酸盐岩，Ca2+和会偏多。水量中等的B1 区、C 区水样在图中向右下角靠近，受到岩石风化作用和大气降水作用的共同控制，其岩石

图4 研究区地下水Gibbs 图Fig.4 Gibbs diagram of underground water in the study area

图5 地下水离子比值图Fig.5 Rates of the selected ions of groundwater

风化作用更显著。C 区以岩浆岩为主，C 区水样化学离子会来自于硅酸盐矿物的风化作用。

结合三大类造岩矿物（碳酸盐矿物、硅酸盐矿物、蒸发盐矿物）风化溶滤作用特征与离子比例关系图（图6）分析（刘永林等，2016；高旭波等，2020；杨芬等，2021），可以看出研究区的水样品大多都分布在硅酸盐矿物到碳酸盐矿物之间，相比离硅酸盐矿物区域更近，表明研究区水化学离子同时受到了硅酸盐矿物和碳酸盐矿物的风化溶滤作用的控制。在B1 型地下水区和B2 型地下水受硅酸盐矿物风化溶滤作用主要控制的同时，受到一定碳酸盐矿物风化溶滤作用控制，其中硅酸盐岩矿物来源于砂岩以及泥质岩；C 型地下水区主要受硅酸盐矿物风化溶滤作用控制，其来源于火山碎屑岩和花岗岩；溶洞暗河发育的D 型地下水区则主要受由碳酸盐矿物风化溶滤作用控制，其碳酸盐岩矿物主要来自石灰岩和白云岩。

图6 研究区水化学离子与矿物风化作用关系图Fig.6 Correlation of hydrochemical ions and mineral weathering

4.3 神经网络水质评价过程

4.3.1 神经网络搭建

笔者基于Python 语言使用Pytorch 开源库搭建BP 神经网络架构，搭建的BP 神经网络模型分4 层结构，分别为输入层，两层隐层以及输出层（图7）。隐层激活函数为Sigmoid 函数，输入输出神经元结构为：输入神经单元18 个，分别为水样品的18 个指标依次是As、Cr、Co、Cu、Pb、Zn、Mn、Ni、Mo、Na+、、Al3+、F-、Cl-、、TDS、总硬度，其他指标由于分异性不大或者低于检测下限，因此剔除；输出层采用线性链接5 个神经元的输出层结果用0、1 表示，水质分五级，Ⅰ类水为 [1，0，0，0，0]，Ⅱ类水为 [0，1，0，0，0]，Ⅲ类水为 [0，0，1，0，0]，Ⅳ类水为 [0，0，0，1，0]，Ⅴ类水为 [0，0，0，0，1]输出，输出表的5 个单元表分别代表Ⅰ-Ⅴ级，若满足则为1，不满足则为0。

图7 BP 神经网络结构图Fig.7 Neural network structure diagram

4.3.2 训练样本设计

训练样本数据依据国家地下水水质标准（GB/T 14848-2017）中各类水质的指标在其各范围内使用计算机插值生成并添加标签（邓大鹏等，2007；袁瑞强等，2021），样本按平均每种类型的样本占总样本的20%分布，共5 000 个训练样本，将样本数据取95%作为训练集，剩余5%作为测试集，以供模型训练过程中对损失函数与预测正确率进行测试。

4.3.3 神经网络训练

训练开始首先读取生成的训练样本为Numpy 格式后进行数据标准化处理，接着将处理后的数据转换成张量（Tensor）数据格式，初始化神经网络后设定基本参数，隐层个数设定在17～9 之间，学习率（lr）设定为0.001～0.000 1，训练轮数（Epochs）为500 轮以保证监测最后训练结果的稳定性，每轮模型迭代次数（Step）为5 000 次，期望误差为1×10-2。在训练过程中通过训练集和测试集的损失函数以及测试集正确率图像监测网络模型训练效果。通过代码自动循环参数比较测试结果，最终确定最优的网络模型。参数为第一隐层14 个神经元，第二隐层9 个神经元，学习率0.000 8，图8 为该参数下训练效果图。在训练结束保存训练好的训练模型，在进行研究区水样评价时调用。

图8 训练效果图Fig.8 Training effects

图中训练集损失在1×10-2收敛，收敛过程平滑，测试集损失低于训练集损失，在2×10-2收敛，测试集正确率稳定在0.988。表明此时的网络模型拟合良好并且泛化能力强，没有出现过拟合现象，在训练过程中学习充足，能够对正常水质样品进行准确且客观的预测。模型整体以及张量数据运算采用GPU 运算，GPU 型号为RTX3080Ti，整体训练时间大幅减少共费时21.7 分钟。

4.4 研究区水质评价结果

调用训练完成后的神经网络模型，将研究区水样品数据进行标准化后导入得到水质评价结果（表3），Ⅰ类水质点2 个占13.3%分别为5-SY、12-SY；Ⅱ类水质点6 个占40%分别为2-SY、4-SY、6-SY、10-SY、11-SY、13-SY；Ⅲ类水质点为7 个占46.6%分别为1-SY、3-SY、7-SY、8-SY、9-SY、14-SY、15-SY。Ⅱ、Ⅲ类水质占多数从研究区数据来看As 含量大多数达到了Ⅲ类标准，含量达到了第Ⅲ类或第Ⅱ类标准。从神经网络训练的权值来看，数据网络的评价主要考虑As、Cr、Mo、、Al3+、F-、以及总硬度等8 个水质指标，表明As、为造成该地区地下水质达到Ⅲ类的指标的重要原因之一。BP 神经网络具有很强的自主学习、自组织、自适应能力，充分学习了训练样本的水质特征，建立起水质指标与水质等级的非线性对应关系。并且权重与偏置都是通过学习得到而非人为给定，很大程度上避免了主观因素影响，从而使其评价结果更加客观合理，然而美中不足之处是缺少真实的训练样本以供模型学习，使预测结果更加贴近现实，并且神经网络评价很难确定主要影响因子。

表3 研究区水质综合排名表Tab.3 Comprehensive ranking of water quality in the study area

结合水样品数据的水文地质背景，以及样品指标在国标中超标部分确定影响因子，对研究区各点水质进行综合评价并由好到坏进行排名。

4.5 地下水安全建议

对于分布在德昌县阿月乡光辉村西番箐、罗家坪子；西昌市安哈镇摆摆顶村、黄水乡观音岩、西溪乡长板桥村；普格县五道菁乡黄坪村五组的Ⅰ类和Ⅱ类地下水地区，水质良好污染少，饲养牲畜能够正常饮用，居民应该经简单处理后饮用，并建议当地政府部门在此建立自来水饮水站。

对于分布在西昌市黄水乡新塘村大庆沟、洼垴村七组、书夫村二组、中坝乡小浸沟、黄联关镇哈土村四组；普格县特尔果乡甲甲沟村、普格县特补乡白庙子等Ⅲ类地下水地区，水质差又受到一定程度污染，牲畜饮用该水需经处理，不建议当地人直接或简单处理后饮用。尤其是As、Pb、Cr 超标的地下水，长期大量饮用会对人体造成伤害，建议：①地方政府及相关部门在Ⅲ类水质地区应该高度重视地下水质是否存在危害人体的元素超标问题。②相关部门以及居民应该加强井水和地下水保护，隔绝农业活动对井水带来的污染，如As、等是大部分农药化肥中带有的元素，应当探明农药化肥是否已对当地地下水造成污染。③建议当地政府及相关部门加强对研究区大气、汇水区及地下水水质监测，查明Ⅲ类地下水污染来源并阻断以保障当地居民的安全生活生产。

5 结论

（1）螺髻山北麓地区地下水水化学结构主要是Mg2+·Ca2+-HCO-3类型，TDS 较低受人类活动干扰较少，岩石风化作用对水化学离子组成控制显著，其次是大气降水。阴阳离子主要是来自于硅酸盐矿物与碳酸盐矿物共同风化溶滤作用。其中硅酸盐矿物主要有花岗岩、长岩、灰岩、砂岩、页岩、泥岩等岩石；碳酸盐矿物主要有泥灰岩、白云岩、泥质灰岩等岩石。

（2）阴阳离子相关性分析和比值分析结果表明，方解石的风化是Ca2+的来源之一，Mg2+主要受到硅酸盐矿物溶解控制，与 Cl-大部分都来自于碳酸盐矿物溶解。

（3）利用BP 神经网络对研究区水质样本进行评价，水质评价结果可为研究区水质分类资源化利用，人民生活生产用水保障提供参考。评价结果表明总体水质较好，其中Ⅰ类水质点占13.3%，Ⅱ类水质点占40%，Ⅲ类水质点为7 个占46.6%。