李祥蓉， 崔東文

(1.云南省水文水資源局文山分局 云南 文山 663000； 2.云南省文山州水務(wù)局， 云南 文山 663000)

1 研究背景

科學(xué)、客觀評價水功能區(qū)水質(zhì)狀況對于區(qū)域河長制考核評價、落實最嚴(yán)格水資源管理制度具有重要意義。開展水質(zhì)綜合評價方法研究一直是水資源評價領(lǐng)域的熱點和難點，目前灰色關(guān)聯(lián)法[1]、TOPSIS法[2]、神經(jīng)網(wǎng)絡(luò)法[3]、模糊法[4]、支持向量機(jī)法[5]、集對分析法[6]等方法已在水質(zhì)綜合評價中得到應(yīng)用。投影尋蹤(projection pursuit，PP)技術(shù)是將高維數(shù)據(jù)投影到低維空間，并在低維空間進(jìn)行數(shù)據(jù)分析來獲取最優(yōu)空間投影向量，目前已在水質(zhì)綜合評價[7-9]中得到廣泛應(yīng)用。研究表明，PP在實際應(yīng)用過程中，科學(xué)確定其最優(yōu)投影向量是客觀評價的基礎(chǔ)和關(guān)鍵。傳統(tǒng)遺傳算法(GA)、粒子群優(yōu)化(PSO)算法等普遍存在收斂速度慢、易陷入局部最優(yōu)等問題，很難搜索得到真正意義上的最優(yōu)投影向量[10]。目前一些新型群體智能算法嘗試用于投影尋蹤最優(yōu)空間投影向量的優(yōu)化，并取得了較好的優(yōu)化效果，如足球聯(lián)賽競爭(SLC)算法[11]、雞群優(yōu)化(CSO)算法[12]、蛾群算法(MSA)[13]、風(fēng)力驅(qū)動優(yōu)化(WDO)[14]算法等。

為能科學(xué)、客觀評價水功能區(qū)水質(zhì)狀況，進(jìn)一步拓展新型群體智能算法與投影尋蹤融合模型在水質(zhì)綜合評價中的應(yīng)用范疇，本文內(nèi)容安排如下：(1)介紹一種新型群體智能算法——靜電放電算法(Electrostatic Discharge Algorithm，ESDA)，選取6個標(biāo)準(zhǔn)測試函數(shù)對ESDA進(jìn)行仿真測試，并與當(dāng)前具有較好尋優(yōu)能力的新型群體智能算法——教學(xué)優(yōu)化(TLBO)算法、灰狼優(yōu)化(GWO)算法、鯨魚優(yōu)化算法(WOA)的仿真結(jié)果進(jìn)行比較。(2)以文山州盤龍河干流4個斷面2015-2017年共30個月水功能區(qū)水質(zhì)綜合評價分析為例，采用實測水質(zhì)數(shù)據(jù)和所選評價指標(biāo)在《地表水環(huán)境質(zhì)量標(biāo)準(zhǔn)》(GB 3838-2002)中Ⅰ～Ⅴ類分級閾值構(gòu)建水質(zhì)綜合評價投影目標(biāo)函數(shù)，利用ESDA求解該目標(biāo)函數(shù)最佳投影方向，構(gòu)建ESDA-PP模型對實例水功能區(qū)水質(zhì)進(jìn)行綜合評價及分析。

2 靜電放電算法及仿真驗證

2.1 靜電放電算法

靜電放電算法(Electrostatic Discharge Algorithm，ESDA)是Houssem[15]于2019年受靜電放電事件啟發(fā)提出的一種新型元啟發(fā)式優(yōu)化算法。ESDA假設(shè)電子設(shè)備適應(yīng)度與搜索空間位置有關(guān)，利用電子設(shè)備個體間直接或間接靜電放電現(xiàn)象，通過具有最低適應(yīng)度值電子設(shè)備向高適應(yīng)度值電子設(shè)備移動來改變其位置以獲得最佳空間位置，即問題最優(yōu)解[15]。與傳統(tǒng)群智能算法相比，ESDA設(shè)置參數(shù)少、收斂速度快、尋優(yōu)精度高，是一種具有較好競爭力的群體智能算法。

ESDA數(shù)學(xué)描述簡述如下[15]：

(1)初始化階段。在搜索空間中隨機(jī)生成群體大小n，即電子設(shè)備數(shù)量。每個電子設(shè)備由類似設(shè)計變量的不同組件組成。電子設(shè)備適應(yīng)度與搜索空間中設(shè)備的位置有關(guān)，適應(yīng)度值越大，設(shè)備的免疫力就越強(qiáng)。每個設(shè)備都有一個計數(shù)器，用于計算該設(shè)備作為損害設(shè)備的次數(shù)。

在搜索空間中，利用公式(1)隨機(jī)生成n個電子設(shè)備初始解：

xi=LB+rand(UB-LB)

(1)

式中:xi為電子設(shè)備群體中第i個電子設(shè)備的空間位置(i=1,2,…,n)；UB、LB分別為搜索空間的上、下界。

(2)在每次迭代過程中，從電子設(shè)備群體中隨機(jī)選擇3個電子設(shè)備的適應(yīng)度值并按降序排序(即最佳適應(yīng)度值個體排序為第1)。生成隨機(jī)數(shù)r1，如果r1>0.5，則靜電放電現(xiàn)象僅在兩個電子設(shè)備間發(fā)生；否則，涉及所有3個電子設(shè)備。每當(dāng)設(shè)備受損害或電子設(shè)備受到靜電放電時，該設(shè)備計數(shù)器將增加1。

①若r1>0.5，假設(shè)具有最低適應(yīng)度值的電子設(shè)備2朝著具有最高適應(yīng)度值的電子設(shè)備1移動。當(dāng)電子設(shè)備2移動到電子設(shè)備1附近時，假設(shè)在兩個個體之間發(fā)生靜電放電現(xiàn)象，且電子設(shè)備2受到損害，這相當(dāng)于直接靜電放電事件。電子設(shè)備2新位置更新公式如下：

x2new=x2+2β1(x1-x2)

(2)

式中：x2new為電子設(shè)備2新空間位置；x1、x2分別為電子設(shè)備1、2當(dāng)前空間位置；β1為服從正態(tài)分布的隨機(jī)數(shù)，其均值μ= 0.7，標(biāo)準(zhǔn)差σ= 0.2。

②若r1≤0.5，靜電放電涉及所有隨機(jī)選擇的3個電子設(shè)備，假設(shè)第3個電子設(shè)備正朝著另外兩個電子設(shè)備移動。當(dāng)電子設(shè)備3移動到電子設(shè)備1和電子設(shè)備2附近時，認(rèn)為發(fā)生靜電放電現(xiàn)象，且電子設(shè)備3受到損害，這相當(dāng)于間接靜電放電事件。電子設(shè)備3新位置更新公式如下：

x3new=x3+2β2(x1-x3)+2β3(x2-x3)

(3)

式中：x3new為電子設(shè)備3新空間位置；x1、x2、x3分別為電子設(shè)備1、2、3當(dāng)前空間位置；β2和β3均為服從正態(tài)分布的隨機(jī)數(shù)，其均值μ= 0.7，標(biāo)準(zhǔn)差σ= 0.2。

(3)假設(shè)電子設(shè)備受靜電放電超過3次，則認(rèn)為該電子設(shè)備已被損壞，必須更新，即在搜索空間中隨機(jī)生成新電子設(shè)備；否則，如果電子設(shè)備受靜電放電事件影響小于或等于3次，則產(chǎn)生隨機(jī)數(shù)r2，如果r2<0.2，則假定該電子設(shè)備的一部分已損壞并且必須改變，否則電子設(shè)備組件是安全的。

(4)將新增電子設(shè)備補(bǔ)充到原電子設(shè)備群體中并保存，再按照電子設(shè)備個體適應(yīng)度值的大小進(jìn)行降序排序，選擇前n個電子設(shè)備作為下一次迭代群體總量。

2.2 仿真驗證

為驗證ESDA尋優(yōu)能力，利用ESDA對文獻(xiàn)[16]中Sphere、Griewank、Quadric、Rastrigin、Schwefel 2.22和Ackley 6個典型測試函數(shù)進(jìn)行仿真驗證，并與當(dāng)前具有較好尋優(yōu)能力的新型群體智能算法TLBO、GWO和WOA的仿真結(jié)果進(jìn)行比較。其中函數(shù)Sphere、Quadric、Schwefel 2.22為單峰函數(shù)，用于測試算法的尋優(yōu)精度；函數(shù)Griewank、Rastrigin、Ackley為多峰函數(shù)，用于測試算法的全局搜索能力。4種算法對6個基準(zhǔn)測試函數(shù)重復(fù)進(jìn)行20次尋優(yōu)計算，并從平均值、標(biāo)準(zhǔn)差2個方面進(jìn)行評估，見表1。

驗證實驗參數(shù)設(shè)置如下：ESDA、TLBO、GWO和WOA 4種算法最大迭代次數(shù)T=100，群體規(guī)模N=100。其中，TLBO算法參數(shù)TF為1～10之間隨機(jī)整數(shù)；WOA對數(shù)螺旋形狀常數(shù)b=2；其他參數(shù)采用各算法默認(rèn)值。

表1 6個典型測試函數(shù)優(yōu)化結(jié)果對比

從表1可以看出，對于Sphere、Griewank、Quadric、Rastrigin、Schwefel 2.22函數(shù)，ESDA 20次尋優(yōu)均獲得了理論最優(yōu)值0，尋優(yōu)效果遠(yuǎn)優(yōu)于TLBO、GWO和WOA 3種算法；對于多峰函數(shù)Ackley，ESDA 20次尋優(yōu)均獲得相對理論最優(yōu)值8.88×10-16，尋優(yōu)精度高于TLBO、GWO、WOA算法6個數(shù)量級以上?？梢?，無論是單峰函數(shù)還是多峰函數(shù)，ESDA均表現(xiàn)出較為理想的尋優(yōu)效果，其具有出色的搜索平衡能力、極值尋優(yōu)能力和跳出局部極值能力。

3 ESDA-PP評價模型

3.1 投影尋蹤(PP)技術(shù)

(1)利用公式(4)、(5)分別對正向、負(fù)向指標(biāo)進(jìn)行處理。

x(i,j)=(x(i,j)-xmin(j))/(xmax(j)-xmin(j))

(4)

x(i,j)=(xmax(j)-x(i,j))/(xmax(j)-xmin(j))

(5)

式中：x(i,j)為指標(biāo)歸一化序列；xmax(j)、xmin(j)分別為第j個指標(biāo)值上、下限值。

(2)構(gòu)造投影值z(i)指標(biāo)函數(shù)：

(6)

式中:a為單位長度向量。

(3)將搜尋最優(yōu)投影向量問題轉(zhuǎn)化為單目標(biāo)非線性最優(yōu)求解問題，即：

(7)

式中：Sz為投影值z(i)的標(biāo)準(zhǔn)差；Dz為投影值z(i)的局部密度;Sz、Dz表達(dá)式參見文獻(xiàn)[12];Q(a)為待優(yōu)化目標(biāo)函數(shù)。

3.2 ESDA-PP模型評價實現(xiàn)步驟

ESDA-PP水質(zhì)綜合評價模型的實現(xiàn)步驟歸納如下：

從資產(chǎn)流轉(zhuǎn)順暢的角度來看，產(chǎn)權(quán)的完整程度主要體現(xiàn)在集體資產(chǎn)是否能夠自由進(jìn)出市場。與此同時，這也能夠更有效地實現(xiàn)要素的優(yōu)化與組合。當(dāng)集體資產(chǎn)順暢流轉(zhuǎn)時，集體資產(chǎn)的保值與增值效果才能達(dá)到最好。隨著產(chǎn)權(quán)的逐步明確，農(nóng)村經(jīng)濟(jì)產(chǎn)權(quán)的改革也得到進(jìn)一步完善，構(gòu)建與之相匹配的經(jīng)濟(jì)產(chǎn)權(quán)交易市場，相關(guān)的交易規(guī)則與程序也能夠逐步得到完善，最終實現(xiàn)公平交易[1]。

Step1選取水質(zhì)綜合評價指標(biāo)、序列長度和各評價指標(biāo)Ⅰ～Ⅴ類閾值，對樣本進(jìn)行一致性處理，利用處理后的樣本構(gòu)建水質(zhì)評價待優(yōu)化目標(biāo)函數(shù)Q(a)，利用ESDA對Q(a)投影向量進(jìn)行優(yōu)化求解。

Step2設(shè)置ESDA群體規(guī)模N、最大迭代次數(shù) 和算法終止條件；利用公式(1)生成電子設(shè)備初始空間位置，令當(dāng)前迭代次數(shù)t=0。

Step3隨機(jī)選擇3個電子設(shè)備，計算待優(yōu)化問題目標(biāo)函數(shù)適應(yīng)度值并按降序進(jìn)行排序。

Step4 在(0,1)之間生成隨機(jī)數(shù)r1，如果r1>0.5，則靜電放電現(xiàn)象僅在兩個電子設(shè)備間發(fā)生，利用公式(2)更新電子設(shè)備2位置；若r1≤0.5，涉及所有3個電子設(shè)備，利用公式(3)更新電子設(shè)備3位置。

Step5 檢查電子設(shè)備空間位置是否超出搜索范圍，若電子設(shè)備空間位置位于搜索空間之外，則利用公式(1)將其強(qiáng)制帶回搜索空間。

Step6 依據(jù)電子設(shè)備受靜電放電次數(shù)更新電子設(shè)備群體，并將新增電子設(shè)備補(bǔ)充到原電子設(shè)備種群中并保存。按電子設(shè)備個體適應(yīng)度值的大小進(jìn)降序排序，選擇前N個電子設(shè)備作為下一次迭代群體規(guī)模。

Step7 令t=t+1，判斷算法是否滿足終止條件，若是，輸出最優(yōu)解xbest，尋優(yōu)過程結(jié)束；否則重復(fù)Step3～Step7。

Step8輸出ESDA最優(yōu)解xbest，即最佳投影向量a。利用a計算各評價對象投影值z及各水質(zhì)類別分級投影值z′，利用水質(zhì)類別分級投影值z′對實例水功能區(qū)水質(zhì)進(jìn)行綜合評價及分析。

4 應(yīng)用實例

4.1 研究區(qū)概況與數(shù)據(jù)來源

盤龍河屬紅河流域，瀘江水系，發(fā)源于紅河州蒙自縣三道灣，流經(jīng)硯山、文山等6縣市，于天保船頭注入越南。盤龍河中越國界以上流域面積6 497 km2，河長247 km，落差1 803 m，平均坡降7.24‰，多年平均水資源量26.7×108m3。依據(jù)《文山州地表水功能區(qū)劃(2015修訂)》，盤龍河干流由上至下依次劃分為盤龍河硯山-文山農(nóng)業(yè)工業(yè)用水區(qū)，水質(zhì)代表斷面天生橋2020和2030年(規(guī)劃年)的水質(zhì)保護(hù)目標(biāo)為Ⅱ類；盤龍河文山飲用水源區(qū)的水質(zhì)代表斷面為灰土寨，規(guī)劃年水質(zhì)保護(hù)目標(biāo)為Ⅱ類；盤龍河文山景觀農(nóng)業(yè)用水區(qū)的水質(zhì)代表斷面為花橋，規(guī)劃年水質(zhì)保護(hù)目標(biāo)為Ⅲ類；盤龍河文山-麻栗坡保留區(qū)的水質(zhì)代表斷面為天保大橋，規(guī)劃年水質(zhì)保護(hù)目標(biāo)為Ⅲ類。開展盤龍河干流水功能區(qū)水質(zhì)綜合評價對于區(qū)域河長制考核評價、落實最嚴(yán)格水資源管理制度具有重要意義。

本文選取對水功能區(qū)水質(zhì)影響較大的DO、NH3—N、CODMn、BOD5和TP作為水質(zhì)綜合評價因子，通過構(gòu)建ESDA-PP模型對盤龍河干流天生橋、灰土寨、花橋、天保大橋4個斷面2015-2017年水功能區(qū)水質(zhì)進(jìn)行綜合評價，其中2015年為雙月監(jiān)測，2016、2017年為逐月監(jiān)測。限于篇幅，水質(zhì)監(jiān)測數(shù)據(jù)從略。

4.2 模型求解

依據(jù)ESDA-PP水質(zhì)綜合評價模型實現(xiàn)步驟，利用天生橋、灰土寨、花橋、天保大橋4個斷面2015-2017年共30個月斷面水質(zhì)實測數(shù)據(jù)和所選5個評價指標(biāo)在《地表水環(huán)境質(zhì)量標(biāo)準(zhǔn)》(GB3838-2002)中Ⅰ～Ⅴ類分級閾值構(gòu)建水質(zhì)綜合評價投影目標(biāo)函數(shù)Q(a)，利用ESDA求解Q(a)最佳投影方向量a。ESDA參數(shù)設(shè)置同上。經(jīng)求解，20次優(yōu)化目標(biāo)Q(a)值均為18 934.447，標(biāo)準(zhǔn)差均為0；20次優(yōu)化投影向量均相同，即a=(-0.0754 0.5530 0.4440 0.5530 0.4307)。該實例優(yōu)化結(jié)果再次驗證了ESDA具有較好的尋優(yōu)精度和極值尋優(yōu)能力。

4.3 結(jié)果與分析

(1)構(gòu)造水質(zhì)評價等級標(biāo)準(zhǔn)。利用上述優(yōu)化獲得的最佳投影向量a和經(jīng)一致性處理后的Ⅰ～Ⅴ類分級閾值，基于公式(6)計算得到Ⅰ～Ⅴ類分級標(biāo)準(zhǔn)投影值z′，z′即為水質(zhì)綜合評價分級標(biāo)準(zhǔn)，即：Ⅰ類<0.1801、Ⅱ類 ∈[0.1801,0.3303)、Ⅲ類 ∈[0.3303,0.5795)、Ⅳ類 ∈[0.5795,0.9411)、Ⅴ類 ∈[0.9411,1.4574)和劣Ⅴ類≥1.4574。

(2)綜合評價。同理，依據(jù)上述優(yōu)化獲得的最佳投影向量a和經(jīng)一致性處理后的天生橋、灰土寨、花橋、天保大橋4個斷面2015-2017年共30個月斷面水質(zhì)實測數(shù)據(jù)，基于公式(6)計算得到各斷面各月份綜合投影值z，利用水質(zhì)綜合評價分級標(biāo)準(zhǔn)z′進(jìn)行評價與分析，見表2及圖1～2。

表2 2015-2017年盤龍河4個斷面水質(zhì)綜合投影值及評價結(jié)果

圖12015-2017年30個月各斷面水質(zhì)評價結(jié)果 圖22015-2017年30個月各斷面枯、豐水期水質(zhì)評價結(jié)果

從表2及圖1、2可以得出以下結(jié)論：

(1)ESDA-PP模型對盤龍河4個斷面評價結(jié)果與單因子評價結(jié)果存在差異，總體評價結(jié)果要優(yōu)于單因子評價結(jié)果Ⅰ～Ⅱ類。從評價數(shù)據(jù)來看，影響盤龍河水質(zhì)主要因子為NH3—N，單因子評價法僅依據(jù)NH3—N單項指標(biāo)評價各斷面水質(zhì)類別，顯然不夠科學(xué)、客觀。

(2)依據(jù)ESDA-PP模型評價結(jié)果，并按頻次法大于等于80%即為達(dá)標(biāo)來看，天生橋、灰土寨、花橋和天保大橋斷面在研究期的30個月內(nèi)水功能區(qū)水質(zhì)達(dá)標(biāo)率分別為93.3%、83.3%、83.3%和100%，能滿足或基本滿足水功能區(qū)水質(zhì)保護(hù)目標(biāo)要求。其中天生橋斷面水質(zhì)處于Ⅰ～Ⅲ類之間，水質(zhì)保護(hù)目標(biāo)Ⅱ類，30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率為93.3%，水功能區(qū)水質(zhì)能滿足Ⅱ類水質(zhì)保護(hù)目標(biāo)要求；灰土寨斷面水質(zhì)處于Ⅰ～Ⅲ類之間，水質(zhì)保護(hù)目標(biāo)Ⅱ類，30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率為83.3%，水功能區(qū)水質(zhì)基本能滿足Ⅱ類水質(zhì)保護(hù)目標(biāo)要求；花橋斷面總體水質(zhì)處于Ⅱ～Ⅴ類，水質(zhì)保護(hù)目標(biāo)Ⅲ類，30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率為83.3%，水功能區(qū)水質(zhì)基本能滿足Ⅱ類水質(zhì)保護(hù)目標(biāo)要求；天保大橋斷面水質(zhì)處于Ⅰ～Ⅱ類，水質(zhì)保護(hù)目標(biāo)Ⅲ類， 30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率為100%，水功能區(qū)水質(zhì)完全能滿足Ⅱ類水質(zhì)保護(hù)目標(biāo)要求。若按單因子評價法進(jìn)行水功能區(qū)水質(zhì)達(dá)標(biāo)評價，則天生橋、灰土寨、花橋、天保大橋30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率分別為46.7%、56.7%、33.3%、100%，除天保大橋外，其余3個斷面均不能滿足水功能區(qū)水質(zhì)保護(hù)目標(biāo)要求。

(3)天生橋、灰土寨、花橋斷面枯水期水質(zhì)總體上要優(yōu)于豐水期，說明農(nóng)田和城市面源污染是影響盤龍河水質(zhì)的主導(dǎo)因素。天保大橋位于盤龍河下游，汛期來水較大，豐水期水質(zhì)總體上要優(yōu)于枯水期，說明盤龍河下游段水質(zhì)總體尚好的主要原因是巨大的水量和流量，也間接說明點源污染是影響盤龍河下游段水質(zhì)的主要因素。

5 結(jié) 論

本文介紹一種新型群體智能算法——靜電放電算法(ESDA)，通過對ESDA進(jìn)行仿真驗證和構(gòu)建ESDA-PP評價模型對盤龍河4個斷面2015-2017年30個月水功能區(qū)水質(zhì)進(jìn)行實例綜合評價，得出以下結(jié)論：

(1)選取6個標(biāo)準(zhǔn)測試函數(shù)對ESDA進(jìn)行仿真驗證，并與TLBO、GWO、WOA仿真結(jié)果進(jìn)行對比，結(jié)果表明ESDA尋優(yōu)精度遠(yuǎn)優(yōu)于TLBO、GWO、WOA，具有出色的尋優(yōu)精度和全局搜索能力，是一種具有較好競爭力的群體智能算法。

(2)ESDA-PP模型對實例評價結(jié)果顯示，盤龍河天生橋、灰土寨、花橋和天保大橋斷面2015-2017年30個月水功能區(qū)水質(zhì)達(dá)標(biāo)率分別為93.3%、83.3%、83.3%和100%，能滿足或基本滿足水功能區(qū)水質(zhì)保護(hù)目標(biāo)要求。天生橋、灰土寨、花橋斷面枯水期水質(zhì)總體上要優(yōu)于豐水期，農(nóng)田和城市面源污染是影響盤龍河水質(zhì)的決定因素；天保大橋豐水期水質(zhì)總體上要優(yōu)于枯水期，點源污染是制約盤龍河下游段水質(zhì)的主要因素。