廢水中鳥糞石回收的機(jī)器學(xué)習(xí)預(yù)測和優(yōu)化

佟 穎, 蔣紹堅(jiān), 康冰艷, 冷立健, 李海龍

(中南大學(xué) 能源科學(xué)與工程學(xué)院, 湖南 長沙 410083)

0 引 言

磷(P)是一種不可再生的基本資源,對生命和農(nóng)業(yè)生產(chǎn)的發(fā)展至關(guān)重要[1]。鑒于P對環(huán)境的不利影響、儲量有限、全球需求增加和單邊依賴,P的可持續(xù)利用和管理日益重要。迄今為止,國內(nèi)外已經(jīng)開發(fā)了許多P重復(fù)利用和回收的工藝,如利用聚磷菌或反硝化聚磷菌等強(qiáng)化生物除磷[2-3]、向富磷溶液中添加試劑生成鳥糞石或磷酸鈣沉淀[4-6]和選用合適的固體吸附材料對磷溶液進(jìn)行選擇性吸附[7]等方法。在氮(N)和P共存的情況下,采用鳥糞石回收磷的方式是最常用的技術(shù)[9]。該方法只需調(diào)節(jié)廢水的pH和投加適當(dāng)鎂源即可進(jìn)行磷回收。含有高濃度N和P的廢水(如城市廢水、工業(yè)廢水和尿液等)是生成鳥糞石的可靠來源。

鳥糞石(磷酸銨鎂,MgNH4PO4·6H2O)是一種富含鎂(Mg)、銨和磷酸鹽的磷酸鹽礦物和緩釋肥料,具有吸磷量大、浸出率低、養(yǎng)分釋放緩慢等優(yōu)點(diǎn),可用于農(nóng)業(yè)生產(chǎn),也用于河流的營養(yǎng)強(qiáng)化[8]。使用鳥糞石結(jié)晶技術(shù)回收磷可以有效緩解礦物磷的資源匱乏問題,并且在減少當(dāng)?shù)氐沫h(huán)境污染方面顯現(xiàn)出巨大潛力。近年來,全球范圍內(nèi)已有近多個(gè)鳥糞石工廠從污水中回收營養(yǎng)物質(zhì)。例如,在比利時(shí)的一個(gè)污水處理廠中,鳥糞石用于回收磷的沉淀和緩釋肥料已經(jīng)商業(yè)化,并被應(yīng)用在OSTARA PEARL、AirPrex?技術(shù)和NuReSys?技術(shù)等流化床反應(yīng)器中[11]。美國某污水處理廠對現(xiàn)有的鳥糞石處理工藝進(jìn)行了質(zhì)量平衡的建模分析,結(jié)果顯示,在符合美國基礎(chǔ)設(shè)施典型設(shè)計(jì)參數(shù)的情況下,可以回收約71%至96%的潛在損失磷[12]。此外,LI等[13]對鳥糞石結(jié)晶的基礎(chǔ)研究進(jìn)展進(jìn)行了總結(jié),其中pH、溫度、Mg∶P、攪拌速率和晶種等因素被證明對除磷效率和產(chǎn)品純度具有重要影響。然而,廢水在不同回收工藝參數(shù)下的鳥糞石回收難以預(yù)測和操作,生產(chǎn)優(yōu)質(zhì)鳥糞石所需的條件尚未得到充分研究。

隨著人工智能的飛速發(fā)展,機(jī)器學(xué)習(xí)(ML,Machine Learning)作為一種新興的技術(shù),可以用來解決復(fù)雜的非線性問題,例如鳥糞石沉淀/結(jié)晶過程中化學(xué)反應(yīng)系統(tǒng)的建立。基于影響目標(biāo)的輸入變量與目標(biāo)本身之間的聯(lián)系,所構(gòu)建的ML模型可以準(zhǔn)確計(jì)算理論目標(biāo)值或預(yù)測目標(biāo)值。利用優(yōu)化算法對構(gòu)建的ML模型進(jìn)行求解,可以優(yōu)化輸入變量,以指導(dǎo)實(shí)驗(yàn)研究。據(jù)知,目前還沒有利用ML模型研究鳥糞石回收的相關(guān)研究。

因此,本文中采用了隨機(jī)森林(RF)和極限梯度提升算法(XGBoost)兩種基于ML算法的經(jīng)典預(yù)測模型,用于輔助從模型廢水中回收鳥糞石的預(yù)測和優(yōu)化。本研究的主要目的是通過ML提供一種省時(shí)省力的工藝設(shè)計(jì)和操作方案,以及相關(guān)的預(yù)測模型,從而促進(jìn)鳥糞石的回收利用。基于機(jī)器學(xué)習(xí)技術(shù)的應(yīng)用,本研究可以更好地指導(dǎo)設(shè)計(jì)和優(yōu)化鳥糞石回收工藝,為相關(guān)領(lǐng)域提供必要的信息和預(yù)測模型。

1 研究方法

1.1 數(shù)據(jù)收集與預(yù)處理

1.1.1 建立數(shù)據(jù)集

在學(xué)術(shù)數(shù)據(jù)庫(如Web of Science、Scopus和中國知網(wǎng)等)中搜索模擬廢水、磷回收、鳥糞石等關(guān)鍵詞。在25篇文獻(xiàn)中,整理了598條有效數(shù)據(jù)[14-38],數(shù)據(jù)集#1的輸入變量為攪拌速率(stir)、反應(yīng)時(shí)間(t)、反應(yīng)溫度(T)、pH、氮磷比(N∶P)、鎂磷比(Mg∶P)和P初始濃度(P_initial),目標(biāo)變量為P回收率(P_recovery);數(shù)據(jù)集#2來自于數(shù)據(jù)集#1中的422條數(shù)據(jù),其中輸入變量將Mg∶P替換成Mg∶N,輸出變量為N回收率(N_recovery);數(shù)據(jù)集#3來自于數(shù)據(jù)集#1中的422條數(shù)據(jù),輸出變量為N回收率和P回收率。

1.1.2 數(shù)據(jù)處理

為了使輸入和輸出參數(shù)服從正態(tài)分布并且比例適當(dāng),在訓(xùn)練ML模型之前,基于式(1)對收集的數(shù)據(jù)集進(jìn)行了標(biāo)準(zhǔn)化:

(1)

1.1.3 變量分析

根據(jù)式(2)所得的皮爾遜相關(guān)系數(shù)(PCC)計(jì)算兩個(gè)變量之間的線性相關(guān)性:

(2)

1.1.4 數(shù)據(jù)集劃分

在進(jìn)行ML建模之前,使用Python編程將整理好的數(shù)據(jù)集拆分為訓(xùn)練集和測試集。隨機(jī)選取80%的數(shù)據(jù)點(diǎn)作為訓(xùn)練集,每個(gè)數(shù)據(jù)集中剩余20%的數(shù)據(jù)點(diǎn)(測試集)用于對訓(xùn)練后超參數(shù)最優(yōu)的ML模型進(jìn)行最終評估[40]。

1.2 機(jī)器學(xué)習(xí)建模

1.2.1 機(jī)器學(xué)習(xí)模型與超參數(shù)調(diào)優(yōu)

本文選用RF和XGBoost兩種算法探究模擬廢水中鳥糞石的回收性能與攪拌速率、反應(yīng)時(shí)間、反應(yīng)溫度、pH、N∶P、Mg∶P和P初始濃度之間的關(guān)系,預(yù)測N回收率和P回收率。RF采用集成學(xué)習(xí)的思想,包含多個(gè)決策回歸樹,每個(gè)回歸樹都盡可能多地增長,通過對原始訓(xùn)練集進(jìn)行裝袋和自舉集成來解決回歸問題[41]。RF模型的最終結(jié)果將是上述所有決策樹的預(yù)測結(jié)果的平均值。在RF模型中,決策樹的數(shù)量(n_estimators)和樹的最大深度(max_depth)是兩個(gè)關(guān)鍵的超參數(shù),其他超參數(shù)是默認(rèn)值。XGBoost算法是基于決策樹的集成學(xué)習(xí)算法,采用梯度提升框架。它使用梯度提升技術(shù)以及幾個(gè)串聯(lián)開發(fā)的決策樹,每個(gè)后續(xù)決策樹都從前一個(gè)樹的錯(cuò)誤中學(xué)習(xí)和成長。當(dāng)后續(xù)決策樹足夠深或前一樹中不再有任何錯(cuò)誤模式時(shí),XGBoost學(xué)習(xí)將停止[42]。學(xué)習(xí)率(learning_rate)、n_estimators以及max_depth是模型訓(xùn)練過程中需要調(diào)整的基本超參數(shù)。這兩種ML模型具有可解釋性強(qiáng)、泛化能力強(qiáng)、訓(xùn)練時(shí)收斂速度快、對多特征數(shù)據(jù)適應(yīng)性強(qiáng)等優(yōu)勢。

交叉驗(yàn)證是提高模型性能和避免過度擬合的常見策略。本文采用5倍交叉驗(yàn)證的方法對RF模型和XGBoost模型的超參數(shù)進(jìn)行優(yōu)化,將訓(xùn)練集隨機(jī)分為五個(gè)部分。其中四部分用于模型訓(xùn)練查找ML模型的最佳超參數(shù),其余一部分將用于評估性能。此過程迭代了五次,以確保五個(gè)部分對于模型都是可行的。然后使用五次驗(yàn)證的平均性能來確定ML模型的最佳超參數(shù)。超參數(shù)調(diào)優(yōu)后,利用訓(xùn)練集對超參數(shù)最優(yōu)的模型進(jìn)行再訓(xùn)練,得到最優(yōu)ML模型。在調(diào)優(yōu)過程中,n_estimators和max_depth從2、4、8、16、32、64、128調(diào)優(yōu)到150。此外,XGBoost模型的超參數(shù)學(xué)習(xí)率值得注意,通過學(xué)習(xí)率可縮小每棵樹的貢獻(xiàn)。學(xué)習(xí)率太小會減慢收斂率,而學(xué)習(xí)率太高會導(dǎo)致梯度在最小值附近振動,甚至無法收斂。

1.2.2 模型評估

使用決定系數(shù)(R2)、平均絕對誤差(MAE)和平均均方根誤差(RMSE)來評估ML模型在訓(xùn)練和測試過程中的預(yù)測性能。R2越大,MAE和RMSE越小,表示訓(xùn)練后的ML模型預(yù)測越準(zhǔn)確,性能越好。R2、RMSE和MAE的計(jì)算方法如下:

(3)

(4)

(5)

1.2.3 特征重要性分析

成功建立具有最優(yōu)超參數(shù)的ML模型后,通過特征分析和部分依賴分析(PDP)對預(yù)測模型進(jìn)行解釋。用偏依賴圖表示必要輸入對輸出的具體平均影響。基于合作博弈論策略的加性解釋模型(SHapley Additive exPlanations,簡稱SHAP)[43]被用于ML模型的特征分析。具體而言,根據(jù)構(gòu)建的ML模型計(jì)算SHAP值,以評估輸入對預(yù)測輸出的相對貢獻(xiàn)和影響。

1.3 正向優(yōu)化

基于最優(yōu)的多目標(biāo)ML模型進(jìn)行優(yōu)化,以獲得在初始P濃度為10 mg/L和1 000 mg/L的情況下最高P回收率和N回收率所對應(yīng)的最佳工藝參數(shù)。在進(jìn)行優(yōu)化之前,先利用Origin軟件畫出數(shù)據(jù)集中各輸入變量的箱線圖,根據(jù)圖表剔除變量中的異常值,確定各變量的范圍。在模型內(nèi)按指定步長對各變量進(jìn)行迭代,得到大量解,篩選出多目標(biāo)ML模型中P回收率和N回收率權(quán)重各占0.5時(shí)的最優(yōu)解。例如,在P回收率和N回收率多目標(biāo)優(yōu)化過程中,N∶P從0.4到1.5,步長為0.1;pH從8到10,步長為0.1;Mg∶P從0.4到2,步長為0.1;t從0到90,步長為20;攪拌速率從40到420,步長為50;初始P濃度為10 mg/L和1 000 mg/L;溫度為25 ℃。通過整改項(xiàng)優(yōu)化篩選出最佳參數(shù)后進(jìn)行氮磷回收實(shí)驗(yàn)。

1.4 實(shí)驗(yàn)驗(yàn)證

1.4.1 儲備液的配制

稱取13.6 g的KH2PO4溶解于100 mL的去離子水中作為含P儲備液,此儲備液中P的濃度為1 mol/L。同時(shí)稱取5.3 g的NH4Cl和9.5 g的MgCl2分別溶解于100 mL的去離子水中,配制1 mol/L的含N儲備液和含Mg儲備液。實(shí)驗(yàn)中所用到的藥劑均為分析純。

1.4.2 實(shí)驗(yàn)過程

使用移液槍將定量的KH2PO4溶液和NH4Cl溶液依次移入100 mL的去離子水中,用NaOH(2 mol/L)和稀鹽酸(2 mol/L)調(diào)節(jié)溶液的pH到實(shí)驗(yàn)的設(shè)定值,之后根據(jù)指定的Mg∶P比立即加入定量MgCl2溶液,放到設(shè)定好參數(shù)的磁力攪拌器上攪拌。沉淀完成后,靜置后使用0.45 μm過濾器對上清液進(jìn)行過濾,采用鉬酸銨分光光度法分析濾液中P的含量,計(jì)算P回收率并與最佳ML模型的預(yù)測值比較。

2 結(jié)果與討論

2.1 數(shù)據(jù)集統(tǒng)計(jì)分析

表1顯示了數(shù)據(jù)集的參數(shù)分布情況。由表可知,數(shù)據(jù)集#1中P回收率在35%～100%之間,由于數(shù)據(jù)集#2和數(shù)據(jù)集#3中的可用數(shù)據(jù)比數(shù)據(jù)集#1的數(shù)據(jù)要少,相較于數(shù)據(jù)集#1大多數(shù)變量的范圍有所縮小。P回收率在58%～100%之間,而N回收率在5%～99%之間。差異是由7個(gè)輸入變量的取值造成的。建立和優(yōu)化以下模型,是識別和調(diào)節(jié)具有良好回收性能的鳥糞石的重要因素。

表1 數(shù)據(jù)集的參數(shù)分布

2.2 皮爾遜相關(guān)系數(shù)

為了評估輸入變量是否可以為ML模型開發(fā)提供有價(jià)值的信息,本研究進(jìn)行了皮爾遜相關(guān)系數(shù)(PCC)分析,分析結(jié)果如圖1所示。

圖1 皮爾遜相關(guān)性矩陣(PCC, **:顯著性水平 p <0.01,*:顯著性水平 p <0.05)Fig. 1 Pearson correlation matrix (PCC, **: significance level p<0.01, *: significance level p<0.05)

PCC的大小表示一個(gè)參數(shù)對另一個(gè)參數(shù)的影響程度,系數(shù)越大輸入與輸出的相關(guān)性越高,負(fù)值表示負(fù)相關(guān)。

如圖1(a)所示,數(shù)據(jù)集#1的PCC中輸入變量攪拌速率、P初始濃度、pH和N∶P與P回收率的相關(guān)性非常顯著(p<0.01),攪拌速率和N∶P呈負(fù)相關(guān);Mg∶P與P回收率的相關(guān)性顯著(p<0.05),呈正相關(guān)。如圖1(b)所示,數(shù)據(jù)集#2輸入變量中攪拌速率、t、N∶P和P初始濃度與N回收率的相關(guān)性非常顯著(p<0.01),其中P初始濃度呈正相關(guān)且最為關(guān)鍵;Mg∶N與N回收率呈顯著相關(guān)(p<0.05)。如圖1(c)所示,數(shù)據(jù)集#3中只有攪拌速率、T、N∶P和P初始濃度與P回收率顯著相關(guān)(p<0.01)。圖1(a)和(c)之間P回收率和其他變量的相關(guān)性差異主要是數(shù)據(jù)集大小和數(shù)據(jù)差異造成的,例如數(shù)據(jù)集#1和#3分別為598和422個(gè)數(shù)據(jù)點(diǎn)。

2.3 機(jī)器學(xué)習(xí)預(yù)測和模型解釋

2.3.1 超參數(shù)調(diào)整

圖2描述了在RF和XGBoost兩種算法下超參數(shù)n_estimators和max_depth的變化對模型預(yù)測效果的影響。模型基于訓(xùn)練集數(shù)據(jù)對超參數(shù)進(jìn)行調(diào)優(yōu),其中以R2作為評估標(biāo)準(zhǔn)。當(dāng)n_estimators高于20而max_depth高于10時(shí),XGBoost模型顯示出令人滿意的結(jié)果。模型的最佳n_estimators和max_depth參數(shù)分別為150、128、64和4、8、8。對于RF模型,當(dāng)max_depth大于8時(shí),隨著n_estimators的增加,R2沒有明顯變化。因此,選擇相對較小的max_depth可以避免預(yù)測模型的過度擬合。

根據(jù)上述超參數(shù)優(yōu)化過程,將RF和XGBoost這兩種模型預(yù)測結(jié)果最好時(shí)的超參數(shù)和模型性能結(jié)果匯總見表2。由此可知,XGBoost模型在單目標(biāo)和多目標(biāo)預(yù)測中表現(xiàn)出優(yōu)于RF的性能,因此本文根據(jù)XGBoost模型對數(shù)據(jù)進(jìn)行特征重要性分析。

表2 超參數(shù)調(diào)優(yōu)和模型性能結(jié)果匯總

2.3.2 P回收率單目標(biāo)預(yù)測

利用整個(gè)訓(xùn)練集對具有最佳超參數(shù)的XGBoost模型進(jìn)行再次訓(xùn)練,使用測試集進(jìn)行測試,結(jié)果如圖3(a)所示。模型XGBoost的最佳參數(shù)所對應(yīng)的測試集R2、RMSE和MAE分別為0.93、6.83和4.19。圖中平分虛線(y=x)表示預(yù)測值與實(shí)際值相等,點(diǎn)越靠近虛線表示預(yù)測性能越好。本次測試幾乎所有的預(yù)測和實(shí)際訓(xùn)練點(diǎn)都集中在平分虛線上。由于ML模型的固有屬性,ML模型對訓(xùn)練數(shù)據(jù)集的預(yù)測性能總是優(yōu)于測試數(shù)據(jù)集[44],所以圖中測試數(shù)據(jù)點(diǎn)的密度分布不如訓(xùn)練數(shù)據(jù)點(diǎn)緊密。

圖3 P回收率單目標(biāo)預(yù)測的預(yù)測性能(a)特征重要性分析(b)和部分依賴圖(c～f)Fig. 3 Prediction performance (a), feature importance analysis (b),and partial dependence graph (c～f) of P recovery single-target prediction

2.3.3 N回收率單目標(biāo)預(yù)測

由表2和圖4(a)可知,模型XGBoost最佳超參數(shù)所對應(yīng)的測試集R2、RMSE和MAE分別為0.91、8.30和3.73。圖4(b)中特征重要性分析表明,Mg∶N是影響N回收率最顯著的變量,其次是P初始濃度和N∶P。P初始濃度、Mg∶N以及N∶P決定了Mg和P對結(jié)晶的有效性,因此從圖4(c～e)中可以看出Mg∶N<1.0時(shí),Mg∶N的增加可以提高M(jìn)g對鳥糞石形成的有效性。高P初始濃度值和低N∶P值可以增強(qiáng)N回收率。同樣的,圖4(f)中pH對N回收率的影響趨勢與圖3(f)中pH對P回收率的影響趨勢相同,pH在9.5左右的N回收率最高,由此可看出P的回收與N的回收主要以鳥糞石結(jié)晶沉淀實(shí)現(xiàn)。

圖4 N回收率單目標(biāo)預(yù)測的預(yù)測性能(a)、特征重要性分析(b)和部分依賴圖(c～f)Fig. 4 Prediction performance (a), feature importance analysis (b),and partial dependence plot (c～f) of N recovery single-target prediction

2.3.4 N、P回收率多目標(biāo)預(yù)測

使用單目標(biāo)預(yù)測模型超參數(shù)的調(diào)優(yōu)方法對多目標(biāo)預(yù)測模型的超參數(shù)進(jìn)行調(diào)整。對P回收率和N回收率同時(shí)預(yù)測,XGBoost模型的最優(yōu)超參數(shù)所對應(yīng)的測試集R2、RMSE和MAE值分別為0.89、6.91和3.89,如圖5(a)所示。多目標(biāo)的預(yù)測結(jié)果遜色于單目標(biāo)預(yù)測。但由于多目標(biāo)預(yù)測數(shù)據(jù)集的變量范圍較窄,泛化能力更好。

圖5 N回收率和P回收率多目標(biāo)預(yù)測的預(yù)測性能(a)和特征重要性分析(b)Fig. 5 Prediction performance (a) and feature importance analysis (b)of multi-objective prediction of N recovery rate and P recovery rate

特征重要性分析(圖5(b))表明,P初始濃度是影響N、P回收率兩個(gè)目標(biāo)最顯著的變量,其次是pH、N∶P、t和Mg∶P,攪拌速率和T的影響最不顯著。

2.4 優(yōu)化與驗(yàn)證

根據(jù)多目標(biāo)預(yù)測的XGBoost模型篩選模型迭代得到的結(jié)果,可對N回收率和P回收率兩個(gè)目標(biāo)進(jìn)行正向優(yōu)化。在鳥糞石回收過程中,Mg含量是最需要考慮的因素,因此根據(jù)圖5(b)中各參數(shù)的SHAP圖,選出除Mg∶P以外的其他參數(shù)最優(yōu)值。根據(jù)Mg∶P的比值從0.4到2,步長為0.1,對迭代結(jié)果進(jìn)行篩選,P回收率和N回收率輸出的和為優(yōu)化目標(biāo),找到最優(yōu)解。對最優(yōu)解進(jìn)行后續(xù)的實(shí)驗(yàn)驗(yàn)證,最優(yōu)方案見表3。

表3 基于ML模型的多目標(biāo)正向優(yōu)化最優(yōu)方案

由表3可知,實(shí)驗(yàn)驗(yàn)證得到的結(jié)果與ML模型由表3可知,實(shí)驗(yàn)驗(yàn)證得到的結(jié)果與ML模型預(yù)測的最優(yōu)解相差不大,這意味著本文所建立的機(jī)器學(xué)習(xí)模型具有良好的應(yīng)用前景。

3 結(jié) 論

本文基于機(jī)器學(xué)習(xí),對從模擬廢水中回收鳥糞石這一問題進(jìn)行了預(yù)測與優(yōu)化。通過XGBoost和RF模型實(shí)現(xiàn)了P回收率和N回收率的單目標(biāo)和多目標(biāo)預(yù)測,其預(yù)測的R2分別為0.89～0.93和0.85～0.90。特征重要性分析表明,P初始濃度、N∶P、pH和Mg∶P對P回收率和N回收率有顯著影響。利用XGBoost多目標(biāo)模型成功預(yù)測和優(yōu)化了鳥糞石結(jié)晶的最佳條件,并進(jìn)行了實(shí)驗(yàn)驗(yàn)證。本研究證明了ML建模和實(shí)驗(yàn)驗(yàn)證相結(jié)合是一種促進(jìn)鳥糞石結(jié)晶研究的有前景的方法。