萊州灣南岸海水入侵特征及基于WA-NARX 混合模型的動態(tài)預測

董凡，張曉影，陳廣泉,，戴振學，王延誠, 3

(1.吉林大學 建設工程學院，吉林 長春 130026；2.自然資源部第一海洋研究所 自然資源部海岸帶科學與綜合管理重點實驗室，山東 青島 266061；3.青島海洋科學與技術試點國家實驗室 海洋地質過程與環(huán)境功能實驗室，山東 青島 266237)

1 引言

近年來，受氣候變化和人類活動的雙重影響，全球海岸地區(qū)都面臨著海水入侵的威脅。海水入侵可引發(fā)地下水水質惡化、土壤退化等生態(tài)環(huán)境問題，嚴重制約社會的可持續(xù)發(fā)展[1–4]?？茖W管理沿海地區(qū)地下水，預測近岸地下水位和電導率，掌握地下水位及鹽度的動態(tài)變化，對合理開發(fā)利用海岸帶地下水資源，防止海水入侵加劇，保護海岸帶生態(tài)環(huán)境有十分重要的指導意義。然而，受自然條件和人為因素的影響，海水入侵在時間序列上表現(xiàn)出時空性、不確定性及動態(tài)變化等特性，多種因素對海水入侵的影響并非多因子影響的線性疊加，從而使預測海水入侵的演化過程變得十分困難。

隨著人工智能的興起，由于人工神經(jīng)網(wǎng)絡（ANN）可以較好地映射非線性變化，因此被廣泛應用在水文變量預測中。其中BP 神經(jīng)網(wǎng)絡模型的應用范圍最廣[5]，但其神經(jīng)網(wǎng)絡訓練的時間長，收斂速度慢，輸入信號由輸入層單向傳遞給輸出層，網(wǎng)絡模型只能使用有限步的歷史信息。而循環(huán)神經(jīng)網(wǎng)絡（RNN）憑借其循環(huán)結構的特點，通過反饋連接可以使用長期的歷史信息，因此更加適合時間序列的預測[6]。具有外部輸入的非線性自回歸神經(jīng)網(wǎng)絡（Nonlinear Autoregressive Models with Exogenous Inputs，NARX）作為循環(huán)神經(jīng)網(wǎng)絡的一種，在地下水位的預測研究領域備受青睞。Guzman 等[7]使用NARX 模型預測密西西比河流域含水層，Wunsch 等[8]運用NARX 模型預測德國西南部數(shù)口井的地下水水位，Di Nunno 和Granata[9]應用NARX網(wǎng)絡預測位于普利亞地區(qū)的76 口監(jiān)測井的日地下水水位。研究結果表明，NARX 神經(jīng)網(wǎng)絡在預測時間序列上具有顯著優(yōu)勢。

盡管ANN 在時間序列預測方面應用廣泛，具有靈活性，但是如果不對輸入和輸出數(shù)據(jù)進行預處理，它可能無法處理非平穩(wěn)數(shù)據(jù)[10–11]。小波變換作為一種數(shù)據(jù)預處理技術，可以將非平穩(wěn)時間序列的小波分解成不同的尺度，提供對序列結構的解釋，并僅用幾個系數(shù)來提取有關其歷史的重要信息。該技術被廣泛應用于非平穩(wěn)信號的時間序列分析[12]，通過將非平穩(wěn)信號分解成不同時間水平的子信號，可以較好地預測非平穩(wěn)信號[13–18]。將小波變換和人工神經(jīng)網(wǎng)絡模型相結合，發(fā)揮這兩種模型的優(yōu)勢，可以明顯提高預測精度，增加有效預測時段。目前，混合模型已廣泛應用于部分水文要素的預測，如降雨、徑流和蒸發(fā)[19]。張建峰等[20]建立小波–BP 神經(jīng)網(wǎng)絡的混合模型預測華北平原城市地下水位。Seo 等[21]將小波分析（WA）分別結合BP 神經(jīng)網(wǎng)絡、自適應神經(jīng)模糊推理系統(tǒng)，預測水庫水位并對結果進行對比。Zhang 等[22]構建了小波分析和NARX 神經(jīng)網(wǎng)絡混合模型對舟山填海區(qū)地下水位進行預測，預測精度較高。

由于受全新世海侵、氣候干旱、海平面上升以及不合理的人類活動的影響，萊州灣地區(qū)是我國海水入侵災害最為嚴重和典型的地區(qū)[23]。近年來，在上述因素綜合影響下，萊州灣地區(qū)地下咸水體分布面積已超過22 000 km2[24]，嚴重污染地下淡水資源，土壤次生鹽漬化面積不斷擴大，對海岸帶區(qū)域的工農業(yè)生產和生態(tài)環(huán)境造成嚴重影響，海水入侵已成為該地區(qū)嚴重的環(huán)境問題之一，嚴重制約經(jīng)濟社會發(fā)展。因此，開展萊州灣地區(qū)海水入侵動態(tài)特征分析，建立預測模型，對提高海岸帶地區(qū)防災減災能力和經(jīng)濟、社會可持續(xù)發(fā)展具有重要意義。在本研究選擇萊州灣南岸典型監(jiān)測剖面，首先分析地下水位和地下水電導率的動態(tài)變化特征，分析對比不同輸入?yún)?shù)對模型預測精度的影響，建立基于WA-NARX 混合的海水入侵預測模型，對近岸地下水位和電導率進行預測模擬，研究結果可為海岸帶海水入侵的實時監(jiān)測、預測、預警提供理論與技術支撐。

2 研究區(qū)與研究方法

2.1 研究區(qū)與資料

2.1.1 研究區(qū)概況

研究區(qū)位于萊州灣南部（圖1），濰坊北部的濱海區(qū)。研究區(qū)內設有6 口監(jiān)測井（BH01–BH06），分別沿垂直于海岸帶的剖面由靠近海岸帶的監(jiān)測井BH06 向內陸的監(jiān)測井BH01 方向分布。研究收集了降雨、蒸發(fā)、潮汐、電導率、地下水位和農業(yè)排灌用電量數(shù)據(jù)，見圖2。其中，潮汐的變化是潮位站監(jiān)測的海水位每小時的變化，降雨量與蒸發(fā)量是日測數(shù)據(jù)，排灌用電量是研究區(qū)每小時灌溉所用電量。研究選 取2011 年10 月1 日 至2012 年10 月31 之 間 的 監(jiān)測數(shù)據(jù)作為研究對象。

圖1 研究區(qū)位置圖及監(jiān)測井分布Fig.1 Location of the study area and monitoring wells

圖2 各影響因素時間序列Fig.2 Time series of influencing factors

2.1.2 研究區(qū)水文地質特征

萊州灣南岸地區(qū)屬于魯北平原水文地質區(qū)、濰北濱海平原水文地質亞區(qū)。該地區(qū)地表水系較發(fā)達，地下水水文地質特征在水平和垂向上均變化明顯。地貌類型由南向北依次為山前沖洪積平原、沖積海積平原和濱海海積平原，海岸類型為粉砂淤泥質海岸。此處的沖洪積平原由多個扇緣相互交接的沖洪積扇群組成，沖洪積平原的前緣被海積物或黃河沖積物所掩埋，含水層厚度大、顆粒粗，具有多層結構。含水層顆粒自上而下呈現(xiàn)由粗變細的趨勢，而在水平方向上則具有沖洪積扇的水文地質特征[25]。

如圖3 所示，研究區(qū)普遍存在一層中更新統(tǒng)砂質黏土弱透水層，埋深一般為 40～70 m，層厚約20 m。有3 層海侵層，分別對應晚更新世以來發(fā)生的3 次大海侵。在3 次海侵中都滯留了濃縮于地層中的古海水（咸–鹵水），與淡水含水層在水平方向上相連，是該區(qū)咸–鹵水入侵的來源，主要出現(xiàn)在晚更新世和全新世地層中。

圖3 研究區(qū)剖面圖（修改自文獻[26]）Fig.3 Hydrogeology of the cross section from south to north in the study area (modified from reference [26])

2.2 研究方法

2.2.1 小波分析

在本研究中，運用小波分析的方法來識別顯著的振蕩周期，從原始數(shù)據(jù)中提取有用的信息，并將這些信息分解成一系列小波函數(shù)，小波函數(shù)由母波函數(shù)生成。本文采用的母波函數(shù)是Morlet 小波：

式中，w0為無量綱頻率；t為時間。通常，選擇w0=6的Morlet 小波[27–28]。連續(xù)小波變換（CWT）用于識別信號在時頻域的主振蕩。連續(xù)小波變換的f(t)如下所示：

式中，a為尺度或頻率參數(shù)；b為時間或位移參數(shù);(t) 是 復共軛函數(shù) ψ(t)。由式（2）可知，小波變換在不同分辨率或尺度下呈現(xiàn)出不規(guī)則和不對稱的特點，將f(t)分解為一系列低分辨率分量。

在實際應用中，連續(xù)小波變換常被離散小波變換（DWT）所代替。令a0>1,b0∈R,k和j為整數(shù)，則f(t)的離散小波變換可表示為

在式（3）中，參數(shù)a0和b0最常用的選擇是使用對數(shù)縮放[29]，即a0=2,b0=1。原始時間序列在頻率（a或j）和時間（b或k）域的特征被連續(xù)小波變換的W(a,b)或離散小波變換的W(j,k)所反映。

將經(jīng)典的連續(xù)小波轉換為離散小波有多種算法，如Mallat 算法和á Trous 算法。本文采用了á Trous 算法[30]。á Trous 的分解算法定義如下：

式中，g(n)為離散低通濾波器；Ci(t) 和Wi(t) (i=1,2,···)分別是分辨率為i級的尺度系數(shù)和小波系數(shù)。W1(t),W2(t),···,Wp(t) 和Cp(t)組成分辨率為P的離散小波變換。

在小波分析中經(jīng)常使用近似特征和細節(jié)分量，近似特征表示信號的低尺度，即低頻信息；細節(jié)分量表示信號的高尺度，即高頻信息。離散小波變換使用低通濾波器g(n)和高通濾波器h(n)將f(t)分解成不同尺度上的低頻成分（C1，C2，C3··· ）和高頻成分（W1，W2，W3···）。因此，f(t)可以分解為包含其近似分量（低頻成分）和細節(jié)分量（高頻成分）。高頻信號捕獲原始時間序列數(shù)據(jù)中的小尺度特征，見圖4 小波分解部分。

圖4 WA-NARX 模型建模流程Fig.4 Flow chart of WA-NARX neural network hybrid model

2.2.2 NARX 神經(jīng)網(wǎng)絡

NARX 神經(jīng)網(wǎng)絡是一種有外部輸入的非線性自回歸網(wǎng)絡模型[31]，表示為

式中，y(t)代表輸出信號；x(t)代表外部輸入信號；f為非線性函數(shù)。即輸出信號y(t)可由輸出信號前ny個值和外部輸入信號x(t)前nu個值通過f非線性映射得到。

NARX 神經(jīng)網(wǎng)絡由網(wǎng)絡的輸出反饋和動態(tài)的神經(jīng)元構成，是一種動態(tài)神經(jīng)網(wǎng)絡。它的結構清晰，通常包含輸入層、輸出層、隱含層及輸入延時和輸出延時這些網(wǎng)絡結構（圖5）。NARX 神經(jīng)網(wǎng)絡具有開環(huán)和閉環(huán)兩種模式，閉環(huán)網(wǎng)絡可以進行多步預測。與一般神經(jīng)網(wǎng)絡相比，NARX 神經(jīng)網(wǎng)絡的反饋結構和梯度下降算法使其擁有更好的收斂速度和學習效率，是一種研究較為廣泛的動態(tài)神經(jīng)網(wǎng)絡[32–33]，因此，NARX 神經(jīng)網(wǎng)絡既可以用于時間序列建模，也可以用于非線性濾波[34]。

圖5 輸入到輸出的NARX 神經(jīng)網(wǎng)絡結構Fig.5 Nonlinear auto-regression with exogenous input NARX

2.2.3 小波與NARX 混合模型

WA-NARX 混合模型是利用多分辨率小波分解得到的子時間序列作為輸入的NARX 模型。利用小波變換將原始信號依次分解成不同尺度上的信號成分，á Trous 算法產生低頻的近似分量（C）和高頻的細節(jié)分量（W）[35]。本研究利用L=int來確定小波的分解層數(shù)[36]，其中L表示小波層的個數(shù)，n表示時間序列中的數(shù)據(jù)點個數(shù)，int[ ]為取整函數(shù)。研究所包含時間序列的數(shù)據(jù)點個數(shù)為9 528，經(jīng)計算需要4 個分辨率級別。

本研究中，利用á Trous 算法首先將水位、潮汐和降雨等參數(shù)的小時時間序列分解為一定數(shù)量的子時間分量。根據(jù)潮汐水位隨時間的變化規(guī)律，運用小波分析的方法來識別顯著的振蕩周期，并從原始數(shù)據(jù)中提取潮汐水位變化的近似特征值和細節(jié)分量，然后將時間分量作為混合模型的輸入。選擇最重要的時間分量作為輸入，這些時間分量與原始數(shù)據(jù)具有很高的相關性。為了判別預測的準確性，本文將均方根誤差（RMSE）和擬合度（R2）作為性能指標來衡量模型的預測精度。同時，利用matlab 中的perform 函數(shù)結合網(wǎng)絡結構與均方誤差（MSE）計算網(wǎng)絡性能來評價網(wǎng)絡的學習及預測能力。

WA-NARX 混合模型在建模過程中，數(shù)據(jù)被分為兩個子集：訓練集和測試集。訓練數(shù)據(jù)集（75%）用于優(yōu)化網(wǎng)絡配置，測試數(shù)據(jù)集（25%）用于測試網(wǎng)絡性能。此外，WA-NARX 混合模型的隱含層神經(jīng)元數(shù)是通過反復試驗確定的。在本次研究中，最終確定10 個隱含層神經(jīng)元。WA-NARX 模型建模流程見圖4，如圖所示：首先是參數(shù)和樣本的選擇，利用小波變換將原始時間序列分解成不同信號成分的小波分量，然后將分解的小波分量輸入到NARX 神經(jīng)網(wǎng)絡模型中進行訓練預測，最終得到預測結果。

3 結果和討論

3.1 海水入侵動態(tài)變化特征

3.1.1 地下水位動態(tài)變化特征

根據(jù)研究區(qū)監(jiān)測井2011–2012 年間的地下水水位監(jiān)測數(shù)據(jù)，結合研究區(qū)內相應時段的降雨量繪制監(jiān)測點的地下水位動態(tài)變化曲線（圖6）。根據(jù)研究區(qū)2011–2012 年各監(jiān)測井采集的地下水位動態(tài)數(shù)據(jù)，計算各監(jiān)測井地下水位的最大值、最小值、平均值以及標準差等空間變化統(tǒng)計特征，計算結果見表1。圖6 顯示，各監(jiān)測井地下水位的變化趨勢基本一致。監(jiān)測井在2011–2012 年間的地下水位的變幅為1.175～4.358 m，各監(jiān)測井地下水位在2011 年12 月、2012 年4 月、2012 年8 月這3 個時間段內均有突變點，與降雨量變化相對應。在降雨后，多數(shù)監(jiān)測井地下水位有明顯上升的趨勢。2011 年的11 月份到2012年的2 月份水位逐漸上升后趨于平緩，因降雨補給存在滯后效應，隨后水位開始下降，2012 年3–4 月份水位隨著降雨上升，在4 月末水位達到全年最高值，因春灌開采影響，5 月份雖有降雨補給，但地下水位也開始下降，7–8 月是降雨最頻繁月份，受強降雨的影響，地下水位隨降雨量波動變化。在降雨過程中，間歇性的小雨只能濕潤土壤表層，不能有效補給地下水；集中式暴雨入滲系數(shù)低，大部分轉化為地表徑流。而連續(xù)小雨可以充分入滲到地下，因此有利于地下水補給[37]?？梢耘袛啾O(jiān)測井2011–2012 年地下水年內動態(tài)變化特征為降雨入滲–開采型。

圖6 監(jiān)測井2011–2012 年間地下水位變化曲線Fig.6 Groundwater level variation in monitoring wells from 2011 to 2012

表1 監(jiān)測井地下水位統(tǒng)計特征分析（單位：m）Table 1 Analysis on statistical characteristics of groundwater level in monitoring wells (unit:m)

地下水超采是海水入侵最直接的原因，地下水位的持續(xù)下降會破壞水壓力平衡，當沿海含水層地下水位低于海平面時，就會發(fā)生海水入侵，因此在分析區(qū)域內海水入侵狀況時，可以根據(jù)地下水位的動態(tài)變化近似推斷海水入侵狀況。對海水入侵的監(jiān)測，通常是對地下水中氯化物和礦化度進行分析[38]，也有相關研究[39–40]用電導率替代Cl–作為監(jiān)測海水入侵的指標，指出水中Cl–含量與電導率成正比。如圖7 所示，地下水位和電導率在一定程度上呈負相關，即地下水位下降，電導率上升。當?shù)叵滤幌陆禃r，淡水入海流量減小，海水入侵加劇，地下水中電導率升高。因此，從圖中可以看出地下水位動態(tài)變化特征在年內時間尺度可較為準確表征海水入侵的變化。

圖7 監(jiān)測井BH4 水位與電導率變化Fig.7 Groundwater level and electrical conductivity variation in BH4 monitoring well

3.1.2 電導率動態(tài)變化特征

根據(jù)研究區(qū)2011–2012 年各監(jiān)測井采集的電導率動態(tài)數(shù)據(jù)，計算各監(jiān)測井電導率的最大值、最小值、平均值、標準差以及標準離散率等空間變化統(tǒng)計特征，結果見表2。

表2 監(jiān)測井電導率統(tǒng)計特征分析（單位：mS/cm）Table 2 Analysis on statistical characteristics of electrical conductivity in monitoring wells (unit:mS/cm)

標準離散率θ的計算公式為

式中，SD為電導率的標準差；K為電導率的平均值。

溫度為20℃時，天然水的電導率介于0.5～1.5 mS/cm。在濱海地區(qū)，電導率可以作為一項評判海水入侵的指標，電導率的值越小說明地下水受海水入侵程度越低[41–42]。由上表可以得出，監(jiān)測井BH1、BH2 和BH3 的電導率值大致在天然水電導率的范圍附近，沒有出現(xiàn)電導率異?，F(xiàn)象，因此這幾口監(jiān)測井的地下水基本未受海水入侵的影響；而BH4、BH5 及BH6 的電導率值遠超出該范圍，說明地下水在很大程度上受到海水入侵的影響。

從電導率的平均值來看，BH5 的值最高，水質相對較差，其主要原因可能與該處的地下水開采所導致的海水入侵有關，距離它較近的BH4 和BH6 的值也較高；BH3 的平均電導率最低，為1.20 mS/cm，說明此處的地下水質相對較好。監(jiān)測井BH4 的標準離散率最高，為47.91%，說明該點地下水電導率相差較大；BH2 的標準離散率最低，為1.27%，說明該點的電導率的穩(wěn)定性較高。

3.1.3 小波相干性分析

為了理解潮汐對地下水位和地下水電導率的影響，引入小波相干性識別潮汐與地下水位和電導率之間的統(tǒng)計關系。潮汐與地下水位和電導率之間的小波相干性分析結果如圖8 所示（圖中以監(jiān)測井BH2和BH6 為例），譜之間的相位關系用箭頭的方向表示。影響錐（白色虛線）表示小波譜的區(qū)域，在該區(qū)域由于基礎數(shù)據(jù)的有限長度特性而產生的邊緣效應不能被忽略。

圖8 小波相干性分析Fig.8 Wavelet coherence analysis

（1）地下水位和潮汐之間的小波相干性分析

小波相干性分析結果表明，每口監(jiān)測井的地下水位與潮汐數(shù)據(jù)之間存在不同的相干性。監(jiān)測井BH2 和BH6，在整個時間序列中頻譜值有相應變化（圖8a 和圖8b），這些頻率代表地下水位對潮汐變化的響應。對于監(jiān)測井BH2 和BH6，地下水位和潮汐在0.25 d 頻率上的相干性很小，0.5 d 相干性較高，1 d的相干性相對較弱，但相關分布在不同時間是間歇性的。在0.5 d 周期內，與監(jiān)測井BH6 同期地下水位和潮汐資料相比，監(jiān)測井BH2 中的地下水位對潮汐的響應時間不同，相干性更弱，這可能與BH2 的位置及地質條件有關。小波相干性分析表明，海岸帶含水層的地下水位受到海洋潮汐變化的影響。潮汐信號對地下水位的影響具有時空差異?？拷０兜牡貐^(qū)更容易受潮汐影響。同時，某些地區(qū)的潮汐效應一般與它們相對于海洋的位置和地質條件的綜合效應相關[43]。

（2）地下水電導率和潮汐之間的小波相干性分析

圖8c 和圖8d 顯示，將小波相干性分析應用于監(jiān)測的地下水電導率數(shù)據(jù)和潮位時間序列時，未檢測到連續(xù)的短期時間相關性。監(jiān)測井BH2 在2012 年7 月1 日至2012 年10 月1 日時間段內電導率和潮汐之間存在16～32 d 的長周期相關性。監(jiān)測井BH6 中電導率和潮汐之間存在短期的相關性，可以看出不同井位的電導率對潮汐變化的響應不同。潮汐對地下水電導率的影響要弱于對地下水位的影響，因為地下水位的波動是受應力變化的影響，而電導率的變化需要含水層儲水量的體積變化。更重要的是，地質條件對控制潮汐對地下水位的影響至關重要。

3.2 模擬結果分析

3.2.1 水位模擬結果分析

WA-NARX 神經(jīng)網(wǎng)絡混合模型在NARX 神經(jīng)網(wǎng)絡模型的基礎上，引入小波分析，對數(shù)據(jù)進行小波分解以減小預測誤差，提高預測精度。在進行海岸帶地下水位預測時，把降雨和潮汐作為輸入?yún)?shù)，經(jīng)WANARX 混合模型模擬得到地下水位預測結果，同水位的觀測值進行對比，得到預測結果的散點圖和時間序列圖。圖9 顯示了BH1–BH6 這6 口監(jiān)測井地下水位的預測結果。

圖9a、圖9c、圖9e、圖9g、圖9i、圖9k 為預測結果散點圖，圖9b、圖9d、圖9f、圖9h、圖9j、圖9l 為預測結果時間序列圖。通過預測結果散點圖可得：訓練數(shù)據(jù)和預測數(shù)據(jù)集中分布一條直線上，模型在水位變幅較大時預測效果有所降低。如圖所示，監(jiān)測井BH3 的散點出現(xiàn)小幅度的偏離。通過預測結果的時間序列圖可得：訓練和測試過程水位變化趨勢與觀測水位變化趨勢一致。混合模型在訓練和測試過程中均表現(xiàn)出較低的誤差，在地下水位波動處誤差相對較高，但模型總體上表現(xiàn)出良好的模擬效果。WANARX 神經(jīng)網(wǎng)絡混合模型在訓練階段和測試階段的RMSE 和R2性能指標值見表3，從表3 可見，所有監(jiān)測井在訓練和測試階段的R2均大于0.99，RMSE 均小于0.03，表明模型性能良好。

圖9 WA-NARX 模型地下水位預測與觀測對比Fig.9 Comparison of groundwater level prediction based on WA-NARX and observation

表3 地下水位預測時WA-NARX 模型的均方根誤差（RMSE）和擬合度（R2）統(tǒng)計表Table 3 The WA-NARX model performance metrics of RMSE and R2 for groundwater level prediction

3.2.2 電導率模擬結果分析

BH1–BH6 6 口監(jiān)測井電導率的預測結果散點圖和時間序列圖如圖10 所示，從混合模型模擬結果可以看出：混合模型的訓練數(shù)據(jù)和測試數(shù)據(jù)基本集中分布在一條直線上，電導率的變化趨勢與觀測電導率的變化趨勢一致，在訓練和測試階段模型的誤差都相對較低。因此WA-NARX 神經(jīng)網(wǎng)絡混合模型適用于地下水電導率的模擬預測。表4 中所有監(jiān)測井在訓練和測試階段的R2均大于0.98，RMSE 均小于0.009，說明該模型在進行電導率的模擬預測時也表現(xiàn)出良好的性能。

表4 電導率預測時WA-NARX 模型的均方根誤差（RMSE）和擬合度（R2）統(tǒng)計Table 4 The WA-NARX model performance metrics of RMSE and R2 for electrical conductivity prediction

圖10 WA-NARX 模型電導率預測的性能圖Fig.10 Performance of WA-NARX model for electrical conductivity prediction

3.3 蒸發(fā)和排灌用電量對混合模型預測的影響

3.3.1 蒸發(fā)對水位預測的影響

在對比不同輸入?yún)?shù)對地下水位預測結果的影響時，在模型訓練完成后，關閉WA-NARX 網(wǎng)絡的回路，輸出反饋進行多步預測。以60 h 為時間步長進行多步預測。表5 為僅有降雨、潮汐參數(shù)的地下水位預測結果與新增蒸發(fā)參數(shù)后（降雨、潮汐和蒸發(fā)）的地下水位預測結果的對比。由預測誤差可知，輸入?yún)?shù)僅為降雨和潮汐時，模型的預測精度較高；當在此基礎上新增蒸發(fā)作為輸入?yún)?shù)時，模型預測結果發(fā)生了變化。與無新增參數(shù)的結果相比，監(jiān)測井BH1、BH3 和BH6 的預測精度有所提高，其余監(jiān)測井的預測精度有所降低。經(jīng)過整體對比發(fā)現(xiàn)，只有降雨和潮汐作為模型的輸入?yún)?shù)時，模型的預測效果略高于新增蒸發(fā)參數(shù)后的預測效果。說明在該地區(qū)進行地下水位預測時，降雨和潮汐對地下水位變化的影響更為顯著。

圖11 是新增蒸發(fā)參數(shù)后得到的預測水位與實測水位的散點圖。圖中實線是觀測水位數(shù)據(jù)的連線，把降雨和潮汐作為WA-NARX 模型的輸入?yún)?shù)，預測60 h 時間步長的地下水位時，模型的預測值基本分布在實測的水位線上，預測效果良好。新增蒸發(fā)參數(shù)后，監(jiān)測井BH1、BH3 和BH6 的預測散點值更集中于實測水位線，表現(xiàn)出更好的預測效果，而其余觀測井的預測效果相比之前有所下降，這與表5 的分析結果一致。

圖11 新增蒸發(fā)參數(shù)的預測和實測地下水位散點圖Fig.11 Scatter diagrams of predicted and observed groundwater level with newly added evaporation parameters

表5 新增蒸發(fā)參數(shù)的地下水位預測誤差及網(wǎng)絡性能統(tǒng)計表Table 5 The prediction error of groundwater level and network performance with newly added evaporation parameters

3.3.2 排灌用電量對水位預測的影響

農業(yè)排灌用電量與農業(yè)抽水灌溉有關，灌溉則與地下水位密切相關，所以把排灌用電量作為一組參數(shù)，研究其對地下水位預測的影響。表6 為新增排灌用電量參數(shù)與只有降雨和潮汐作為參數(shù)時地下水位預測結果的對比。由預測誤差可知，只有降雨和潮汐作為輸入?yún)?shù)時模型有較高的預測精度；當新增排灌用電量作為輸入?yún)?shù)后，模型的預測精度同之前相差甚小。其中監(jiān)測井BH5、BH6 的預測效果有細微的提高，其余監(jiān)測井的預測效果有所降低。說明在研究區(qū)內農業(yè)排灌用電量數(shù)據(jù)可以作為預測地下水位的一項參數(shù)，在使用時要結合研究區(qū)監(jiān)測井的實際情況，在農業(yè)抽水灌溉量較大處的影響會更顯著。同時從電量時間序列可以看出，因統(tǒng)計關系，每個月的灌溉用電量基本為定值，而地下水位的預測周期是小時變化，所以，用電量的觀測與地下水預測頻率不完全一致，對模擬結果也有較大影響。如果用電量的統(tǒng)計頻率與預測相符，模型的預測精度應有所提高。

表6 新增用電量參數(shù)的地下水位預測誤差及網(wǎng)絡性能統(tǒng)計Table 6 The prediction error of groundwater level and network performance with newly added power consumption parameters

圖12 是新增排灌用電量參數(shù)得到的預測水位與實測水位的散點圖。從圖中可以看出，每口井對新增排灌用電量參數(shù)后的預測效果響應并不盡相同。BH1 與BH6 的預測結果相對近似；而BH2、BH3 及BH4 的預測水位分布及表6 中RMSE 共同顯示其預測值較只應用降水與潮汐時有所下降。對于BH5 來說，加入用電量后預測結果整體向上偏移，精度有一定提高。對比結果說明井位地理位置分布及其區(qū)域地質條件對用電量的響應有較大關聯(lián)。

圖12 新增用電量參數(shù)的預測和實測地下水位散點圖Fig.12 Scatter diagrams of predicted and observed groundwater level with newly added power consumption parameters

3.3.3 蒸發(fā)對電導率預測的影響

圖13 是新增蒸發(fā)參數(shù)后得到的預測電導率與觀測電導率的曲線圖，表7 為不同輸入?yún)?shù)地下水電導率預測結果的對比統(tǒng)計表。由分析結果可知，在降雨、潮汐基礎上增加蒸發(fā)數(shù)據(jù)作為輸入?yún)?shù)時，無論在變化平緩區(qū)還是劇烈區(qū)，監(jiān)測井中電導率的預測精度都有顯著提高。特別對監(jiān)測井BH5，未加入蒸發(fā)數(shù)據(jù)時預測值略有偏高，加入蒸發(fā)對該井有明顯的校正作用。預測結果說明在該地區(qū)預測地下水電導率時，蒸發(fā)是重要影響因素。

圖13 新增蒸發(fā)參數(shù)的預測和實測電導率曲線圖Fig.13 Curve diagrams of predicted and observed electrical conductivity with newly added evaporation parameters

3.3.4 排灌用電量對電導率的影響

為了對比不同輸入?yún)?shù)對地下水電導率預測結果的影響，關閉WA-NARX 網(wǎng)絡的回路，輸出反饋進行多步預測。把收集到的農業(yè)抽水排灌用電量數(shù)據(jù)作為一項參數(shù)，對比其對地下水電導率預測效果的影響。由表7 可知，只有降雨和潮汐作為輸入?yún)?shù)時，模型有較好預測效果；新增農業(yè)排灌用電量參數(shù)，監(jiān)測井BH1、BH2 和BH5 的預測效果有所提高，其余監(jiān)測井的預測效果基本保持一致，總體上新增排灌用電量參數(shù)后混合模型的預測效果有所提升。

表7 不同輸入?yún)?shù)電導率預測誤差及網(wǎng)絡性能評價統(tǒng)計Table 7 The prediction error of electrical conductivity and network performance with different input parameters

圖14 是新增排灌用電量參數(shù)得到的預測電導率與觀測電導率的曲線圖。預測值的變化趨勢和實測電導率的變化趨勢基本吻合，表現(xiàn)出較好的預測效果；當新增農業(yè)排灌用電量參數(shù)時，監(jiān)測井BH1、BH2 和BH5 的預測曲線與觀測電導率曲線的變化趨勢更加貼合，與表7 的分析結果一致。

圖14 新增用電量參數(shù)的預測和實測電導率曲線圖Fig.14 Curve diagrams of predicted and observed electrical conductivity with newly added power consumption parameters

4 結論與展望

萊州灣南岸是受海水入侵影響嚴重地區(qū)。通過對沿海典型剖面含水層地下水位和電導率的變化分析可知，該區(qū)海水入侵季節(jié)性變化趨勢明顯。主要表現(xiàn)在地下水電導率與水位存在一定的負相關關系，在降雨補給充足的情況下，地下水水位有較為明顯的上升趨勢，地下水電導率降低，海水入侵程度減弱。但由于受地下水開采的影響，水位和電導率均存在突變現(xiàn)象，水位降低打破原有的咸淡水平衡，使海水入侵加劇。同時，本文建立了基于小波分析的非線性自回歸神經(jīng)網(wǎng)絡混合模型（WA-NARX），分別對地下水位和地下水中的電導率進行預測，通過增加模型的輸入?yún)?shù)，對比分析不同輸入?yún)?shù)對模型預測效果的影響。研究發(fā)現(xiàn)，地下水位和電導率的變化與降雨和潮汐的變化密切相關。當降雨和潮汐作為模型的輸入?yún)?shù)預測地下水位與電導率時，模型預測效果良好。蒸發(fā)參數(shù)和農業(yè)排灌用電量對地下水位與電導率也有一定影響。模型預測依然具有較高精度，個別觀測點精度有所提高。通過上述研究表明，WA-NARX 模型可以作為海水入侵預測的有效工具，不僅可以為水資源的高效利用提供參考，還可以掌握海岸帶地下水鹽度的變化，為海水入侵的監(jiān)測預測預警、保護地下淡水資源及海岸帶環(huán)境保護提供科學依據(jù)。