朱湘臨,陳 威,丁煜函,王 博,朱 莉,姜哲宇2,宋 彥
(1.江蘇大學(xué) 電氣信息工程學(xué)院,江蘇 鎮(zhèn)江 212013; 2.無(wú)錫太湖水務(wù)有限公司,江蘇 無(wú)錫 214000)
光合細(xì)菌(PSB)是具有光能合成體系的原核生物,在新能源開(kāi)發(fā)、環(huán)境保護(hù)、養(yǎng)殖、醫(yī)藥等方面有廣泛應(yīng)用[1-5]。然而其發(fā)酵過(guò)程中菌液紅色會(huì)不斷加深,越來(lái)越渾濁,是一個(gè)高度時(shí)變性、強(qiáng)耦合復(fù)雜過(guò)程。由此導(dǎo)致直接反應(yīng)發(fā)酵品質(zhì)的關(guān)鍵參量活菌濃度難以在線(xiàn)測(cè)量,進(jìn)而影響發(fā)酵過(guò)程的優(yōu)化控制。目前主要采用的是人工離線(xiàn)檢測(cè)的手段,但存在著較大的測(cè)量延時(shí)及誤差,且易染菌,會(huì)降低發(fā)酵品質(zhì)。軟測(cè)量技術(shù)[6-8]是解決此問(wèn)題的有效手段。
最小二乘支持向量機(jī)(LSSVM)[9-12]是一種基于統(tǒng)計(jì)學(xué)習(xí)理論的機(jī)器學(xué)習(xí)方法,具有優(yōu)秀的小樣本數(shù)據(jù)學(xué)習(xí)能力和預(yù)測(cè)能力,很適合樣本數(shù)據(jù)較少的發(fā)酵過(guò)程[13-15]。實(shí)踐表明,LSSVM的核參數(shù)與懲罰參數(shù)的取值對(duì)模型的泛化能力與精度起關(guān)鍵作用,不恰當(dāng)?shù)倪x擇可能會(huì)使LSSVM預(yù)測(cè)模型容易出現(xiàn)過(guò)度擬合和泛化能力不佳的問(wèn)題。針對(duì)這個(gè)問(wèn)題,大量學(xué)者針對(duì)LSSVM模型參數(shù)的選擇研究了各種優(yōu)化算法。文獻(xiàn)[16]將混沌優(yōu)化算法與LSSVM進(jìn)行有機(jī)結(jié)合并用于發(fā)酵過(guò)程中,得到誤差較小的發(fā)酵預(yù)測(cè)曲線(xiàn);文獻(xiàn)[17]采用粒子群算法優(yōu)化LSSVM的模型參數(shù),提高了模型的函數(shù)逼近能力,取得了良好的預(yù)測(cè)效果;文獻(xiàn)[18]以遺傳算法完成對(duì)LSSVM的模型參數(shù)選擇,能夠有效地跟蹤因變量和預(yù)測(cè)變量之間復(fù)雜的非線(xiàn)性關(guān)系,顯示出很好的預(yù)測(cè)與泛化能力。
蝙蝠算法(BA)是由楊新社根據(jù)蝙蝠的生物特性于2010年提出的一種新型智能算法[19],相較于粒子群算法、遺傳算法等傳統(tǒng)智能算法,BA算法具有模型簡(jiǎn)單、搜索能力強(qiáng)、收斂速度快等特點(diǎn),有效性和準(zhǔn)確性方面也有明顯的提高。目前已在電機(jī)運(yùn)行監(jiān)測(cè)[20]、數(shù)據(jù)挖掘[21]、風(fēng)速預(yù)測(cè)[22]、水力發(fā)電[23]等領(lǐng)域得到廣泛應(yīng)用。但BA算法仍有著易陷入局部極小的問(wèn)題[24-25],基于此,本文對(duì)BA算法的速度更新公式進(jìn)行改進(jìn),并將差分進(jìn)化算法(DE)的變異機(jī)制引入BA,提升了BA算法的全局及局部搜索能力。然后將改進(jìn)的蝙蝠算法(Improve Bat Algorithm,IBA)對(duì)LSSVM的模型參數(shù)進(jìn)行組合尋優(yōu),建立光合細(xì)菌發(fā)酵過(guò)程活菌濃度的IBA-LSSVM模型。仿真結(jié)果表明,該軟測(cè)量模型有不錯(cuò)的泛化能力和預(yù)測(cè)精度。
最小二乘支持向量機(jī)是Suykens[26]在支持向量機(jī)(SVM)的基礎(chǔ)上提出的新型支持向量機(jī),用于解決模型分解和函數(shù)估計(jì)問(wèn)題。其建模原理如下:
設(shè)有l(wèi)個(gè)訓(xùn)練樣本{(xi,yi)|i=1,2,…,l},其中樣本為n維向量,xi∈Rn為樣本輸入,yi∈Rn為樣本輸出,對(duì)樣本數(shù)據(jù)逼近的LSSVM可表示為:
(1)
s.t.yi=ωTφ(xi)+b+ξi,(i=1,2…,l)
(2)
式中,ω為權(quán)矢量;g∈R+是懲罰參數(shù);ξi是誤差變量;b為偏差量;φ(·)是一個(gè)非線(xiàn)性映射,能把xi從輸入空間映射到高維(甚至無(wú)限維)的特征空間。
用拉格朗日法對(duì)上述問(wèn)題進(jìn)行優(yōu)化:
(3)
式中,αi是拉格朗日乘子。根據(jù)KKT條件將求解的優(yōu)化問(wèn)題最終可轉(zhuǎn)換為求解線(xiàn)性方程:
(4)
式中,y=[y1,y2,…,yl]T;1l=[1,…,1]T;Il為l階單位矩陣;α=[α1,…,αl]T;K為滿(mǎn)足Mercer條件的核函數(shù)矩陣,K=K(xi,xj)=φT(xi)φ(xj),i,j=1,…,l。
本文選用徑向基函數(shù)作為核函數(shù):
(5)
式中,σ為核函數(shù)寬度。
最終得到LSSVM的函數(shù)估計(jì)為:
(6)
fi=fmin×(fmax-fmin)θ
(9)
(10)
(11)
式中,θ∈[0,1]是均勻分布的隨機(jī)數(shù);fi∈[fmin,fmax]是第i只蝙蝠的搜索脈沖頻率。
蝙蝠算法的局部搜索是通過(guò)隨機(jī)擾動(dòng)實(shí)現(xiàn)的,隨機(jī)在當(dāng)前最優(yōu)解集中選擇一個(gè)解zr,并在其附近產(chǎn)生一個(gè)局部新解znew,如式(12)所示:
znew=zr+μAt
(12)
式中,At為當(dāng)前蝙蝠種群的平均脈沖音強(qiáng);μ是[-1,1]上均勻分布的隨機(jī)數(shù)。
蝙蝠i的搜索脈沖頻度和脈沖音強(qiáng)更新公式為:
Ai(t+1)=βAi(t)
(13)
Ri(t+1)=R0[1-exp(-γt)]
(14)
式中,R0表示蝙蝠群體的最大脈沖頻度;γ表示頻度增加系數(shù);β表示脈沖音強(qiáng)衰減系數(shù);Ai(t)表示t時(shí)刻的脈沖音強(qiáng);Ai(t+1)表示t+1時(shí)刻的脈沖音強(qiáng);Ri(t+1)表示t+1時(shí)刻的脈沖頻度。
(15)
(16)
λ1=1-λ2
(17)
(18)
式中,M是當(dāng)前迭代次數(shù),Mmax是最大迭代次數(shù),τ是[0,1]的1個(gè)實(shí)數(shù),隨著迭代次數(shù)的增大,λ2線(xiàn)性遞減,而λ1則從1-τ線(xiàn)性遞增到1,qmax和qmin別為慣性權(quán)重的最大和最小值。
該自適應(yīng)方法使得算法在迭代早期具有較強(qiáng)的全局搜索能力,跳出局部極小,而在后期以當(dāng)前全局最優(yōu)解為主導(dǎo)進(jìn)行局部搜索。針對(duì)BA算法因缺乏變異機(jī)制而導(dǎo)致局部搜索能力不強(qiáng)的缺陷,將DE算法中的變異機(jī)制引入BA算法中,以豐富種群多樣性,強(qiáng)化算法局部搜索能力。本文采用標(biāo)準(zhǔn)DE/rand/1/bin[27]進(jìn)行變異操作,公式如下:
η=e1-Mmax/(Mmax-M+1)
(19)
F=F0×2η
(20)
znew=zr1+F(zr2-zr3)
(21)
式中,η為變異算子;F為縮放因子,F(xiàn)0為縮放因子的初始值;zr1,zr2,zr3是從當(dāng)前局部最優(yōu)解集中隨機(jī)選出的3個(gè)解。
綜上所述,IBA步驟如下:
Step 1:初始化種群數(shù)目N、最大迭代次數(shù)Mmax或搜索精度ε、搜索脈沖頻率范圍[fmin,fmax]、最大脈沖音強(qiáng)A、音強(qiáng)衰減系數(shù)β、最大脈沖頻度R0、頻度增加系數(shù)γ、τ、權(quán)重范圍[qmin,qmax]、初始縮放因子F0,初始化種群速度和位置并計(jì)算每只蝙蝠的適應(yīng)度,找到當(dāng)前全局最優(yōu)解z*。
Step 4:生成隨機(jī)數(shù)rand2,若rand2>Ai且更新位置后的蝙蝠i優(yōu)于當(dāng)前全局最佳位置,則接受此新位置為當(dāng)前全局最佳位置,并根據(jù)式(13)、式(14)更新脈沖頻度Ri和響度Ai。
Step 5:當(dāng)滿(mǎn)足搜索精度或達(dá)到最大迭代次數(shù)時(shí),則算法結(jié)束,輸出結(jié)果;否則轉(zhuǎn)回Step 2繼續(xù)迭代。
rand1與rand2是[0,1]上均勻分布的隨機(jī)數(shù)。
在光合細(xì)菌發(fā)酵過(guò)程中,影響到關(guān)鍵參量活菌濃度,且可實(shí)時(shí)在線(xiàn)檢測(cè)的潛在輔助變量有很多,如:光照強(qiáng)度E、空氣流量H、發(fā)酵罐壓力p、發(fā)酵液溫度T、發(fā)酵液體積V、氨水流加速率S、葡萄糖流加速率C、電機(jī)攪拌速度U、發(fā)酵液酸堿度pH等環(huán)境參量。軟測(cè)量模型的精準(zhǔn)度離不開(kāi)對(duì)可測(cè)輔助變量的合理選取[28-29],故采用一致關(guān)聯(lián)度法[30]獲得各參量與活菌濃度的關(guān)聯(lián)度,選擇關(guān)聯(lián)度最高的前5個(gè)作為輔助變量。表1為各潛在輔助變量與關(guān)鍵參量活菌濃度的關(guān)聯(lián)度。
表1 環(huán)境變量和主導(dǎo)變量的關(guān)聯(lián)度
由表1可知,光照強(qiáng)度E、發(fā)酵罐溫度T、空氣流量H、葡萄糖流加速率C、發(fā)酵液酸堿度pH與活菌濃度的相關(guān)性更密切,因此選擇這5個(gè)變量作為輔助變量。
采用離線(xiàn)訓(xùn)練構(gòu)建基于IBA-LSSVM的光合細(xì)菌發(fā)酵過(guò)程活菌濃度軟測(cè)量模型,在訓(xùn)練過(guò)程中采用IBA對(duì)LSSVM的懲罰參數(shù)g與核函數(shù)寬度σ進(jìn)行組合尋優(yōu),得到最佳的模型參數(shù),然后將訓(xùn)練好的軟測(cè)量模型進(jìn)行測(cè)試,驗(yàn)證IBA的優(yōu)化性能。整個(gè)軟測(cè)量模型構(gòu)建流程如圖1所示。
圖1 光合細(xì)菌發(fā)酵軟測(cè)量模型構(gòu)建流程圖
算法具體步驟如下:
1)采集樣本數(shù)據(jù),劃分為測(cè)試樣本和訓(xùn)練樣本兩部分,并進(jìn)行預(yù)處理。
2)設(shè)定參數(shù)g和σ的搜索范圍,初始化IBA算法的相關(guān)參數(shù)。
3)隨機(jī)產(chǎn)生蝙蝠種群,每個(gè)蝙蝠位置由g和σ構(gòu)成,使用訓(xùn)練樣本進(jìn)行訓(xùn)練,取訓(xùn)練樣本集的均方差作為適應(yīng)度值。
4)依據(jù)適應(yīng)度值最小法則,使用IBA算法進(jìn)行迭代搜索以及必要的越界處理,保留適應(yīng)度值最優(yōu)的位置。
5)判斷是否滿(mǎn)足搜索精度或達(dá)到最大迭代次數(shù),若滿(mǎn)足,結(jié)束迭代,轉(zhuǎn)至6);否則返回4)循環(huán)進(jìn)行迭代。
6)選擇全局適應(yīng)度最佳的模型參數(shù)gbest和σbest建立LSSVM軟測(cè)量模型,最終得到基于IBA-LSSVM的光合細(xì)菌發(fā)酵過(guò)程活菌濃度軟測(cè)量模型。
依據(jù)光合細(xì)菌發(fā)酵工藝進(jìn)行分批發(fā)酵試驗(yàn)。發(fā)酵罐進(jìn)行高溫消毒后,發(fā)酵過(guò)程溫度控制在25~34℃,電機(jī)攪拌轉(zhuǎn)速為350~450 r·min-1,罐壓為0.04~0.06 Mpa,通氣量為0.2~0.6 V·(V·min)-1,pH控制在6~8,光照強(qiáng)度控制在2000~5000 Lux,葡萄糖流加速率控制在0.01~0.015 g·(L·h)-1發(fā)酵周期為48小時(shí)。由數(shù)字系統(tǒng)每隔1分鐘采集發(fā)酵過(guò)程中光照強(qiáng)度、溫度、發(fā)酵液pH值、葡萄糖流加速率,由下位機(jī)傳送到上位機(jī),形成數(shù)據(jù)庫(kù)。在發(fā)酵正常工況下,每2 h取樣一次發(fā)酵液,經(jīng)離心、洗滌后,于105℃烘干后得到菌體含量。匯集現(xiàn)場(chǎng)數(shù)據(jù)進(jìn)行數(shù)據(jù)預(yù)處理,一共采集10個(gè)發(fā)酵批次的數(shù)據(jù),從上述批次中取出6批數(shù)據(jù)(包含480個(gè)樣本)作為訓(xùn)練樣本,2批(包含120個(gè)樣本)作為驗(yàn)證樣本,其余2批數(shù)據(jù)(包含480個(gè)樣本)作為測(cè)試樣本。
為了驗(yàn)證以上方法對(duì)光合細(xì)菌發(fā)酵過(guò)程軟測(cè)量建模的可行性采用IBA-LSSVM軟測(cè)量方法建立了光合細(xì)菌發(fā)酵過(guò)程活菌濃度的軟測(cè)量模型,選用標(biāo)準(zhǔn)BA-LSSVM與其做比較分析。模型參數(shù)設(shè)置為:種群規(guī)模N=50,最大迭代次數(shù)Mmax=300、搜索精度ε=0.05、搜索脈沖頻率范圍[0,10]、最大脈沖音強(qiáng)A=0.5、音強(qiáng)衰減系數(shù)β=0.95、最大脈沖頻度R0=0.5、頻度增加系數(shù)γ=0.9、τ=0.6、權(quán)重范圍[0.2,0.9]、初始縮放因子F0=1,參數(shù)g和σ的取值范圍均設(shè)置為[0.01,1000]。兩種優(yōu)化模型的預(yù)測(cè)效果分別如圖2和圖3所示。
圖2 BA-LSSVM預(yù)測(cè)效果圖
圖3 IBA-LSSVM預(yù)測(cè)效果圖
對(duì)比圖2和圖3可以發(fā)現(xiàn),圖3中IBA-LSSVM模型的預(yù)測(cè)曲線(xiàn)比圖2中BA-LSSVM模型的預(yù)測(cè)曲線(xiàn)更逼近于離線(xiàn)化驗(yàn)值,即IBA-LSSVM模型的預(yù)測(cè)效果明顯優(yōu)于BA-LSSVM模型的預(yù)測(cè)效果。將兩者對(duì)光合細(xì)菌活菌濃度的軟測(cè)量誤差在圖4中體現(xiàn)出來(lái),并計(jì)算兩種模型的均方誤差(RMSE)與最大相對(duì)誤差(MRE)如表2所示,以便能更直觀體現(xiàn)兩種模型的預(yù)測(cè)性能。
圖4 軟測(cè)量誤差
模型RMSEMRE/%BA-LSSVM0.31749.7IBA-LSSVM0.13584.1
從圖4和表2可以看出相較于BA-LSSVM,IBA-LSSVM的軟測(cè)量效果更優(yōu)異,MSE僅為4.1%表明模型的預(yù)測(cè)精度高,RMRE為0.1358表明模型的可行性佳。
針對(duì)光合細(xì)菌發(fā)酵過(guò)程中活菌濃度難以在線(xiàn)測(cè)量的問(wèn)題,以BA算法為基礎(chǔ),對(duì)BA的速度更新公式進(jìn)行改進(jìn),并引入DE的變異性,使得改進(jìn)的IBA算法有更好的局部收斂能力及收斂精度,將該方法確定LSSVM的模型參數(shù),建立了基于IBA-LSSVM的光合細(xì)菌發(fā)酵過(guò)程活菌濃度的軟測(cè)量模型。仿真實(shí)驗(yàn)結(jié)果表明,該軟測(cè)量方法的學(xué)習(xí)能力、預(yù)測(cè)性能優(yōu)于BA-LSSVM,為某些工業(yè)中存在難以在線(xiàn)測(cè)量的參量提供了一種可行的測(cè)量方法。