于天鑫 彭 璇
(北京化工大學(xué) 信息科學(xué)與技術(shù)學(xué)院, 北京 100029)
近年來(lái),甲烷作為一種清潔燃料越來(lái)越被人們所重視,而采用金屬有機(jī)框架材料(MOFs)實(shí)現(xiàn)甲烷的吸附[1-3]和儲(chǔ)存也引起了較為廣泛的關(guān)注。隨著實(shí)驗(yàn)室制備的MOFs以及計(jì)算機(jī)虛擬合成的MOFs的數(shù)量呈現(xiàn)爆發(fā)式的增長(zhǎng),僅僅利用巨正則系綜蒙特卡洛模擬(GCMC)方法[4-5]實(shí)現(xiàn)高性能吸附材料的高通量計(jì)算篩選已經(jīng)無(wú)法滿足要求。
基于GCMC的高通量篩選方法往往受限于龐大的MOFs數(shù)據(jù)庫(kù)和有限的計(jì)算資源,因此,具有強(qiáng)大數(shù)據(jù)分析和挖掘能力的機(jī)器學(xué)習(xí)方法被研究者們用來(lái)進(jìn)行高效的MOFs高通量篩選研究[6-8]?;诖?,本文采用機(jī)器學(xué)習(xí)建模的方法,通過(guò)決策樹(shù)(DT)模型及其衍生的隨機(jī)森林(RF)模型、極端隨機(jī)樹(shù)(ET)模型和梯度提升樹(shù)(GBDT)模型這4種模型對(duì)吸附甲烷的MOFs材料進(jìn)行高通量的計(jì)算篩選以選擇出最佳性能材料;對(duì)兩種較優(yōu)模型(RF模型和GBDT模型)的參數(shù)優(yōu)化進(jìn)行了探究,并推薦了合適的材料結(jié)構(gòu)特征參數(shù)。
目前,MOFs數(shù)據(jù)庫(kù)基本上可劃分為兩類,即由實(shí)驗(yàn)合成的MOFs(eMOFs)所組成的數(shù)據(jù)庫(kù)和由計(jì)算機(jī)合成的MOFs(hMOFs)所組成的數(shù)據(jù)庫(kù)。盡管通過(guò)計(jì)算機(jī)合成的hMOFs為MOFs的種類提供了無(wú)限的可能,但是hMOFs數(shù)據(jù)庫(kù)中的材料僅有一小部分能夠在實(shí)驗(yàn)中合成,絕大部分hMOFs設(shè)計(jì)的合理性和可行性存在著很大問(wèn)題,導(dǎo)致無(wú)法通過(guò)實(shí)驗(yàn)合成相應(yīng)的材料。
本文采用eMOFs數(shù)據(jù)庫(kù)[9-10],實(shí)驗(yàn)數(shù)據(jù)集中包含1 800個(gè)真實(shí)的MOFs數(shù)據(jù)樣本,其中每一種MOFs由9種特征描述符來(lái)表征,即表1中的前6種結(jié)構(gòu)描述符和后3種化學(xué)信息描述符。通過(guò)GCMC模擬計(jì)算每種材料在溫度298 K和壓力35 bar (1 bar=0.1 MPa)下的甲烷吸附量,根據(jù)美國(guó)能源局對(duì)吸附甲烷的金屬有機(jī)框架材料在該條件下的劃分標(biāo)準(zhǔn),將吸附量高于180(單位氣體吸附量與單位材料的體積比)的數(shù)據(jù)樣本標(biāo)記為高性能材料,反之,則標(biāo)記為低性能材料。
表1 每種材料特征向量的描述符表示Table 1 Descriptors used to construct a feature vector for each material
1.2.1相關(guān)性分析
本文計(jì)算了每個(gè)描述特征之間的相關(guān)性,如圖1所示。
圖1 特征向量的相關(guān)性Fig.1 Correlation of feature vectors
從圖1可以看出,材料的最大孔徑(MPD)和主導(dǎo)孔徑(DPD)的相關(guān)性非常強(qiáng),達(dá)到了99%。由此可見(jiàn),絕大多數(shù)材料的最大孔徑和主導(dǎo)孔徑是一致的。其次,可以看出每種材料的孔體積(Pv)和比表面積(Sa)的相關(guān)性也比較強(qiáng),達(dá)到91%,實(shí)際上,當(dāng)材料的孔徑較大時(shí),其相應(yīng)的比表面積也會(huì)增大,以支撐MOFs的有機(jī)骨架結(jié)構(gòu),從而更好地實(shí)現(xiàn)對(duì)甲烷的吸附。與此同時(shí),對(duì)于化學(xué)信息描述特征來(lái)說(shuō),它們之間的相關(guān)性都不高,而且與結(jié)構(gòu)描述特征的相關(guān)性也不強(qiáng)。鑒于兩者是從不同的角度對(duì)材料信息的提取,因此應(yīng)該結(jié)合結(jié)構(gòu)特征與化學(xué)信息特征共同完成材料的篩選。
1.2.2重要度分析
基于構(gòu)造決策樹(shù)時(shí)分裂節(jié)點(diǎn)的原理[11],進(jìn)一步計(jì)算每個(gè)特征描述符對(duì)甲烷吸附能力的重要度。在每棵樹(shù)的節(jié)點(diǎn)分裂時(shí)需要選擇該節(jié)點(diǎn)的分裂特征,通過(guò)計(jì)算基尼系數(shù)來(lái)確定節(jié)點(diǎn)特征,基尼系數(shù)越小,劃分的純度越高,則節(jié)點(diǎn)特征越好,特征的重要度就越高。樹(shù)的節(jié)點(diǎn)特征的順序就是重要度的順序。從圖2可以看出,MOFs材料的孔體積(Pv)對(duì)材料的吸附能力的重要度最高,這是因?yàn)椴牧系目左w積增大,甲烷的吸附量也會(huì)相應(yīng)增加。除此之外,結(jié)構(gòu)特征描述符對(duì)甲烷吸附的重要度較高,影響較大,而由于甲烷是非極性分子,材料的化學(xué)信息描述符對(duì)于甲烷吸附的重要度較小。因此,結(jié)構(gòu)特征對(duì)于甲烷吸附材料性能的影響更大。
圖2 特征向量對(duì)甲烷吸附的重要度Fig.2 Importance of feature vectors for methane adsorption
數(shù)據(jù)庫(kù)中的很多材料由于結(jié)構(gòu)原因?qū)е履承┨卣鳠o(wú)法測(cè)量,存在有缺省值問(wèn)題,此外當(dāng)按照分類標(biāo)準(zhǔn)劃分時(shí),存在高、低性能材料數(shù)量不平衡的問(wèn)題,極有可能造成數(shù)學(xué)模型的不穩(wěn)定。相比于其他機(jī)器學(xué)習(xí)的算法,由單棵決策樹(shù)衍生出的多棵決策樹(shù)是采用集成的學(xué)習(xí)方法,利用該方法建立模型對(duì)數(shù)據(jù)的要求相對(duì)較低,輸出的結(jié)果更加可靠。為了比較不同機(jī)器學(xué)習(xí)算法的篩選能力,本文選擇了決策樹(shù)基礎(chǔ)模型,及由它改進(jìn)而來(lái)的隨機(jī)森林、極端隨機(jī)樹(shù)和梯度提升樹(shù)3種樹(shù)模型,隨機(jī)地將數(shù)據(jù)集劃分為訓(xùn)練集和測(cè)試集兩組,采用普遍的7∶3的劃分方式,即訓(xùn)練集和測(cè)試集的材料數(shù)分別為1 260種和540種。利用不同的機(jī)器學(xué)習(xí)方法對(duì)訓(xùn)練集進(jìn)行學(xué)習(xí),并使用建立的模型對(duì)測(cè)試集的數(shù)據(jù)進(jìn)行篩選預(yù)測(cè)。
2.1.1混淆矩陣計(jì)算
通過(guò)模型對(duì)材料的測(cè)試集進(jìn)行篩選,計(jì)算各個(gè)模型的混淆矩陣[12-13]。從表2中各模型混淆矩陣的計(jì)算結(jié)果可以看出其分類效果,例如,GBDT模型在低性能材料的分類結(jié)果中,有375種材料分類正確,21種材料分類錯(cuò)誤;而在高性能材料的分類結(jié)果中,有135種材料分類正確,9種材料分類錯(cuò)誤。比較4個(gè)模型的混淆矩陣,發(fā)現(xiàn)它們的錯(cuò)誤分類數(shù)量大小順序?yàn)镈T>ET>RF>GBDT,GBDT模型的誤分個(gè)數(shù)明顯低于其他模型。
表2 4種模型的混淆矩陣Table 2 Confusion matrix for four models
2.1.2接收者操作特征(ROC)曲線
圖3給出了各個(gè)模型的ROC曲線,該曲線可以用來(lái)衡量模型的擬合程度[14]。由圖3可以看出,隨著誤診率的增加,召回率也逐漸增加。召回率T與誤診率F的計(jì)算公式如式(1)、(2)所示。
圖3 4種模型的ROC曲線Fig.3 ROC curves of four models
(1)
(2)
式中,TP表示樣本的真實(shí)類別是正例,并且模型將其預(yù)測(cè)成為正例的數(shù)量;FN表示樣本的真實(shí)類別是負(fù)例,并且模型將其預(yù)測(cè)成為負(fù)例的數(shù)量;TN表示樣本的真實(shí)類別是正例,模型將其預(yù)測(cè)成為負(fù)例的數(shù)量;FP表示樣本的真實(shí)類別是負(fù)例,模型將其預(yù)測(cè)成為正例的數(shù)量。對(duì)于每一個(gè)模型,我們希望其有一個(gè)較高的召回率以及較低的誤診率,所以圖3中每一個(gè)圖形的拐點(diǎn)越接近左上方則模型的效果越好,也即曲線與橫坐標(biāo)軸圍成的面積越大越好。DT、ET、RF以及GBDT這4個(gè)模型曲線與橫坐標(biāo)軸所圍成的面積分別為0.81、0.93、0.95和0.98。從面積上看,GBDT模型曲線的拐點(diǎn)更加靠近左上方,所圍成的面積最大,表明GBDT模型比其他模型的擬合效果更好。
2.1.3查準(zhǔn)率-查全率(PR)曲線
由于材料數(shù)據(jù)庫(kù)中低性能的材料較多,高性能的材料較少,這種較差的樣本均衡性會(huì)對(duì)模型的篩選造成一定的影響。因此,可以通過(guò)PR曲線來(lái)反映樣本均衡性對(duì)模型的影響[15]。4種模型的查準(zhǔn)率- 查全率曲線如圖4所示,查全率R以及查準(zhǔn)率P的計(jì)算公式如(3)、(4)所示。
圖4 4種模型的PR曲線Fig.4 PR curves of four models
(3)
(4)
可以看出,隨著查全率的不斷增加,查準(zhǔn)率則在不斷下降。對(duì)于一個(gè)較好的模型而言,應(yīng)該有較高的查全率及查準(zhǔn)率,即PR曲線的拐點(diǎn)盡量靠近右上方,使曲線與橫坐標(biāo)軸及左邊框圍成的面積越大越好。4種模型的PR曲線所圍成的面積大小順序?yàn)镈T 2.2.1測(cè)試集 基于DT、RF、ET和GBDT這4種機(jī)器學(xué)習(xí)模型對(duì)540種材料構(gòu)成的測(cè)試集進(jìn)行高性能甲烷吸附材料的篩選。從表3可以看出,利用4種機(jī)器學(xué)習(xí)模型篩選的類別為0的低性能材料,其各項(xiàng)指標(biāo)普遍比篩選出的類別為1的高性能材料要高,原因在于在訓(xùn)練集中進(jìn)行高低性能的分類時(shí),低性能材料的數(shù)量遠(yuǎn)多于高性能材料的數(shù)量,導(dǎo)致4種模型對(duì)于高性能材料的學(xué)習(xí)不充分,故而針對(duì)高性能材料篩選的效果不明顯。4種模型篩選的準(zhǔn)確度大小順序?yàn)镈T 2.2.2學(xué)習(xí)曲線 RF是基于套袋(bagging)的思想,有放回地均勻取樣,而GBDT則是基于梯度提升(boosting)的思想,根據(jù)訓(xùn)練錯(cuò)誤率對(duì)樣本賦予不同的權(quán)重。實(shí)驗(yàn)所選取的驗(yàn)證集是在數(shù)據(jù)訓(xùn)練進(jìn)行有放回抽取時(shí)未被抽取的數(shù)據(jù)的集合,這些未被抽到的材料數(shù)據(jù)稱作袋外數(shù)據(jù)[16]。繪制RF和GBDT這兩種較優(yōu)模型的學(xué)習(xí)曲線,如圖5所示。由圖可知,GBDT模型相對(duì)于RF模型的學(xué)習(xí)效果更好。在RF模型中,訓(xùn)練集的準(zhǔn)確度在訓(xùn)練過(guò)程中基本保持不變,說(shuō)明該模型在訓(xùn)練過(guò)程中擬合程度較好;而交叉驗(yàn)證集的準(zhǔn)確度則是從較低的數(shù)值逐漸上升的,且并沒(méi)有無(wú)限接近訓(xùn)練集的準(zhǔn)確度,兩者之間的間距較大,導(dǎo)致誤差比較大。也即在訓(xùn)練過(guò)程中,RF模型的擬合準(zhǔn)確度非常高,達(dá)到100%,但是在交叉驗(yàn)證過(guò)程中僅達(dá)到90%左右。這說(shuō)明RF模型對(duì)于新的數(shù)據(jù)集適應(yīng)性較差,存在過(guò)擬合的問(wèn)題。而對(duì)于GBDT模型,訓(xùn)練集的準(zhǔn)確度在訓(xùn)練過(guò)程中有微小的下降,而交叉驗(yàn)證集的準(zhǔn)確度則有所上升,且兩者有向同一準(zhǔn)確度值靠近的趨勢(shì)(兩條數(shù)據(jù)線趨近的準(zhǔn)確度值在95%左右)。由此可見(jiàn),GBDT模型能夠改善RF模型中存在的過(guò)擬合現(xiàn)象。 圖5 RF與GBDT模型的學(xué)習(xí)曲線Fig.5 Learning curves of RF and GBDT models 2.3.1RF模型參數(shù)曲線 圖6 不同參數(shù)對(duì)RF模型的影響Fig.6 Effect of different parameters on the RF model 2.3.2GBDT模型參數(shù)曲線 影響GBDT模型擬合效果的兩個(gè)最重要的因素分別為迭代次數(shù)n和每棵回歸樹(shù)的學(xué)習(xí)速率l,因此本文考察了這兩個(gè)因素對(duì)測(cè)試集偏差e的影響。從圖7可以看出,不同學(xué)習(xí)速率下曲線的變化趨勢(shì)大致相同,即隨著n的增加,e值是逐漸減小的。當(dāng)n小于20時(shí),3條曲線的e值下降得非???;而且l=0.5時(shí)對(duì)應(yīng)的偏差值是最低的,說(shuō)明在有限的迭代次數(shù)內(nèi),l越高,所達(dá)到的測(cè)試效果越好。而在n大于40時(shí),l=0.5曲線對(duì)應(yīng)的e值一直維持在較高的水平,而l=0.1和l=0.2曲線仍然有下降的趨勢(shì)。l=0.2曲線在迭代100次左右后偏差達(dá)到最低,而l=0.1曲線在迭代140次左右時(shí)偏差達(dá)到最低。這說(shuō)明當(dāng)回歸樹(shù)的l較低時(shí),要增加n的值才能保證e值降低。因此,在訓(xùn)練模型時(shí),可以調(diào)節(jié)回歸樹(shù)的學(xué)習(xí)速率l和迭代次數(shù)n兩個(gè)參數(shù)來(lái)改善GBDT模型。在本文測(cè)試中,當(dāng)n=100,l=0.2時(shí),可以達(dá)到有效改善模型性能的效果。 圖7 不同參數(shù)對(duì)GBDT模型的影響Fig.7 Effect of different parameters on the GBDT model 2.3.3高性能吸附材料的特征向量 對(duì)影響甲烷吸附量的重要度進(jìn)行分析發(fā)現(xiàn),影響甲烷氣體吸附的主要因素為材料的孔體積、密度、限制孔徑及最大孔徑。利用GBDT模型篩選測(cè)試集內(nèi)的高性能材料,分析前20種高性能材料的特征向量與甲烷吸附量之間的關(guān)系,結(jié)果如圖8所示。從圖中可以看出,當(dāng)孔體積為0.5~0.75 cm3/g,限制密度為2~3 g/cm3,材料孔徑在4 ?左右,最大孔徑在6~10 ?時(shí),甲烷的吸附量較高。 圖8 高性能材料的特征向量與甲烷吸附量的關(guān)系Fig.8 Relationship between the feature vectors and methane adsorption of high-performance materials 本文采用DT模型及其衍生的RF、ET、GBDT模型對(duì)金屬有機(jī)框架材料進(jìn)行分類篩選,通過(guò)對(duì)模型的篩選性能進(jìn)行比較,得出GBDT模型的篩選效果最好。當(dāng)?shù)螖?shù)為100,學(xué)習(xí)速率為0.2時(shí),GBDT的模型性能最佳。利用GBDT模型篩選出的前20種金屬有機(jī)框架材料進(jìn)行構(gòu)效關(guān)系分析,得出當(dāng)孔體積為0.5~0.75 cm3/g,材料密度為2~3 g/cm3,材料限制孔徑在4 ?左右,最大孔徑在6~10 ?時(shí),甲烷的吸附量較高。所得結(jié)果可望為用于甲烷吸附的金屬有機(jī)框材料的設(shè)計(jì)提出合理化建議。2.2 模型測(cè)試結(jié)果
2.3 模型參數(shù)討論
3 結(jié)論