一種多策略改進鯨魚優化算法的混沌系統參數辨識

潘悅悅，吳立飛，楊曉忠

（1.華北電力大學 控制與計算機工程學院, 北京 102206; 2.華北電力大學 數理學院信息與計算研究所, 北京 102206）

在非線性科學領域，混沌系統的控制與同步是研究熱點，其在信息科學、保密通信等方面應用廣泛[1-2]，準確的混沌系統模型是實現同步控制的基礎。然而，實際應用中，混沌系統結構復雜，部分參數不可知，因此精確辨識混沌系統參數具有重要的現實意義[3-4]。

近年來，許多學者開始關注和研究混沌系統參數辨識問題，并應用智能優化算法辨識系統參數。Chen等(2019年)[5]使用改進花授粉算法估計混沌和超混沌系統參數。數值模擬證明了新算法的有效性和魯棒性。Ahandani等(2020年)[6]提出一種混洗復雜進化算法解決混沌系統參數估計問題，數值結果表明改進算法提高了參數估計的收斂速度和精度。Turgut等(2021年)[7]使用單元拓撲改進鯨魚算法和正余弦算法，辨識混沌系統參數，測試結果證明，新算法成功地解決了參數識別問題。Ebrahimi等(2021年)[8]提出了一種改進洛茲映射混沌優化算法，數值結果表明，改進算法能夠高效、準確地估計混沌系統參數。

鯨魚優化算法(whale optimization algorithm，WOA)是一種新型群智能算法[9]，具有搜索能力強、計算穩定等特點，在最優控制、光伏系統等方面應用廣泛[10-11]。然而WOA存在全局和局部搜索不平衡，易陷入局部最優等問題，因此國內外學者對WOA進行了有效的改進。Chen等(2020年)[12]提出一種準對立混沌WOA(whale optimization algorithm with chaos mechanism based on quasi-opposition, OBCWOA)，數值仿真結果證明了改進算法的全局搜索能力。Chen等(2020年)[13]為提高WOA的收斂精度和速度，提出一種雙自適應隨機備用增強WOA(reinforced whale optimization algorithm, RDWOA)。仿真結果驗證了RDWOA的有效性。Chakraborty等(2021年)[14]為解決高維問題，結合多種策略，提出一種增強WOA(enhanced whale optimization algorithm, eWOA)，測試結果表明所提算法能有效解決高維問題。Shen等(2023年)[15]提出一種多策略進化WOA(whale optimization algorithm based on multi-population evolution,MEWOA)，試驗結果證明了MEWOA求解優化問題的有效性。Deng等(2023年)[16]提出了一種具有多策略混合算法的改進WOA(improved whale optimization algorithm, IWOA)，仿真結果表明，IWOA有較好的收斂速度和穩定性。

以上改進算法雖然有一定提升，但在收斂精度和速度方面仍有待提高，本文在已有工作的基礎上，以WOA為基礎，使用Chebyshev混沌映射產生均勻分布種群，增加初始種群多樣性。將收斂因子非線性化，增加自適應權重，兼顧全局搜索和局部挖掘。動態選擇自適應t分布或蟻獅優化算法[17-18]，跳出局部最優。通過對10個基準函數和高維測試函數進行仿真試驗，驗證了MIWOA的優越性。將MIWOA應用于Ro¨ssler和Lu¨混沌系統參數辨識，仿真結果證明了MIWOA辨識混沌系統參數的有效性。

1 多策略改進鯨魚優化算法的建立

1.1 基本WOA

WOA模擬了鯨魚捕食行為，其仿生學原理描述如下：

1) 包圍捕食階段。

鯨魚尋找食物時，當前鯨魚根據最佳鯨魚位置更新自身位置，公式為

式中：a=2-2t/T，t和T分別為當前迭代次數和最大迭代次數，X*(t)和X(t)分別為最佳鯨魚位置和當前位置，r1、r2為[0,1]的隨機數。

2) 螺旋更新階段。

在螺旋更新位置時，鯨魚通過螺旋向上的方式，靠近群體中最優位置，更新公式為

式中：D′為第i只鯨魚和獵物之間的距離，b為改變螺旋形狀的常數，l為[-1,1]中的隨機數。

在鯨魚搜索獵物時，收縮包圍和螺旋更新同步進行，位置更新公式為

式中p為[0,1]之間的隨機數。

3) 隨機搜尋階段。

鯨魚通過|A|大小，選擇隨機搜尋或者包圍捕食。當|A|＞1時，鯨魚在包圍圈外，通過隨機搜尋獲得獵物信息，位置更新公式為

式中Xrand為當前的一個隨機鯨魚位置。

1.2 WOA的改進

WOA初始化種群時采用隨機搜索策略，全局搜索能力弱。包圍捕食和螺旋更新階段，收斂因子線性遞減不能平衡全局和局部搜索。算法后期容易陷入局部最優。針對以上缺點，本文通過Chebyshev混沌映射初始化種群。包圍捕食和螺旋更新階段，非線性化收斂因子，加入自適應權重。算法后期，動態使用自適應t分布或蟻獅優化算法更新鯨魚位置，對WOA進行改進。

1) Chebyshev混沌映射。

Chebyshev混沌映射對初值敏感[19]，可產生大量無周期實值序列。本文采用Chebyshev混沌映射對WOA進行種群初始化，迭代方式為

式中k為階次。Chebyshev映射迭代1 000次的分布如圖1所示。

圖1 Chebyshev映射分布Fig.1 Chebyshev map distribution

由圖1可知，映射值分布于[-1,1]，Chebyshev混沌映射能更均衡地選取初始種群。將Chebyshev混沌序列映射到WOA解空間中，規定鯨魚種群數為N，隨機產生第1只鯨魚個體向量，Y=(y1,y2,···,yd)。使用式(10)對Y的各維進行n-1次迭代，產生其余的n-1只鯨魚。映射生成n只鯨魚個體為

式中：ud和ld分別為搜索空間第d維的上下界，yid和xid分別為第d維的第i只鯨魚和第i只鯨魚的坐標值。

2) 非線性收斂因子和自適應權重。

WOA中，參數A∈[-a,a]調節全局搜索和局部挖掘。當|A|≥1時，算法進行全局搜索，當|A|＜1時，算法進行局部挖掘。隨著迭代的進行，a=2-2t/T線性遞減不能體現實際優化過程，因此改進收斂因子，公式為[20]

式中：u為大于零的常數，調節a的衰減程度。T=500時，圖2為a隨u的變化曲線。

圖2 a的衰減曲線Fig.2 Attenuation curves of a

如圖2，與WOA相比，u越大，a＞1所占的迭代次數比例越大，算法全局搜索能力越強，局部挖掘能力越弱。反之，u越小，a＜1所占的比例越大，算法局部挖掘能力越強，全局搜索能力越弱。所以a的變化趨勢很大程度地影響優化求解。由于本文提出的Chebyshev混沌映射初始化種群策略，可以有效提高算法全局搜索能力，故選取u=0.6，增強局部挖掘能力，加快收斂速度，提高精度。

針對WOA后期局部挖掘時，權重固定，不利于算法尋優，本文提出一種自適應權重策略，公式為

3) 自適應t分布。

t分布又稱學生分布，概率密度函數為[17]

圖3 函數分布對比Fig.3 Function distribution comparison

由圖3可知，t分布的自由度n越小，曲線的雙尾翹越高，中間峰值越小，整體越平滑。反之，自由度n越大，中間峰值越大，整體越陡峭。自適應t分布對鯨魚位置的更新為

式中：xti為更新后的鯨魚位置，n和t(n)分別為迭代次數和以迭代次數為自由度的t分布函數。

4) 蟻獅優化算法。

在蟻獅優化算法中，螞蟻通過隨機游走更新位置，公式為[18]

式中：cumsum為螞蟻游走位置累加和，T為最大迭代次數，t為當前迭代次數，r(t)為0或1的隨機數。螞蟻隨機游走公式的標準化形式為

式中：ai和bi分別為第i個變量的最小值和最大值，和分別為第t代第i個變量的最小值和最大值。

2 混沌系統參數辨識的MIWOA

2.1 混沌系統的參數辨識原理

混沌系統參數辨識原理如圖4所示。

圖4 混沌系統參數辨識原理Fig.4 Principle of parameter identification for chaotic system

考慮如下n維混沌動力學系統：

式中：x=[x1x2···xn]T∈Rn為系統的n維狀態變量，x0為初始狀態量，θ0=[θ10θ20···θm0]T為參數真實值。當辨識系統參數時，假設結構為

式中：y=[y1y2···yn]T∈Rn為辨識系統的狀態變量，θ=[θ1θ2···θm]T為參數估計值。系統參數辨識是尋找一組最優未知參數，使系統真值x和估計值y的誤差值最?。?/p>

式中：xk和yk分別為k時刻的系統真值和估計值，m為參數辨識狀態變量序列的長度。目標泛函J(θ)的值越小，參數辨識的精度越高。MIWOA流程如圖5所示。

圖5 MIWOA流程Fig.5 Flow chart of MIWOA

2.2 MIWOA的時間復雜度分析

設種群規模為N，搜索空間維度為D，最大迭代次數為T，WOA的時間復雜度為O(NDT)。本文MIWOA以WOA為基礎進行改進，Chebyshev混沌映射初始化種群的時間復雜度為O(ND)，非線性收斂因子和自適應權重未在基本WOA基礎上增加循環嵌套，時間復雜度為O(NDT)，動態選擇自適應t分布或蟻獅優化算法的時間為t1，則其時間復雜度為O(NDT+t1)，雖然MIWOA的時間復雜度相對WOA有所增加，但在可接受范圍內，下面通過數值試驗說明MIWOA的卓越性。

3 MIWOA性能測試與分析

3.1 測試函數的選取

仿真試驗基于AMD Ryzen R7 5700U CPU@1.80 GHz，在Matlab R2019a 環境下運行。為了驗證MIWOA有更好的尋優性能，選取表1的10個基準函數進行測試，ε表示絕對誤差精度。f1～f7為單峰函數，f8～f10為多峰函數，10個測試函數的最優值均為0。

表1 基準函數Table 1 Benchmark functions

3.2 MIWOA與其他智能算法對比

采用10個基準函數檢驗MIWOA的性能，與WOA、粒子群算法(particleswarmoptimization，PSO)[21]、人工蜂群算法(artificialbee colonyalgorithm,ABC)[22]、灰狼算法(greywolf optimizer，GWO)[23]和哈里斯鷹算法(harris hawksoptimization，HHO)[24]進行測試結果比較。令N=30，D=30，T=1000，運行30次，表2為6種算法的最優值、最差值、平均值和標準差(黑色粗體為最優結果)。

表2 測試函數優化結果Table 2 Test function optimization results

由表2數據可知，對于單峰函數f1、f2，MIWOA的尋優結果均為0， WOA的計算結果優于PSO、ABC和GWO算法。在f3、f4測試函數上，MIWOA的計算結果優于HHO，可取得理論最小值。對于f5～f7函數，雖然MIWOA沒有達到最小值，但優于WOA。對于多峰函數f8～f10，MIWOA可收斂到最優值附近，其中f9函數可收斂到最優值，說明MIWOA相對于其他算法有更好的收斂精度和穩定性。

基準測試函數的收斂曲線可以清晰地展現算法的收斂速度，圖6(a)～(d)和(e)～(f)分別為單峰和多峰函數的平均收斂曲線。由圖6可以看出，對于6個測試函數，在收斂精度相同的情況下，MIWOA收斂速度更快，說明Chebyshev混沌映射初始化種群策略提高了種群中高質量個體的比例，非線性收斂因子和自適應權重平衡了全局和局部搜索，加速收斂。MIWOA收斂曲線波動下降，說明自適應t分布或蟻獅優化算法的動態選擇策略有助于算法跳出局部最優，提高收斂精度。

圖6 30維測試函數的平均收斂曲線Fig.6 Average convergence curves of 30 dimensional test functions

3.3 Wilcoxon秩和檢驗

Wilcoxon秩和檢驗能夠檢測更為復雜的數據分布，并與算法多次運行的數據對比，公平地體現MIWOA的優越性[24]。實驗設定顯著性差異為5% ，當p＜5%判定兩算法有明顯差異，反之無明顯差異。符號“＋”“-”和“＝”分別表示MIWOA的性能優于、劣于和相當于對比算法，N/A表示無法進行顯著性判斷。選取MIWOA在10個測試函數的運行結果與WOA、PSO、ABC、GWO以及HHO運行結果進行Wilcoxon秩和檢驗，計算p值。表3結果顯示大部分p＜5%，說明MIWOA的尋優能力優于5種對比算法。

表3 Wilcoxon秩和檢驗結果Table 3 Results of Wilcoxon rank sum test

3.4 MIWOA不同改進策略的有效性分析

為比較3種改進策略對MIWOA性能的影響，令N=30，D=30，T=500，對f1～f10進行尋優計算，將MIWOA與WOA、WOA1(采用Chebyshev混沌映射)、WOA2(采用非線性收斂因子和自適應權重)和WOA3(采用自適應t分布或蟻獅優化算法)比較，計算結果如表4。

表4 不同改進策略算法性能對比Table 4 Performance comparison of algorithms for different improved strategies

由表4可知，3種改進策略對WOA均有不同程度的提升，WOA3的改進效果最好，WOA2次之，WOA1的計算結果低于其他2種算法的計算結果，但對WOA仍有明顯的改進效果。將3種改進策略相結合時，算法搜索最精確，明顯優于WOA1、WOA2和WOA3，表明3種改進策略是有效的。

3.5 MIWOA與不同策略改進WOA對比

為了對比MIWOA與其他改進WOA的改進效果，令N=30，D=500，T=500，引用3種算法(OBCWOA[12]、eWOA[14]和MEWOA[15])的數據，在相同測試條件下，對f1～f10測試，結果如表5。

表5 MIWOA與其他改進WOA性能對比Table 5 Performance comparison between MIWOA and other improved WOAs

表5結果顯示，對于函數f1～f4，MIWOA與MEWOA可以在有限次迭代收斂到最優值，而OBCWOA和eWOA相對收斂精度較低。對于函數f5～f10，MIWOA可以使用更少的迭代次數收斂到最優值，說明MIWOA比其他改進WOA尋優效果更好。

平均絕對誤差(meanabsoluteerror, MAE)是評價算法有效可行的重要指標，計算公式為[23]

式中：mi為算法求解結果平均值，oi為各測試函數理論值，Nf為測試函數個數。表6為除理論值非零的f8函數外，4種算法的MAE排序，由計算結果可知，MIWOA的MAE值最小，排名第1，證明了本文改進策略的有效性。

表6 4種改進WOA的MAE排名Table 6 MAE ranking of four improved WOAs

3.6 MIWOA求解高維函數的實驗分析

由上述計算結果可知，對于低維測試函數，本文MIWOA尋優效果良好，但實際應用中高維大規模問題普遍存在，為了檢驗MIWOA求解高維問題的可行性，將其在10個基準函數，N=30，T=500，D=200、500、1000情況下，運行30次取平均值，根據文獻[25]求解算法尋優成功率。絕對誤差精度如表1所示，計算結果如表7。

表7 高維測試函數優化對比Table 7 High dimensional test functions optimization comparison

從表7可以看出，兩算法在求解高維測試函數時，尋優精度比求低維函數略有下降，但整體仍可達到較高精度，獲得理想結果。MIWOA在7個基準函數上尋優成功率達到了100%，相同維度下，MIWOA的尋優效果比WOA更好，說明MIWOA對于求解高維優化問題具有較高的計算精度和穩定性。

3.7 MIWOA平均運行時間分析

為了驗證MIWOA的計算速度，將其與PSO、ABC、GWO、HHO、WOA、WOA1、WOA2、WOA3運行30次的平均時間對比，N=30，D=1000，T=500，計算結果如表8。

表8 基準函數尋優平均時間對比Table 8 Comparison of average time for optimization of benchmark functions

由表8可知，對于元啟發式算法，WOA與ABC、HHO運行時間相近，可以快速收斂到最優值。對于改進算法，WOA、WOA1和WOA2的平均時長相當，說明WOA1和WOA2的改進策略未增加WOA運行時間。由于算法后期動態選擇自適應t分布或蟻獅優化算法，增加了循環嵌套，WOA3耗時多于WOA。融合多種改進策略的MIWOA搜索范圍變廣，尋優時間增加，但均在合理范圍內，與理論分析相符。

4 MIWOA辨識混沌系統參數

4.1 MIWOA對Ro¨ssler 混沌系統的參數辨識

Ro¨ssler 系統具有簡單、非對稱吸引子結構，在保密通信中承擔非常重要的作用。以典型Ro¨ssler混沌系統為例，驗證MIWOA可以精確辨識混沌系統參數。Ro¨ssler系統表達式為[7]

當參數a=0.2，b=0.2，c=5.7時，系統為混沌狀態，其演化過程如圖7所示。令初值向量為[-101]T，步長h=0.01，使用四階Runge-Kutta法解式(21)，利用前100個解及PSO、WOA、MIWOA辨識待定參數[abc]T，各算法獨立運行20次。

圖7 Ro¨ssler 混沌系統動力學軌跡Fig.7 Dynamic track of Ro¨ssler chaotic system

首先測試3種算法在具有一個未知參數的Ro¨ssler系統參數辨識中的搜索性能，即每次只辨識a、b、c3個參數中的一個，測試結果如表9。由表中數據可知，對于參數a，MIWOA可收斂到理論最優值，其最差值優于其他2種算法的最優值。對于參數b，WOA和MIWOA均可達到理論最優值，但MIWOA尋優效果更好。對于參數c，MIWOA的適應度值小于PSO和WOA的適應度值，估計值很接近真實值，說明MIWOA具有良好的全局搜索能力和計算魯棒性。

表9 不同算法的一維參數估計結果比較Table 9 Comparison of one-dimensional parameter estimation results for different algorithms

然后測試各算法在具有2個未知參數的混沌系統參數辨識中的搜索性能，計算結果如表10。從表10估計a和b可以看出，各算法搜索精度相似，相對于一個未知參數的情況，搜索精度略有下降，但MIWOA的計算結果仍優于PSO和WOA。另外2種情況也得出同樣的結論，驗證了MIWOA辨識2個未知參數Ro¨ssler系統的可行性。

表10 不同算法的二維參數估計結果比較Table 10 Comparison of two-dimensional parameter estimation results for different algorithms

最后，使用3種算法辨識具有3個未知參數的混沌系統，表11為各算法運行20次的統計結果。由表11中數據可知，各算法的搜索精度下降，但MIWOA的最差值優于PSO和WOA的最優值，在3種算法中搜索精度最高，表明MIWOA對Ro¨ssler混沌系統參數辨識更具有效性。

表11 不同算法的三維參數估計結果比較Table 11 Comparison of three-dimensional parameter estimation results for different algorithms

圖8和圖9分別是WOA和MIWOA對Ro¨ssler混沌系統的參數辨識曲線和適應度值曲線，從圖中可以看到MIWOA的參數估計值可以快速收斂到真實值，適應度值迅速逼近0，表明MIWOA高效的全局搜索能力和快速收斂速度。

圖8 WOA和MIWOA參數辨識曲線Fig.8 Parameter identification curves of WOA and MIWOA

圖9 WOA和MIWOA適應度值曲線Fig.9 Fitness value curves of WOA and MIWOA

4.2 MIWOA對Lu¨ 混沌系統的參數辨識

為了驗證MIWOA對混沌系統參數辨識的普適性，以Lu¨系統為例進行仿真，表達式為[26]

其中，系統參數真實值為a=30，b=22.2，c=8.8/3時，系統為混沌狀態，其演化過程如圖10所示。

圖10 Lu¨混沌系統動力學軌跡Fig.10 Dynamic track of Lu¨ chaotic system

令初值向量為[-101]T，使用PSO、WOA和MIWOA辨識Lu¨混沌系統3個參數未知時的情況，3種算法獨立運行20次，辨識結果如表12所示。

表12 各算法的三維參數估計結果比較Table 12 Comparison of three-dimensional parameter estimation results for each algorithm

由表12中數據可知，除了MIWOA的最優適應度值大于WOA的最優適應度值，MIWOA的計算結果均優于其他2種算法，說明MIWOA辨識混沌系統具有普適性。

5 結束語

本文以WOA為基礎，提出一種MIWOA。通過分析初始種群分布，使用Chebyshev混沌映射提高了初始種群質量。采用非線性收斂因子和自適應權重，提高了算法全局和局部搜索能力。動態選擇自適應t分布或蟻獅優化算法更新鯨魚位置，避免過早收斂。通過對10個基準函數和高維測試函數進行測試，以及Wilcoxon秩和檢驗，證明了MIWOA的優越性。

將MIWOA應用于Ro¨ssler和Lu¨ 混沌系統的參數辨識，仿真結果優于PSO和WOA，驗證了MIWOA辨識混沌系統參數的高效性。今后將繼續研究WOA的優化策略，將其應用到混沌系統控制與同步等其他領域。

致謝本文作者衷心感謝德國阿爾弗雷德韋格納研究所 Sergey Danilov 教授和王強博士在寫作過程中提出的寶貴建議和意見！