基于分子指紋和拓撲指數(shù)的工質臨界溫度理論預測

引 言

臨界溫度(Tc)作為工質能維持液相的最高溫度，是建立狀態(tài)方程的基礎，也可以用于計算工質其他物性如焓、熵、比熱容、黏度、熱導率等。同時，臨界溫度是超臨界萃取過程中的重要參數(shù)。因此，獲取工質準確的臨界溫度具有重要的科學意義和工程價值[1-5]。實驗是獲取臨界溫度最有效的方式。然而由于實驗研究代價高昂、復雜性高，無法僅依靠實驗手段獲得工質的臨界溫度。因此，有必要提出一種能夠準確預測工質臨界溫度的理論模型。

臨界溫度的預測方法主要包括經(jīng)驗公式法、狀態(tài)方程法和定量結構-性質關系法(quantitative structure-property relationship, QSPR)。經(jīng)驗公式法采用一些易于測量的參數(shù)，如沸點、密度等，建立相應的關聯(lián)式得到臨界溫度。Reid等[6]最早提出了臨界溫度與沸點的關聯(lián)式Tc=1.5Tb。周傳光等[7]基于沸點與對比密度，提出了部分化合物臨界溫度的關聯(lián)式，平均偏差為1.35%。王新紅等[8]以沸點、對比密度、相對分子質量為參數(shù)，提出了一個新的有機物臨界溫度計算模型，平均偏差為2.36%。經(jīng)驗公式法形式簡單、計算精度較高，但缺乏理論基礎。狀態(tài)方程法可以基于pVT數(shù)據(jù)，擬合獲得工質狀態(tài)方程中相應參數(shù)，而后反推得到物質的臨界溫度。例如，Kontogeorgis等[9]采用狀態(tài)方程法估算了6種烷烴的Tc，絕對平均偏差均在2%以內(nèi)。Hsieh等[10]依據(jù)同樣的思路，首先獲得Peng-Robinson(PR)狀態(tài)方程的參數(shù)，進而得到392種純物質的臨界溫度，平均偏差為5.4%。狀態(tài)方程法需要已知工質pVT數(shù)據(jù)，且計算流程復雜，適用于密度數(shù)據(jù)較為豐富的物質。定量結構-性質關系法(QSPR)根據(jù)分子結構-物質性質之間的構效關系，對物質相關性質進行建模和預測?；鶊F貢獻法是QSPR中最常用的一種方法，包括經(jīng)典的Lydersen法[11]、Joback法[12]等。這些方法假設分子性質為各基團貢獻的線性加和，而基團貢獻度在不同分子中保持不變。這種線性加和的方法使用較方便，但沒有考慮不同基團的位置信息，導致該方法不能有效區(qū)分同分異構體。盡管后續(xù)的一些方法如Constantinou-Gani法[13]、Marrero-Pardillo法[14]等，通過引入多級基團和鍵貢獻在一定程度上緩解了上述缺陷，但適用范圍依然有限。綜合分析以上方法可知，現(xiàn)有模型無法對常見工質進行準確估算，須采用新的思路，以實現(xiàn)對包括同分異構體工質在內(nèi)的常見工質臨界溫度的精準預測。

分子結構描述符[如分子指紋(molecular fingerprints, MF)[15]、拓撲指數(shù)(topological index, TI)[16]等]作為一種將分子結構編碼為結構化數(shù)據(jù)的方法，可以將一種物質與其他物質進行明確區(qū)分。將分子描述符的概念引入QSPR模型，有望解決工質同分異構體的區(qū)分問題。在實際使用中，分子描述符通常與機器學習方法(machine learning, ML)相結合，以構建物質特性預測模型[17-19]。近年來，隨著計算機性能的不斷提高，有學者將分子描述符和機器學習應用于工質物性[20-24]的預測，預測效果良好。

本研究受上述分子描述符工作的啟發(fā)，首先以分子指紋表征分子結構，并借助機器學習算法建立16種臨界溫度的QSPR預測模型。此外，為了進一步提升本文模型的預測精度，本研究還將分子指紋與拓撲指數(shù)相結合，得到新的MF+TI-ML模型(以分子指紋和拓撲指數(shù)表達分子結構，結合機器學習算法建立模型)，以期取得良好的預測效果。

1 方 法

1.1 數(shù)據(jù)庫的搭建

本研究中工質的臨界溫度實驗數(shù)據(jù)取自物理性質設計研究所(DIPPR?801)[25]及相關文獻[26]。根據(jù)實驗數(shù)據(jù)不確定度對其進行篩選后，獲得了155種工質的Tc (本文所涉及工質的詳細信息，參見文末附錄)。搭建的數(shù)據(jù)庫中，臨界溫度的范圍分布在190.56~583.00 K。數(shù)據(jù)庫中工質可分為五種：烷烴、烯烴、鹵代烷烴、鹵代烯烴、醚類。為提升模型泛化能力，從每種類型工質中選取其中70%的數(shù)據(jù)點進入訓練集，剩下的30%作為測試集。訓練集用于建立工質臨界溫度的模型，測試集用于評估模型的預測性能。

1.2 分子指紋的生成

通過ChemDraw程序獲得工質分子的線性輸入規(guī)范(simplified molecular input line entry system, SMILES)，隨后利用在線轉換工具ChemDes [27]將SMILES字符串轉換為相同長度的二進制位串(即分子指紋)。為了研究不同長度/類型的指紋對QSPR模型性能的影響，本文選擇了計算四種分子指紋，包括兩種Key型指紋：MACCS(166位)和Pubchem(881位)，一種Path型指紋：Extended(1024位)和一種Circular型指紋：Morgan(2048位)。

1.3 回歸算法的選擇

本文使用了四種機器學習算法，包括支持向量回歸(support vector regression, SVR)、回歸樹(regression tree, RT)、隨機森林(random forest, RF)以及多層感知機(multi-layer perceptron, MLP)。

支持向量回歸通過核技巧將非線性數(shù)據(jù)映射到高維空間中，將非線性關系轉換為線性的形式，其精度取決于參數(shù)的選擇，例如核函數(shù)、寬度系數(shù)γ、不敏感損失系數(shù)ε、懲罰系數(shù)C等[28]。在本文中，將采用5折交叉驗證和網(wǎng)格搜索確定參數(shù)的最佳組合。決策樹(decision tree, DT)利用多節(jié)點的樹結構來描述各變量與目標之間的非線性關系，回歸樹是決策樹的回歸版本。由于樹模型具有較高的方差，可能導致結果不穩(wěn)定，基于樹模型的集成學習算法隨機森林相對樹模型有較大的改進[29]，在物性預測中應用較多。人工神經(jīng)網(wǎng)絡(artificial neural network, ANN)模擬神經(jīng)系統(tǒng)的結構，通過不斷調整神經(jīng)元間的權重和偏差，使整個網(wǎng)絡能更好地擬合數(shù)據(jù)[30-34]。多層感知機(MLP)是一種前饋神經(jīng)網(wǎng)絡，通過相互連接的人工神經(jīng)元和復雜的拓撲結構來模擬非線性關系[35]。本文利用深度學習庫Keras搭建了具有雙隱層的MLP，并通過試錯法確定了神經(jīng)元數(shù)、激活函數(shù)、學習率的最優(yōu)組合。

1.4 評估指標的選擇

本文選用均方根偏差(RMSE)、絕對平均偏差(AAD)、決定系數(shù)(R2)評估模型的預測性能， 相關定義式如下。

AAD=1m∑i=1m(100×∣∣yi?y?iyi∣∣)$AAD=\frac{1}{m}\sum _{i=1}^m\left(100\times \left(\frac{y_i-{\stackrel{?}{y}}_i}{y_i}\right)\right)$(1)RMSE=∑i=1m(yi?y?i)2m????????√$RMSE=\sqrt[]{\frac{\sum _{i=1}^m(y_i-{\stackrel{?}{y}}_i{)}^2}{m}}$(2)R2=1?????∑i=1m(yi?y?i)2∑i=1m(yi?yˉ)2????$R^2=1-\left(\frac{\sum _{i=1}^m(y_i-{\stackrel{?}{y}}_i{)}^2}{\sum _{i=1}^m(y_i{-\bar{y})}^2}\right)$(3)

式中，m表示樣本個數(shù)；yi$y_i$和y?$\stackrel{?}{y}$i 分別表示工質i臨界溫度的實驗值和預測值；yˉ$\bar{y}$表示臨界溫度實驗數(shù)據(jù)平均值。將評估指標應用于測試集時，RMSE、AAD越低，R2越高，模型的表現(xiàn)越好。

2 實驗結果與討論

2.1 模型的建立與評估

將四種分子指紋(MACCS、Pubchem、Extended、Morgan)分別用作四種機器學習算法(SVR、RT、RF、MLP)的輸入特征，得到16種臨界溫度的QSPR模型。各模型在測試集中的預測性能(以絕對平均偏差AAD為評價指標)如圖1所示。

圖1

圖1   以不同指紋為輸入的各QSPR模型的預測精度

Fig.1   Prediction accuracy of QSPR models with different fingerprints as inputs

從圖1可以看出，以MACCS指紋為輸入特征的模型預測性能較差，其中表現(xiàn)最好的MACCS-SVR(以MACCS指紋為輸入，結合SVR建立的模型)在測試集中的絕對平均偏差(AAD)也僅達到了7.49%。其原因是MACCS指紋長度過短，包含的結構信息有限，導致工質某些結構片段并不包含于MACCS指紋中。因此，以短位數(shù)的MACCS為輸入，模型預測精度并不高。

Extended指紋結合SVR算法建立的模型在測試集的AAD為7.36%。這是因為在ChemDes中，Extended指紋最大路徑長度默認設置為5 (即結構片段包含的最大鍵數(shù)為5)，導致許多線性路徑大于5的分子具有相同的Extended指紋。由于缺乏碳鏈長于5的工質Tc的實驗數(shù)據(jù)，如果提高路徑最大長度，特征維度會急劇增加，從而造成模型過擬合。因此目前來看Path型指紋不是建立工質QSPR模型的最優(yōu)選擇。

Circular型指紋Morgan作為一種立體型指紋長度最長，包含的結構信息也最多，因此可以有效地表征分子結構，進而有效區(qū)分工質同分異構體。綜合來看，雖然以Morgan指紋為輸入特征的模型預測性能要比上述兩種類型的指紋好，但仍不理想。其原因可能是位數(shù)過長導致了模型過擬合，因而Morgan指紋也不適用于搭建樣本數(shù)較少的QSPR模型。

Pubchem-MLP模型(Pubchem指紋結合MLP算法建立的模型)在訓練集、測試集的AAD分別為1.12%、4.76%。相比其他分子指紋而言，基于Pubchem指紋的QSPR模型預測表現(xiàn)最好。這說明Pubchem指紋可以合理表征工質的結構信息，從而在有限的訓練樣本中有效建立分子結構與臨界溫度之間的構效關系，準確預測工質臨界溫度。針對本文所研究的155種工質，Pubchem-MLP模型在工質臨界溫度實驗值和計算值的比較如圖2所示。

圖2

圖2   Pubchem-MLP模型在工質臨界溫度實驗值和計算值的比較

Fig.2   Comparison between experimental and caculated values of Pubchem-MLP model

從結果來看，四種ML算法建立的模型對工質臨界溫度的綜合預測效果排序如下：SVR > MLP > RF > RT。SVR模型預測精度最高且表現(xiàn)穩(wěn)定?；诩伤惴≧F的模型相比RT，在預測精度上有了明顯的提高，但和SVR仍有較大差距。

2.2 模型的優(yōu)化

Pubchem指紋可以很好地表達工質結構。但由于該類型指紋需要預先指定子結構，可能會造成工質中極少數(shù)同分異構體(如順反異構體)無法區(qū)分的問題。因此本文考慮在分子指紋的基礎上添加拓撲指數(shù)，以“分子指紋+拓撲指數(shù)”(MF+TI)作為新型分子結構描述符，采用效果較好的SVR和MLP算法，以期完全解決區(qū)分工質中同分異構體的問題。

拓撲指數(shù)是一種量化分子結構的指標，通過對表征分子圖的矩陣執(zhí)行數(shù)值運算獲得。這里引入拓撲指數(shù)(molecular topological index, MTI′)，在MTI′的基礎上添加幾何校正數(shù)(geometric modification, GM)區(qū)分工質中的同分異構體，拓撲指數(shù)S的計算公式[16]如下：

DvVw=DvDVDw$D_{vVw}=D_v{D}_V{D}_w$(4)MTI'=∑i=1N(vDvVw)i$MTI\text{'}=\sum _{i=1}^N{\left(v{D}_{vVw}\right)}_i$(5)GM=∑i=1N[MGF(DvVw+DTvVw)]i$GM={\sum _{i=1}^N\left(M_{GF}\left(D_{vVw}+D_{vVw}^T\right)\right)}_i$(6)S=MTI'+GM$S=MTI\text{'}+GM$(7)

式中， Dv、 DV、 Dw分別表示工質結構的價矩陣、頂點權重矩陣、鄰接矩陣；N表示分子的原子總數(shù)； v 表示價向量； MGF是用以區(qū)分異構體的對角矩陣。文末附錄給出了拓撲指數(shù)的具體計算流程和案例。

采用新型描述符后兩種模型的回歸和預測性能如圖3、圖4所示。可以看出引入拓撲指數(shù)S后，模型的預測精度明顯提升。Pubchem+TI-SVR模型(新型描述符輸入SVR算法建立的模型)在測試集的決定系數(shù)R2提高到0.8426，而Pubchem+TI-MLP模型(新型描述符輸入MLP算法建立的模型)在測試集的AAD降低至3.99%，R2提高到0.9143。對比圖2、圖4可以發(fā)現(xiàn)，相比Pubchem-MLP模型，Pubchem+TI-MLP模型預測性能明顯提高。這表明引入拓撲指數(shù)得到的新型描述符可以很好地解決區(qū)分工質中同分異構體的問題，提升模型的預測性能。

圖3

圖3   Pubchem+TI-SVR模型在工質臨界溫度實驗值和計算值的比較

Fig.3   Comparison between experimental and caculated values of Pubchem+TI-SVR model

圖4

圖4   Pubchem+TI-MLP模型在工質臨界溫度實驗值和計算值的比較

Fig.4   Comparison between experimental and caculated values of Pubchem+TI-MLP model

表1給出了本文搭建的Pubchem+TI-MLP模型在工質各數(shù)據(jù)集、各物質體系預測值和實驗值的AAD。從表中可以看出，新提出模型對烷烴類工質臨界溫度的回歸和預測都具有很高的精度，分別達到了0.90%和1.65%。模型對烯烴類工質的擬合回歸效果很好，但預測效果較差。醚類、鹵代烷烴類、鹵代烯烴類工質的計算精度相比上述兩類更低。從整個數(shù)據(jù)集來看，五種類型工質的絕對平均偏差均低于3%，取得了很好的計算效果。

表1   本文模型在各數(shù)據(jù)集、各物質體系的AAD

Table 1  AAD for each dataset and category of working fluids

圖5給出了155種工質的相對偏差(ARD)分布情況，相對偏差的定義如式(8)所示：

ARD=∣∣yi?y?iyi∣∣×100%$ARD=\left(\frac{y_i-{\stackrel{?}{y}}_i}{y_i}\right)\times 100\%$(8)

其中ARD<3%的工質有113種，占比72.9%，ARD>9%的工質僅7種，最大偏差為15.98%。

圖5

圖5   工質臨界溫度ARD分布情況

Fig.5   Distribution of ARD for Tc of working fluids

2.3 模型的對比

將新提出模型的計算結果與現(xiàn)有其他方法進行對比，用于對比的經(jīng)典方法列在表2中，模型的比較結果如表3所示。從表3中可以看出，本文模型的計算精度最高，Lydersen法和Joback法次之，C-G法精度最低?；诜悬c實驗值的Joback法計算工質Tc精度很高。但必須注意的是，并非所有工質都具有準確的沸點數(shù)據(jù)。當使用估算的沸點值(Testb$T_b^{est}$)計算時，Joback法的計算效果明顯降低。C-G法不需要使用沸點值，但該方法對工質臨界溫度的預測精度較低。K-R法將臨界溫度與分子量(Mw)、沸點關聯(lián)，并給出了一個簡單的線性方程，估算效果較好。然而，這種關聯(lián)缺乏理論基礎，普適性較差。綜合來看，本文提出的Pubchem+TI-MLP模型基于分子結構計算工質的臨界溫度，不僅無須沸點值，還獲得了最高的計算精度。

表2   現(xiàn)有預測臨界溫度的方法

Table 2  Existing method for estimation of critical temperature

表3   提出模型與以往方法計算效果的對比

Table 3  The comparisons between proposed model and previous methods

為了進一步驗證本文新提出模型和C-G法在區(qū)分同分異構體性能上的差異，表4給出了C-G二級基團貢獻法和本文模型在區(qū)分各類同分異構體(包括順反異構、位置異構和碳架異構)上的案例，其中T1cal$T_{cal}^1$表示本文模型計算值，T2cal$T_{cal}^2$表示二級C-G法計算值。

表4   C-G法和本文模型對同分異構體的區(qū)分案例

Table 4  Samples of isomers for the comparison of C-G method and the proposed model

從表4中可以看出，本文模型對于各類同分異構體的臨界溫度都取得了良好的預測精度。而C-G法對于所有順反異構體的臨界溫度預測結果完全一致，這表明C-G法無法區(qū)分順反異構體。

表5給出了本文模型和C-G二級基團貢獻法對155種工質中三類同分異構體計算結果的統(tǒng)計結果。從表5中可以看出，C-G法在碳架異構體的計算表現(xiàn)良好，37種碳架異構體臨界溫度的AAD為1.87%，但是由于不能辨別順反異構，對10種順反異構體的計算精度較差。在位置異構體臨界溫度的計算上，C-G法的精度也較低。而本文提出的Pubchem+TI-MLP模型不僅可以有效區(qū)分工質中存在的各類同分異構體，在計算精度上也遠高于C-G法。本文模型對順反異構體、位置異構體、碳架異構體臨界溫度計算值和實驗值的AAD分別為2.35%、2.51%、0.87%。

表5   C-G法和本文模型計算同分異構體的統(tǒng)計結果

Table 5  Statistical parameters of C-G method and the proposed model for the isomers

3 結 論

本文基于分子指紋和拓撲指數(shù)，采用機器學習算法建立了工質臨界溫度的Pubchem+TI-MLP模型。將新模型應用于155種常見工質的臨界溫度預測中，取得了良好的計算精度，針對測試集預測的絕對平均偏差為3.99%。通過與經(jīng)典模型的比較可以得出，新模型不僅可以有效區(qū)分工質中各類同分異構體，其計算精度相比現(xiàn)有模型也更高。通過對模型進一步分析還可看出，對指紋長度的選擇，必須綜合考慮樣本總數(shù)以及數(shù)據(jù)集包含的物質種類。在指紋類型的選擇上，Key型指紋Pubchem雖然在本文工質的臨界溫度預測上表現(xiàn)最好，但其自身不能區(qū)分少數(shù)順反異構體，需要引入拓撲指數(shù)以提高區(qū)分能力。長度更長的Path型和Circular型指紋對同分異構體的區(qū)分能力更好，但不適用于樣本數(shù)少的數(shù)據(jù)集。隨著以后工質實驗數(shù)據(jù)的不斷補充，可考慮使用更長的分子指紋搭建性能更加優(yōu)異的QSPR模型。