李　匡，胡宇豐，梁犁麗，徐海卿

PSO算法在馬斯京根法參數率定中的應用

李匡，胡宇豐，梁犁麗，徐海卿

（北京中水科水電科技開發(fā)有限公司，北京 100038）

將PSO算法用于馬斯京根法參數的自動率定，開發(fā)了基于PSO算法的馬斯京根法參數自動率定系統(tǒng)，該系統(tǒng)可以對馬斯京根法參數進行快速的率定和檢驗；是分析率定結果穩(wěn)定性的一種新方法；此系統(tǒng)具有較高的實用性，已在多個流域中得到了應用。

馬斯京根法；參數率定；PSO算法；參數自動率定系統(tǒng)；數據穩(wěn)定性

1　引言

馬斯京根法是一種廣泛應用于河道洪水流量演算的方法，是美國工程師G.T.Mccarthy于1938年提出的流量演算法，此法最早在美國馬斯京根河流域上使用，因而稱為馬斯京根法。該法建立馬斯京根槽蓄曲線，并與水量平衡方程聯解，進行河段洪水演算[1]。該方法在實際應用中的一個重要問題是模型的參數估計。在傳統(tǒng)方法中利用試算法對參數K、x進行率定，再計算出決策變量C0、C1、C2，但是試算法具有盲目性、不確定性、不易程序化且精度不高等缺陷[1-2]。因此有學者選擇了直接優(yōu)化決策變量C0、C1、C2，得到較好的效果，主要方法有最小二乘法[3]，遺傳算法[4-5]，粒子群算法[6]等。

粒子群優(yōu)化算法（ParticleSwarmOptimization Algorithm，PSO）是進化算法的一種，系統(tǒng)初始化為一組隨機解，通過迭代搜尋最優(yōu)值。傳統(tǒng)PSO算法在應用中和其它全局優(yōu)化算法一樣可能會陷入局部最優(yōu)，導致收斂精度不高，后期收斂速度較慢，經過學者的研究對粒子群算法進行了改進，使得粒子群算法的效率得到了較大的提高[7-8]。另外，同其它的進化算法，如遺傳算法相比，PSO模型具有參數較少，原理簡單，易于實現，計算效率高等優(yōu)點。鑒于以上優(yōu)點，本文選擇PSO算法率定馬斯京根法參數。

馬斯京根法參數率定是一個非線性優(yōu)化問題，利用PSO算法進行率定是合適的。在應用PSO算法進行馬斯京根法參數率定時，需要給出馬斯京根法各參數的取值范圍，然而在已有的研究文獻中，作者并沒有明確給出各參數的取值范圍。本文根據嚴格的數學證明給出了馬斯京根法參數的取值范圍，為廣大學者在應用各類優(yōu)化算法進行馬斯京根法參數率定提供參考。

作者開發(fā)了基于PSO算法的馬斯京根法參數自動率定系統(tǒng)，按照給定的數據庫表結構輸入資料（流量，斷面信息），系統(tǒng)可以對流域中各河道的馬斯京根法參數進行自動率定，并對率定出的馬斯京根法參數進行檢驗。

2　馬斯京根法簡介

馬斯京根法基本原理如下：對水量平衡方程（1）和馬斯京根槽蓄方程（2）進行求解：

式中：I上斷面的流量；I1，I2時段始末上斷面的流量；Q下斷面的流量；Q1，Q2時段始末下斷面的流量；△t計算時段長；S河段蓄水量；S1，S2時段始末河段蓄水量；K穩(wěn)定流情況下的河段洪水傳播時間；x槽蓄系數。C0，C1，C2是 K，△t，x的函數，并且C0+C1+C2=1

3　粒子群（PSO）算法及其計算流程

PSO算法是Kennedy和Eberhart于1995年提出的一種基于對鳥群捕食行為模擬的智能群集優(yōu)化算法。其基本思想可以簡單地表述為：在一個優(yōu)化問題的解空間中，每一個可行解被看做一個“粒子”，這些粒子在解空間內不停地飛行，在飛行的過程中根據自身和種群中其他粒子積累的經驗不斷調整自己的飛行策略，最終這些粒子都趨于解空間中的最優(yōu)區(qū)域，也即所謂的“食物”[10]。

PSO算法的數學描述如下：在一個n維的搜索空間中，由m個粒子組成一個種群，即X={x1，x2，x3，…，xi，…xm}T，第i個粒子的位置為xi=（xi1，xi2，xi3，…，xin），其速度為Vi=（Vi1，Vi2，Vi3，…，Vin），粒子的個體極值為：Pi=（Pi1，Pi2，Pi3，…，Pin），種群的全局極值為：Pg=（Pg1，Pg2，Pg3，…，Pgn），基本PSO算法的粒子在搜索過程中通過式（3）和式（4）不斷地進行位置和飛行速度的更新：

其中：d=1，2，3，…，n；i=1，2，3，…，m。

式中：n為粒子的維數；m為種群的規(guī)模；t為當前的迭代步數；r1，r2分別為0和1之間的隨機數，c1，c2為加速度常數。

粒子的維數n與具體的研究問題有關，例如在率定馬斯京根法參數時，有三個參數需要率定，則n=3。

種群規(guī)模m影響算法的穩(wěn)定性和效率，m取值越大，PSO算法的穩(wěn)定性越好，但同時計算量也較大，耗時長，計算效率較低。因此在對精度要求較高時可選擇較大的種群規(guī)模，對計算效率要求較高時選擇較小的種群規(guī)模，種群規(guī)模m的取值范圍一般在[50，100]之間。

從式（3）可以看出，粒子飛行速度更新公式包括三部分：第一部分是粒子的歷史飛行速度，說明了粒子目前的狀態(tài)，起到平衡全局和局部搜索能力的作用；第二部分是粒子自身的認識，表示粒子在飛行過程中自身的思考；第三部分是粒子的社會認識，表示粒子群中各粒子之間信息的交流。三個部分共同作用，決定了粒子的空間搜索能力。

傳統(tǒng)的粒子群算法容易陷入局部極小值，影響PSO算法找到全局最優(yōu)值。為了提高PSO算法的優(yōu)化性能，平衡算法的全局搜索能力和局部搜索能力，文獻[11]對基本PSO算法進行了改進，在基本PSO算法的基礎上引進了慣性權重項ω，提出了標準PSO算法。式（3）也相應地修正為：

式中ω的取值方式有兩種。一種是遞減策略，慣性權重ω滿足ω（t）=0.9-t×0.5，其中t為當前迭T 代次數，T為最大迭代次數；另外一種是固定策略，通常取慣性權重為[0.4，0.6]中的某一固定值。

根據PSO算法基本思想，可以設計標準PSO算法的流程如下：

（1）將優(yōu)化的問題數學模型化，選定優(yōu)化問題的目標函數（即適應度函數），選定粒子的維數；

（2）初始化算法。對粒子群中的粒子位置和速度進行初始化設定，即在一定的范圍內隨機產生出每一個粒子的位置和速度；

（3）根據優(yōu)化問題的目標函數計算每個粒子的適應度值；

（4）對每個粒子，將其當前適應度值與其所經歷的最優(yōu)適應度值進行比較，如果該粒子的當前適應度值更優(yōu)，那么將當前位置記錄為該粒子的局部最優(yōu)位置；

（5）對每個粒子所經歷的最優(yōu)適應度值與全局最優(yōu)適應度值進行比較，如果個體粒子的最優(yōu)適應度值較全局最優(yōu)適應度值為優(yōu)，則將其作為全局最優(yōu)位置；

（6）對每個粒子的位置和飛行速度進行更新；

（7）判斷是否達到優(yōu)化的終止條件。如果滿足終止條件，結束循環(huán)，否則返回第三步。

4　馬斯京根法參數率定系統(tǒng)

本文利用Vb.net開發(fā)了馬斯京根法參數率定系統(tǒng)，該系統(tǒng)利用PSO算法對馬斯京根法參數進行率定。將各水文站的測站信息、流量資料等按照系統(tǒng)設計的數據庫表結構輸入到數據庫中，系統(tǒng)即可按照默認設置的PSO算法參數對馬斯京根法參數進行率定，并且可以對率定的參數進行檢驗及分析。

在用PSO算法進行馬斯京根法參數率定時，還需要解決以下幾個問題：

4.1目標函數

采用PSO算法率定模型參數時必須構建目標函數，一方面用以評價模擬效果的優(yōu)劣，另一方面作為粒子飛行策略調整的標準。在對馬斯京根法參數進行率定時，一般有河槽蓄量誤差最小和流量誤差最小兩種目標函數，經學者研究表明，選用流量誤差最小作為目標函數率定的馬斯京根法參數效果更好[12]。

本系統(tǒng)將計算流量和實測流量的誤差平方和最小作為判斷馬斯京根法參數尋優(yōu)的目標函數：

式中：Qi為第i時段計算出流量，QiC為第i時段實測出流量，n為總時段數。

4.2參數率定結果的穩(wěn)定性分析

在利用PSO算法率定馬斯京根法參數時，由于PSO算法的初始解是隨機給定的，因此可能會造成率定的參數不穩(wěn)定。該問題可以通過多次計算，得到率定參數結果序列，然后檢驗參數結果序列是否穩(wěn)定，如果穩(wěn)定，則取參數結果的平均值作為最終率定參數結果；如果未通過穩(wěn)定性檢驗，則調整PSO參數重新率定，直到率定的參數結果穩(wěn)定為止。

可以通過以下算法檢驗率定的參數結果是否穩(wěn)定：

該算法的核心是計算聚集在平均值附近的數值的概率，當該概率大于穩(wěn)定指標時，則確定該數據序列是穩(wěn)定的，反之數據不穩(wěn)定。判斷過程可按照以下步驟進行：

（1）分析數據的特點，確定合適的穩(wěn)定性指標S以及穩(wěn)定性判別變化幅度Adec；

（2）計算該段數據序列的平均值Vavg；

（3）給定變化幅度A（0≤A≤1）圍繞平均值Vavg上、下波動，得到區(qū)間[Vavg-Vavg×A，Vavg+Vavg×A]，稱此區(qū)間為穩(wěn)定區(qū)域F。計算穩(wěn)定區(qū)域內數據出現的概率P，即穩(wěn)定區(qū)域內出現的數據個數與整個時間段內數據的總個數的比值。如果P大于穩(wěn)定指標F并且abs（P-S）≥ε，則轉至（4）繼續(xù)；如果P小于穩(wěn)定指標S并且abs（P-S）≥ε，則轉至（5）繼續(xù)；如果abs（P-S）≥ε，則轉至（6）；

（4）減小變化幅度A，轉至（3）重新計算；

（5）增大變化幅度A，轉至（3）重新計算；

（6）如果A≤Adec，則認為數據序列穩(wěn)定，否則認為數據序列不穩(wěn)定。

此法的優(yōu)點是：不僅可以判斷數據是否穩(wěn)定，而且可以判別出數據在多大的范圍內穩(wěn)定。

4.3有支流斷面的參數率定問題

對有支流匯入的斷面，系統(tǒng)中采用“先合后演”的方法進行參數率定，即將干、支流所有上斷面相同時刻的流量相加，作為河段總入流，按無支流河段方法進行演算[2]。

4.4系統(tǒng)中PSO算法默認參數

（1）計算次數：ComputeCount=50；

（2）種群規(guī)模：m=70；

（3）粒子維數：n=70；

（5）限制條件C2=1-C0-C1，當C2超出其取值范圍，則調整至其取值范圍；

（6）加速因子：設定為C1=C2=2.0；

（7）慣性權重：從0.9～0.4遞減變化或者取固定值ω=0.5，默認選擇后一種；

（8）速度限制系數：k=0.7；

（9）終止條件：最大迭代次數T=200，或者目標函數之差小于給定精度ε=10-8。

表1　率定誤差與檢驗誤差比較

在實際應用中用戶可以通過對以上參數進行修改，找出率定效果較好的馬斯京根法參數。

5　實例

采用雅礱江流域上三個測站的資料進行馬斯京根法參數率定和檢驗。三個測站的地理位置示意圖如圖1所示：

圖1　測站地理位置示意圖

將麥地龍站和列瓦站同時刻的資料相加作為入流量，將錦屏站的資料作為出流量。利用2006～2008年的資料進行參數率定，再利用2009年的資料進行檢驗，我們得到如下結果：

率定結果C0=0.056，C1=0.001，C2=0.942。

資料的率定誤差和檢驗誤差分別如表1。

先對每次計算所得的最優(yōu)目標函數值序列以及參數序列分別進行歸一化處理；再計算目標函數值的方差為0.1466，率定參數的方差分別為0.0396，0.0218，0.0849，證明PSO算法具有較高的穩(wěn)定性。

每次計算當目標函數達到終止條件所需的平均迭代次數為27.66次，最大次數為54次，最小次數為13次，因此PSO算法也有較高的效率。

由于在參數率定時，選用的資料序列較長，計算量較大（因為計算50次需要的時間較長）。鑒于PSO算法的穩(wěn)定性和效率，在實際應用中可以不用計算50次，只需計算5次左右即可，以期節(jié)省計算時間。另外，如果對精度要求不高，可以將終止條件中的目標函數之差設為ε=10-5，此舉也可節(jié)省計算時間。

6　結束語

本文將PSO算法用于馬斯京根法參數的自動率定中，開發(fā)了基于PSO算法的馬斯京根法參數自動率定系統(tǒng)，此系統(tǒng)具有較高的實用價值，適合內嵌入洪水預報系統(tǒng)進行參數率定和檢驗；該系統(tǒng)在雅礱江流域和豐滿流域等的洪水預報過程中均得到了應用，取得了較好的參數率定和檢驗效果。但是當率定參數采用的資料較長時，需要的時間也較長，下一步我們將再對PSO算法進行改進，以提高系統(tǒng)的運算效率。

[1]詹道江，葉守澤.工程水文學[M].北京:中國水利水電出版社，2000.

[2]包為民.水文預報[M].3版.北京:中國水利水電出版社，2006.

[3]翟國靜.馬斯京根模型參數估計方法探討[J].水文，1997，3: 40-43.

[4]魯帆，蔣云鐘，王浩，等.多智能體遺傳算法用于馬斯京根模型參數估計[J].水利學報，2007，38（3）：289-294.

[5]董純，康玲，侯國祥.基于遺傳模擬退火法的馬斯京根方程參數估計[J].華中科技大學學報（自然科學版），2002，30

（10）:20-22.

[6]李明明，李承軍，張銘.改進PSO法在馬斯京根模型參數估計中的應用[J].人民長江，2008，39（3）：60-62.

[7]李建勇.粒子群優(yōu)化算法研究[D].杭州:浙江大學，2004.

[8]熊鷹.微粒群算法的若干改進及應用[D].武漢:武漢理工大學，2006.

[9]李匡，付力，胡宇豐，等.馬斯京根法參數C0，C1，C2取值范圍的確定[J].南水北調與水利科技，2012，10（5）:43-45.

[10]劉蘇寧，甘泓，魏國孝.粒子群算法在新安江模型參數率定中的應用[J].水利學報，2010，41（5）:537-544.

[11]ShiY，EberhartRC．EmpiricalStudyofParticleSwarm Optimization[C]//Proceedingofthe1999CongressonEvolutionaryComputation．Piscataway，NJIEEEServiceCenter，1999．

[12]董曉華，劉冀，徐剛，等.目標函數對馬斯京根法河道流量演算精度的影響[J].水文，2009，29（6）：9-13.

TV734

B

1672-5387（2016）08-0088-04

10.13599/j.cnki.11-5130.2016.08.027

2016-06-29

李匡（1983-），男，高級工程師，從事洪水預報調度系統(tǒng)開發(fā)工作。