基于蜻蜓-高斯過程回歸耦合的居民社區(qū)時用水量動態(tài)實時區(qū)間預測

劉龍龍，李文竹，劉 心

(1.河北工程大學，河北 邯鄲 056038；2.邯鄲河務局，河北 邯鄲 056001)

我國水資源存在的時空分布不均、用水矛盾突出等現實問題，越來越成為新時期經濟社會發(fā)展的制約性因素[1]。水資源日益短缺，對區(qū)域水資源進行合理規(guī)劃、利用和調度勢在必行。居民社區(qū)用水是水資源消耗的重要組成部分[2]。科學的管網調度能節(jié)省大量的供水消耗，全面提高管網的安全性和可靠性，但調度方案是否可行很大程度上取決于用水量的預測精度[3]。

對用水量預測方法的研究一直是國內外學者研究的熱點，主要有以時間序列法為代表的傳統方法、以人工神經網絡(ANN)為代表的人工智能方法和以支持向量機(SVM)為代表的機器學習算法[4]。然而這些方法都只能得到確定的點預測結果和未來某一時刻的預測結果，而居民實際用水量蘊含各種不確定性因素，點預測結果必然存在不同程度的誤差，進而使調度決策工作面臨一定程度的風險。另外，如果預測的結果不是實時的，僅是未來某一時刻的值，那么調度人員無法根據準確的時用水量值進行實時調度。如果能夠給出實時預測結果和其變化區(qū)間，使調度決策人員隨時了解未來用水量的變化波動范圍，有利于其做出更合理的決策。

高斯過程回歸(GPR)是通過尋找訓練數據之間的關系來進行系統辨識的一種非參數黑箱模型,是用于概率問題預測建模的一種強大工具[5]。與ANN和SVM方法相比，GPR不僅能夠預測未知量的期望值，還能給出其分布狀況,同時也能進行實時預測。此外，GPR模型參數較少、易實現和泛化能力強等優(yōu)點[6]，可以直接方便的用于區(qū)間預測和實時預測[7]。通常采用共軛梯度法求解最優(yōu)超參數，但存在優(yōu)化效果對初始值依賴性強和容易陷入局部最優(yōu)的缺點。

蜻蜓算法是美國學者MIRJALILI S在2015年提出的一種基于種群的啟發(fā)式智能優(yōu)化算法[8]。自然界中的蜻蜓生活方式包括靜態(tài)和動態(tài)群體行為，這和啟發(fā)式優(yōu)化算法中的兩個重要階段探索和開發(fā)非常類似[9]。蜻蜓算法通過模擬蜻蜓群體航行、捕食及躲避外敵等行為進行全局和局部搜索，尋找獵物的過程就是算法尋優(yōu)的過程[10,11]。本文提出用擅長全局搜索和對初始值沒有依賴的蜻蜓算法進行超參數尋優(yōu)，形成蜻蜓-高斯過程回歸耦合算法(DA-GPR)，對居民社區(qū)時用水量進行動態(tài)實時區(qū)間預測。

1 蜻蜓-高斯過程回歸耦合預測模型

本文通過高斯過程回歸與蜻蜓算法的有機結合，利用擅長全局搜索和對初始值沒有依賴的蜻蜓算法進行最優(yōu)超參數尋優(yōu)，并利用高斯過程回歸(GPR)尋找時用水量訓練數據之間的關系來進行系統辨識,給出其分布狀況,對居民社區(qū)時用水量進行動態(tài)實時區(qū)間預測。

1.1 DA-GPR基本原理

從函數空間角度來看，高斯過程回歸其性質由均值函數m(x)和協方差函數k(x,x′)確定[12]，定義為：

(1)

式中：函數f(x)的高斯過程(Gaussian process,GP)數學表達式為f(x)～GP[m(x),k(x,x′)],x,x′∈Rd是任意的隨機變量。

(2)

式中：K(X,X)=Kn=(kij)為n×n階對稱正定協方差矩陣。

預測值的后驗分布為：

(3)

(4)

(5)

本文的協方差函數采用平方指數協方差函數，為：

(6)

1.2 DA-GPR目標函數設計

蜻蜓-高斯過程回歸耦合時用水量預測模型的最終目標是對未來的時用水量進行較準確的動態(tài)實時區(qū)間預測，并具有較好的泛化能力。因此在使用蜻蜓算法優(yōu)化高斯過程回歸的參數時，選取能充分代表泛化能力的目標函數是非常重要的。常規(guī)的蜻蜓優(yōu)化算法將訓練集上獲得最小平均相對誤差(MAPE)的參數作為最優(yōu)解。即，將訓練集上獲得最小平均相對誤差(MAPE)作為目標函數，對蜻蜓個體進行選著，目標函數值越小，蜻蜓個體越優(yōu)，對應的參數越優(yōu)，目標函數的選著至關重要。

居民社區(qū)每小時的用水量比較少，基本上在幾立方米到十幾立方米，僅僅采用平均相對誤差(MAPE)作為目標函數值進行參數尋優(yōu)的效果不理想。在此基礎上，本文對DA-GPR進行改進，將體現絕對誤差大小的平均絕對誤差(MAE)和均方誤差(MSE)引入到目標函數中，將平均相對誤差(MAPE)、平均絕對誤差(MAE)和均方誤差(MSE)進行線性組合，構造出新的目標函數，如下式所示。

F(x)=3MAPE+2MAE+MSE

(7)

1.3 DA-GPR局部搜索公式設計

蜻蜓群體在一個小的區(qū)域內來回飛翔尋找獵物以及躲避外敵。通過模擬蜻蜓群體之間尋找獵物以及躲避外敵的社會互動行為來實現優(yōu)化的過程，就是算法的局部搜索過程。

蜻蜓個體在群體運動中可以分為分離、對齊、內聚、食物吸引和天敵排斥等5種行為模式。該行為的具體意義和數學表達方法如下[13]。

(1)分離：指蜻蜓與相鄰個體之間避免碰撞。該行為的計算公式如下:

(8)

(2)對齊：指相鄰個體之間傾向于保持相同的速度。該行為的計算公式如下:

(9)

(3)內聚：指蜻蜓個體試圖向自己認為所屬的群體靠近。該行為的計算公式如下:

(10)

(4)食物吸引力：指食物對蜻蜓的吸引力。該行為的計算公式如下:

Fi=X+-X

(11)

(5)天敵排斥力：指蜻蜓對天敵的排斥力。該行為的計算公式如下:

Ei=X-+X

(12)

式中：Si為第i個個體的分離量；X為當前個體的位置；Xj為相鄰個體j的位置；N為相鄰個體的數量；Ai為第i個個體的對齊量；Vj為第j個鄰近個體的速度；Ci為第i個個體的內聚量；Fi為第i個個體對獵物的吸引力；X+為獵物所在的位置；Ei為第i個個體需逃離天敵的距離；X-為天敵所在的位置。

1.4 DA-GPR全局搜索公式設計

在蜻蜓尋優(yōu)的過程中，大量的蜻蜓群體朝著同一個方向進行長距離遷移，這個群體活動就是對算法的全局搜索過程。

根據5種蜻蜓行為，下一代蜻蜓的步長計算如下：

ΔXt+1=(sSi+aAi+cCi+fFi+eEi)+ωΔXt

(13)

有鄰近蜻蜓時，下一代蜻蜓的位置為：

Xt+1=Xt+ΔXt+1

(14)

無鄰近蜻蜓時，設置為隨機游走(Le′xy飛行)行為，下一代蜻蜓的個體位置為:

Xt+1=Xt+Le′xy(d)×Xt

(15)

式中：t為當前迭代次數；i為第i個蜻蜓個體；Xi為當前第t代種群個體位置；ΔXt+1為下一代種群位置更新步長；Xt+1為下一代種群個體位置；s，a，c，f，e分別對應于5種行為的權重；ω為慣性權重；d為個體位置向量的維數。

Le′vy函數計算公式如下：

(16)

(17)

式中：r1，r2為[0,1]隨機數；Γ(x)=(x-1)!；β是一個常數(這里取為0.5)。

蜻蜓算法在尋優(yōu)的過程中，對每個個體的鄰近個體數量的計算是非常重要的，因此這里假定一個鄰域半徑，該半徑隨迭代次數的增加而成比例地增長。同時，為了達到算法的局部搜索和全局搜索的平衡，權重(s，a，c，f，e和ω)在優(yōu)化過程中自適應地調整。

1.5 DA-GPR算法步驟

流程圖如圖1所示。

圖1 DA-GPR模型區(qū)間預測流程圖Fig.1 DA-GPR model interval prediction flow chart

具體實現步驟如下：

Step1：讀取樣本數據，產生訓練集和測試集。

Step2：初始參數設置：設置蜻蜓種群規(guī)模、最大迭代次數和參數l,σf,σn的取值范圍。

Step3：初始化蜻蜓：隨機初始化蜻蜓的位置H、步長ΔH、鄰域半徑R，慣性權重ω，分離權重s，對齊權重a，內聚權重c，食物吸引權重f，天敵排斥權重e。

Step4：將蜻蜓個體位置H的信息依次賦值給l,σf,σn。其中位置H矩陣的第一行存放參數l的值，第二行存放參數σf的值，第三行存放參數σn的值，每個蜻蜓個體對應一組參數值。

Step5：創(chuàng)建GPR模型，對訓練樣本進行訓練，求出每組參數對應的目標函數值，并將其作為蜻蜓算法的適應度函數值，判斷當前的適應度函數值是否為最佳適應度值，若是則將對應的超參數保存為最優(yōu)超參數值,否則仍保存原適應度值及其所對應的參數值。

Step6：判斷是否達到最大迭代次數，如果達到最大迭代次數，則輸出最優(yōu)超參數并創(chuàng)建最優(yōu)的GPR模型；如果達不到最大迭代次數，則依次更新蜻蜓最優(yōu)個體和最差個體、更新鄰域半徑、更新個體位置，然后返回步驟Step4繼續(xù)迭代。

Step7：將測試集數據輸入創(chuàng)建的最優(yōu)GPR模型，輸出預測值的均值和方差。

Step8：根據給定置信度，得到時用水量實時區(qū)間預測結果。

2 仿真實驗

2.1 仿真實驗平臺

本次實驗的運行環(huán)境為：操作系統：Windows 7；處理器：Intel(R) Core(TM)i7-2640M；主頻：2.80GHz；內存：4.00GB；編程工具：MATLAB R2015(b)。

2.2 樣本數據

本文所用的數據來自河北工程大學在線水量檢測平臺，采用居民家屬院2016年6月17日到7月7日，共21 d每天24 h的實際社區(qū)時用水量資料，如圖2所示。將樣本集21 d/h用水量中的前20 d 480組時用水量資料作為訓練集，最后一天24 h用水量資料作為測試集。

圖2 單日時用水量變化趨勢Fig.2 Change in water consumption on a single day

2.3 模型輸入和輸出

從圖3可以明顯看出居民社區(qū)時用水量具有明顯的24 h周期性，為預測下一時刻用水量，模型的輸入取該時刻前24 h的時用水量。模型輸出是居民社區(qū)下一時刻用水量，輸入是該時刻前24 h每小時用水量(數據來源于2016年6月17日到7月6日居民社區(qū)的用水數據)。

圖3 不同模型時用水量預測結果Fig.3 Water consumption prediction results of different models

2.4 參數設置

慣性權重ω、分離權重s、對齊權重a、內聚權重c、食物吸引權重f和天敵排斥權重e均為自適應線性遞減權重，其最小值為0.4，最大值為0.9；設種群規(guī)模為10，最大迭代次數為20，參數M的取值范圍為[-1，1]，參數σf的取值范圍為[-1，1]，參數σn的取值范圍為[-1，1]。

2.5 對比模型及評價標準

為了使DA-GPR模型的預測結果具有可比性，選著BP神經網絡、最小二乘支持向量機和高斯過程回歸模型作為對比模型，模型性能采用平均相對誤差(MAPE)、平均絕對誤差(MAE)和均方誤差(MSE)進行衡量，它們的定義為：

(18)

(19)

(20)

2.6 結果與分析

為了驗證本文創(chuàng)建的蜻蜓-高斯過程回歸耦合模型的預測可靠性，本文將其與BP神經網絡、最小二乘支持向量機和高斯過程回歸預測模型進行試驗對比，預測結果如圖3所示。在95%置信度下，高斯過程回歸和蜻蜓-高斯過程回歸兩個模型的區(qū)間預測結果如圖4和圖5所示。從圖3可知，BP神經網絡和最小二乘支持向量機模型預測結果誤差較大，蜻蜓-高斯過程回歸模型預測的結果誤差較小。從圖4可知，高斯過程回歸模型預測的結果誤差較大，極個別點落在了預測區(qū)間外。從圖5可知，蜻蜓-高斯過程回歸預測模型預測的結果走勢與實測值一致，預測值與實測值基本上重合，誤差較小。

圖4 高斯過程回歸(GPR)時用水量預測結果Fig.4 Water consumption prediction results in Gaussian process regression (GPR)

圖5 蜻蜓-高斯過程回歸(DA-GPR)時用水量預測結果Fig.5 Prediction of water consumption in the dragonfly Gauss process regression (DA-GPR)

不同預測模型的預測相對誤差如圖6所示。從圖6可以看出，BP神經網絡預測模型相對誤差最大為0.5，最小二乘支持向量機預測模型相對誤差最大為0.9，高斯過程回歸預測模型相對誤差最大為1.653，蜻蜓-高斯過程回歸預測模型相對誤差最大為0.019。顯然，蜻蜓-高斯過程回歸預測模型具有更高的預測精度。

圖6 不同預測模型預測的相對誤差Fig.6 Relative error of prediction model of different prediction models

采用平均相對誤差(MAPE)、平均絕對誤差(MAE)和均方誤差(MSE)對不同模型的預測結果進行分析，結果如表1所示。通過表1可以看出，蜻蜓-高斯過程回歸模型的居民社區(qū)時用水量預測精度優(yōu)于BP神經網絡、最小二乘支持向量機和高斯過程回歸模型，取得了較好的預測效果。

表1 不同預測模型誤差分析Tab.1 Error analysis of different prediction models

3 結 語

居民社區(qū)時用水量具有較強的不確定性和隨機性，傳統的ANN、SVM等預測方法只能得到確定的點預測結果和未來某一時刻的預測結果，無法給出預測的區(qū)間，也不能進行實時預測。為了克服這些缺點，本文提出了一種基于蜻蜓-高斯過程回歸耦合的居民社區(qū)時用水量動態(tài)實時區(qū)間預測方法。為進一步提高預測精度，進行了改進，最終得到了一定置信水平下的區(qū)間預測結果。仿真結果表明，本文構建的區(qū)間預測方法與常規(guī)方法相比，不僅能夠預測未知量的期望值，還能給出其分布狀況,同時也能進行實時預測。而且預測精度較高，最大的相對誤差為僅0.019，具有較強的實用價值，為未來水資源實時調度提供理論依據。