0 引言

地下水位作为区域水资源管理的重要依据，对水资源的合理开发及可持续利用有着重要意义；但由于受一系列自然条件与人为活动影响，地下水水位在时间序列上呈现出高度的随机性与滞后性。随着生态文明建设的推进，探索出简单有效的地下水位预测方法迫在眉睫。传统的地下水水位的预测方法往往采用简单的线性函数描述其动态变化特征，忽略了外部条件的影响，结果往往不尽人意。随着科学技术的发展，地下水水位预测方法逐渐完善，一些学者将一些非线性算法、模型应用于地下水水位预测研究：刘博等^[1]利用SOM-RBFN模型对吉林市地下水位进行预测，预测结果较为理想；马新双等^[2]采用前馈神经网络模型对北方浅层地下水水位情况进行预测，获得了较高预测精度；周振民等^[3]将灰色模型(grey model，GM(1, 1))应用到滦河下游地区地下水位预测中，结果表明该模型预测精度较高，对地下水的合理开发利用具有一定的实用价值。但这些模型的建立往往需要大量高质量的数据作为支撑且调参较为复杂，故利用率较低。同时，时间序列模型也被广泛应用于地下水水位预测：Ahn^[4]、陈汉军等^[5]将时间序列模型用于地下水水位估算中，均得到了较为理想的结果；张展羽等^[6]建立主成分-时间序列模型对济南市陡沟灌区地下水水位进行预测，结果表明，其预测精度要高于逐步回归模型与BP神经网络模型。但是由于传统时间序列模型过于突出时间序列，忽略了外部条件影响，当遇到外界因素发生较大变化时往往会产生较大偏差，因此有较大的预测缺陷。

近些年，随着深度学习技术的高速发展，深度学习模型逐渐被应用于时间序列数据的研究中。其中，循环神经网络(recurrent neural network，RNN)将时间序列引入神经网络中，使模型在时序数据分析上表现出更强的适应性。但基于结构特点，循环神经网络仍存在一些问题，其中较为严重的是容易出现梯度消失或者梯度爆炸的问题，即时间过长而造成记忆值较小的现象。在众多循环神经网络变体中，长短期记忆(long short-term memory，LSTM)神经网络通过遗忘门(forget gate)和更新门(update gate)来调控梯度，将反向传播累乘作用变成累加作用，解决了梯度消失或爆炸的问题，使得循环神经网络能够更加有效地处理分析长距离时序信息。LSTM神经网络在不同领域的时序数据处理分析方面已有较为成功的案例，Yuting等^[7]将LSTM神经网络应用于海平面温度的预测、Tax等^[8]将其用于业务流程监控的预测均取得了不错的结果。因此，为克服传统时间序列模型过于突出时间序列的缺陷，同时考虑外部条件对地下水水位的影响，本文将LSTM神经网络应用于地下水水位预测研究，并引入多个影响变量，以期获得合理有效的预测结果。

1 研究方法及相关理论 1.1 RNN循环神经网络

对于传统神经网络而言，从输入层到隐含层，再到输出层，层与层之间为全连接，但是每层之间的节点却是无连接的；因此传统神经网络在处理序列数据时，结果往往有很大的偏差。不同于传统神经网络，RNN循环神经网络会对之前的信息进行记忆并将其应用到当前的输出计算中，即隐含层之间的节点为连接的，并且隐含层的输入包括输入层的输出与上一时刻隐含层的输出。令x为输入单元，s为隐含层单元，o为输出单元，L^*为损失函数，y为训练集标签，t为时刻，V、W、U分为单元共享权重系数矩阵，则RNN前向输出如图 1所示。

故对于时刻t，隐含状态为

式中：f为激活函数(比如tanh函数)；x_t为t时刻输入；b_s为隐含层单元偏置向量，b为偏置向量，下标为对应单元。

输出量为

式中，soft max为归一化指数函数。

损失函数为

基于处理时序数据的特点，RNN常见训练方法为BPTT(back-propagation through time)算法，即梯度下降法，其特点为沿需要优化参数(V、W、U)的负梯度方向不断寻找更优的点直至收敛。但在寻优过程中，求解某时刻参数的偏导需追溯该时刻之前所有时刻的信息，整体偏导函数为所有时刻的累加; 嵌套激活函数后，部分累乘会导致激活函数(tanh或sigmoid函数)导数的累乘，进而导致“梯度消失”或“梯度爆炸”。^[9]。

1.2 多变量LSTM神经网络模型

为解决RNN“梯度消失”或“梯度爆炸”的问题，Hochreiter等^[10]于1997年提出了LSTM，通过门控制将短期与长期记忆结合起来，从而在一定程度上解决了上述问题。令Γ_f、Γ_i、c、Γ_o、分别表示遗忘门、输入门、单元状态、输出门、新记忆，权重系数矩阵下标f、i、o分别表示遗忘门、输入门、输出门，h_t表示t时刻隐含状态，σ表示激活函数sigmoid，则LSTM结构可用图 2表示。

LSTM网络通过3个“门”来实现信息的保护和控制：

1) 遗忘门，决定单元c_t-1需要丢弃或保留什么信息，通过读取h_t-1与x_t，输出[0, 1]赋值给c_t-1，其中“1”表示完全保留，“0”表示完全丢弃。遗忘门输出为

2) 输入门，决定着让多少信息加入到单元中，通过sigmoid的信息决定与tanh的信息生成，结合遗忘门联合对单元状态进行更新。输入门步骤为：

3) 输出门，决定着当前单元状态哪一部分信息被输出，依然通过sigmoid与tanh完成。输出门步骤^[11]为：

基于时间序列数据特点以及循环神经网络从简设计的原则，多变量LSTM预测模型框架如图 3所示，该框架大致可分为输入层、隐含层、模型训练、输出层4个部分，其中X=[X₁, X₂, …，X_L]为实际输入，P=[P₁, P₂, …，P_L]为实际输出。各模块功能如下。

输入层：对原始变量时间序列集进行分割、标准化处理，以满足网络输入要求。

隐含层：采用图 2所示的LSTM单元结构对参数进行权重更新、优化。

输出层：输出预测结果、反标准化处理、验证误差。

模型训练：模型采用Adam优化算法^[12-13]进行网络权重更新，步骤为：

1) 初始化参数向量θ₀、一阶矩向量λ₀=0、二阶矩向量v₀=0、时间步t=0。

2) 当参数θ没有收敛时，令t=t+1，循环迭代更新各部分。

3) 更新随机目标函数在时间t上对参数θ所求梯度：

式中：g_t为时间t上的目标函数梯度；▽_θ为关于θ的偏导；f_t(θ_t-1)为时间t-1上的随机目标函数。更新偏差的一阶矩估计λ_t和二阶原始矩估计v_t：

式中，β₁、β₂为指数衰减速率。

4) 计算偏差修正的一阶矩估计、二阶矩估计：

5) 更新模型参数θ：

式中：α为步长；ε为常数。

2 研究区概况与数据来源 2.1 研究区概况

泰安市位于山东省中部，地处116°04′E—118°00′E，35°38′N—36°29′N，总面积7 762 km²。年平均气温12.8 ℃，多年平均降水量为680 mm，多集中在6—9月，占全年降水量的75%，多年平均蒸发量1 813 mm，属暖温带半湿润季风气候。境内水系主要属大汶河水系和淮河水系，地形东高西低，为侵蚀或剥蚀构造地形，海拔37.5~1 532.0 m，受地形地貌影响，东部降水多于西部，山区降水量多于平原，总趋势为自东北向西南逐渐减少。区内地下水类型以裸露岩溶孔洞水为主，含水层厚度为70.0~80.0 m，富水性强，单井涌水量达1 000~5 000 m³/d。

2.2 研究数据来源

本研究所选长观井为泰安市岱岳区的J1号井，分析该井地下水位主要受气象因素和人工开采影响，因此选取2001—2016年逐月平均地下水位观测数据和该地蒸发量、降雨量、气温、气压、相对湿度、日照时长等6项气象数据，以及人工开采量数据^[14-15]，共192组时序数据作为研究数据。气象数据则根据中国气象局地面气候资料泰安气象站(站号54827)日值数据集统计而来，其中2007年7—12月地下水水位监测数据有缺失情况，为模型训练需要，取每年该月平均值进行填补，历年地下水位变化如图 4所示，地下水开采量见表 1。

3 研究区地下水水位预测 3.1 实验准备

为构建基于多变量LSTM神经网络的地下水水位预测模型，选取2001—2013年数据作为模型训练集，2014年数据作为测试集，2015—2016年数据作为预测集，以地下水位、蒸发量、降雨量、气温、气压、相对湿度、日照时长、开采量作为模型变量输入，以测试集的均方根误差(root mean squared error，RMSE)作为模型训练预测精度评价指标。

为消除量纲及数量级的影响，提高模型预测精度，加快模型训练收敛速度，本文采用以下方法对原始变量时序集进行标准化处理，公式可表示为

式中：x^*表示标准化后的值；x表示标准化前的值；x表示原始数据集每个变量的平均值；x_δ表示原始数据集中每个变量的标准差。以地下水月均开采量为例，选取该变量时序数据中平均值为5 075.081 3，标准差为545.193 2，计算得到地下水月均开采量标准化后数值(表 2)。

3.2 模型训练

根据LSTM模型结构，多变量LSTM神经网络预测模型训练步骤如下：

1) 在输入层中，定义2001—2014年数据集为

式中：x_t={a_t, b_t, c_t, d_t, e_t, f_t, g_t, h_t}, a_t, b_t, c_t, d_t, e_t, f_t, g_t, h_t分别表示t时刻对应的地下水位、蒸发量、降雨量、气温、气压、相对湿度、日照时长、开采量。

2) 划分标准化处理后的训练集与测试集，可表示为：

式中：F_train^′，F_test^′分别表示标准化后的训练集与测试集；训练集截止点m=156，测试集截止点n=168。

3) 对数据F_train^′进行分割，设分割窗口长度为L，L=12，故分割后输入隐含层的数据集为：

式中，X为输入到隐含层L个细胞单元的数据集。

4) 隐含层理论输出为：

实际输出为：

式中，LSTM为式(4)—(9)的计算过程。

5) 训练过程中均方误差(E_ms)作为损失计算公式^[16]：

式中：y_i为理论输出；p_i为实际输出。

6) 以损失最小为目标，应用Adam优化算法不断更新网络权重，直至网络损失最小。

3.3 模型调参

在模型训练的同时，需对模型进行调参，调参的部分主要有批次大小、隐含层节点数和学习率的选取。由于本次研究所取的数据集较小，故可采用全数据集形式，其优点为由全数据集所计算的梯度可以更好地代表数据集总体，从而更准确地朝向极值所在的方向。对于隐含层节点数以及学习率的调参，本研究则采取控制变量的方法确定最优参数。首先固定非关键参数取值，选择批次大小为156，确定重复训练次数为2 000；其次确定重要参数的范围以简化调参过程，隐含层节点数经验公式如下式所示，学习率一般按3的倍数来调整^[17]。

式中：p为隐含层节点数；q为输入层节点数；l为输出层节点数；a∈[1,10], a∈N。因此，隐含层节点数p∈[11,20], p∈N。令学习率η={0.001, 0.003，0.010, 0.030, 0.100, 0.300}，在非关键参数条件下，依次选择不同的隐含层节点数，判断不同的学习率对测试集RMSE的影响，从而确定最优参数。其中，多变量LSTM模型相关组合的RMSE如表 3所示。

由表 3可得，当η=0.001、p=19时，测试集RMSE最小，为0. 770 4，因此选取该参数对模型训练。

3.4 模型预测

模型采用单点预测的方式，确定参数后，将2001—2013年时序数据输入模型中训练。为预测2015—2016年地下水位，设定预测次数为24，理论输出的最后L项数据为

将x_m^′输入LSTM模型中，预测输出结果为

将n+1时刻模型变量输入模型中，得到P_n+2^′，依此类推，得到预测时序集为

对P_*^′进行反标准化处理得到实际预测序列集P_*={P_n+1, P_n+2, …, P_n+24}，并计算其与预测集F_FCT的偏差以定量解释模型的预测精度，整个训练、预测过程在python环境中进行。

3.5 实验结果分析

为验证多变量LSTM神经网络预测地下水水位的优点，使用单变量LSTM神经网络和BP神经网络对地下水水位进行预测，并作对比分析。其中，BP神经网络以7个影响因子的时序数据作为输入值、地下水位时序数据作为输出值训练，采用回归方式进行预测。已有大量学者对BP神经网络^[18]进行过研究，故本文对BP神经网络的训练、预测过程不再赘述。单变量LSTM模型的数据集的构建过程与多变量LSTM时间序列数据集的构建过程大致相同，区别在于多变量LSTM时间序列每一个时刻的数据不是一个数值，而是一个多维向量，向量的维度由变量的个数决定，本研究单变量LSTM模型仅以地下水水位作为变量输入。三者预测集及实测水位拟合曲线如图 5所示，测试集、预测集拟合RMSE见表 4。

由表 4可知，单、多变量LSTM模型以及BP模型测试集的RMSE低于预测集；这是因为模型调参过程是一个寻找靶心的过程，选择测试集误差最低的参数作为最优参数，说明靶心对于测试集而言处于最佳位置，而预测新数据时，原有靶心则会偏离，导致预测精度下降。由图 5可知，BP拟合曲线较实测值偏离较大，多变量LSTM模型预测曲线与实测值较为接近，单变量LSTM模型的预测曲线呈现出较为平稳的周期性变化，这是由于：BP模型在训练过程中未考虑时序变化规律，仅考虑了相关变量的影响，总体预测精度要显著低于多变量和单变量LSTM模型；多变量LSTM考虑地下水位时间序列变化前后依赖关系的同时，引入了其他相关影响变量的相互作用，故预测效果较好；单变量LSTM模型仅以地下水位作为变量进行训练，根据地下水位的时序变化规律作出预测，因此无法考虑相关变量的影响，在实测值波动较大时误差较大，所以预测精度介于二者之间。综上可以认为，在三种模型中，基于多变量的LSTM预测模型拟合曲线拟合最好，仅在某些峰值处误差较大，但总体预测效果较为理想。

4 结论与建议

1) 基于多变量的LSTM神经网络模型可以作为预测地下水水位动态变化简单有效的工具，与单变量LSTM神经网络模型和BP神经网络模型相比，拥有更高的预测精度。

2) 基于多变量的LSTM神经网络模型预测精度对参数设置依赖较大，不同的参数组合对模型预测精度影响较为显著。本文通过确定非关键影响参数，对重要参数进行组合调参，相比于全参数调参，在弱化模型训练复杂程度的同时确保了预测精度。

3) 多变量LSTM神经网络预测模型在拟合曲线个别峰值段存在误差偏大的现象，因此如何提高模型在地下水位突变处的预测精度将是日后研究的重点，并且在今后的工作中将引入不同地区监测的资料以验证模型的适应能力。