天然苦味分子識別及苦味閾值預測模型

馮寶龍，任海斌，段佳慧，張厚森，溫春輝，白曉森，高 飛，王玉堂,,4,

（1.東北農業大學現代教育技術中心，黑龍江哈爾濱 150030；2.乳品科學教育部重點實驗室（東北農業大學），黑龍江哈爾濱 150030；3.東北農業大學食品學院，黑龍江哈爾濱 150030；4.中國農業科學院農產品加工研究所，北京 100193）

苦味是人能感知的6種味覺之一，既能導致食品的風味變差，也能給食品帶來更加豐富的味感。啤酒、茶飲和咖啡等帶苦味的產品深受人們歡迎。人們對食物中的苦味組分知之甚少，如藍靛果、柚子、核桃等食物中苦味的呈味分子并不十分清楚[1]。這一方面給掩蓋苦味、改善風味帶來困難，也使這些苦味成分很少被食品研發人員用于開發新的產品。

目前鑒定苦味的方法都來源于藥物研究。主要通過志愿者來品嘗，配體驗證和構效關系模型（Structure Activity Relationships, SAR）的方法進行苦味驗證。志愿者品嘗通常要進行志愿者篩選及苦味標準化培訓，培訓周期一般較長[2]；配體驗證一般通過將待鑒定物質與苦味受體作用靶點進行分子對接，需要消耗較大的計算機資源[3]；而SAR的方法因其準確度高、簡便成為應用最多的方法。Rodgers等[4]使用包括649個苦味分子和13530個從MDL藥物數據倉庫（MDL Drug Data Repository, MDDR）中隨機選擇的分子作為數據集，利用圓形指紋（MOLPRINT 2D）和樸素貝葉斯算法對苦味化合物進行分類，在五折交叉驗證中，預測模型的準確度、精密度、特異性和靈敏度分別為88%、24%、89%和72%；Banerjee等使用從SuperSweet和BitterDB數據集中獲得的包含517種人工和天然甜味劑、685種苦味化合物的數據集，利用分子指紋和隨機森林分類器預測苦味，模型準確率為95%，ROC曲線下面積（Area Under Curve, AUC）為0.98[5-7]；Margulis等[8]使用包括極致梯度提升（eXtreme Gradient Boosting, XGBoost）在由169個非常苦的化合物和324個表現不是非常苦的化合物組成的493個化合物的訓練集上建立苦味程度識別模型，其準確率達到了87%。這些研究中，模型使用的苦味數據集包含較多非天然化合物，對食品中天然苦味分子的鑒定適用性尚未可知。到目前為止，還沒有針對性預測食品中苦味分子的數據集和SAR研究。

本研究以苦味分子研究文獻、專利及公開數據集為數據源，采用人工交叉驗證的方法搜集清洗數據，在Mysql數據集中建立最大的人工修正的苦味、非苦味數據集和苦味分子苦味閾值公開數據集。利用機器學習算法，首先建立定性構效關系模型，鑒別食品中的分子是否呈現苦味，然后，建立定量構效關系模型，對分子的苦味閾值做出預測，最后聯用上述模型預測FooDB數據集中潛在的苦味分子及其苦味閾值，并結合感官評價進行結果驗證。本研究建立了一種食品中苦味分子的快速鑒定方法，對于改善食品風味和開發新的苦味回甘的食品具有實際意義，以及對解釋苦味分子的結構特征、了解苦味受體及苦味信號傳遞奠定了數據基礎。

1 材料與方法

1.1 苦味-非苦味數據集及苦味閾值數據集建立

苦味和非苦味分子的數據都來自公開的數據集、文獻及專利。其中苦味分子及閾值來自Bitter-DB[6]中的天然化合物，并在此基礎上，擴展了從文獻、專利中獲得的天然苦味分子的閾值；非苦味分子來自FooDB[9]、FlavorDB[10]及Fenaroli的風味成分手冊（第5版）[11]。

非苦味分子數據集的選擇至關重要。非苦味分子數據集作為負集，既要有明確的文獻證明沒有苦味，也要在結構上盡量與苦味分子相似，這樣才能讓模型更加科學有效。根據苦味分子結構和物化特性，本文建立了如下規則來建立非苦味分子數據集：首先根據苦味分子的分子量范圍，篩選分子量范圍在85～1224之間的分子建立苦味分子數據集；其次苦味分子有較強的疏水性，根據苦味分子的疏水特性，選擇疏水性范圍與苦味分子相近的分子建立苦味分子數據集。

將不同來源的分子收集起來，通過Python語言，調用Pubchem API接口，查閱Pubchem數據集[12]獲得所有分子的簡化分子線性輸入規范（SMILES,Simplified molecular input line entry specification）格式，通過開源軟件（R語言（version 3.5.1），Python語言（version 3.6））和人工篩查方式去除重復、錯誤的分子，形成有意義的化學空間來訓練和評估機器學習模型。因苦味分子閾值的取值范圍較大，本研究進行適當縮放處理，將其取對數后，存入數據集。

1.2 分子描述符篩選

本研究利用分子操作環境（Molecular Operating Environment, MOE）[13]、ChemoPy[14]及Mordred[15]生成描述符。分別將分子的SMILES作為輸入，加載到描述符生成軟件中，輸出分子描述符。剔除含有缺失值的描述符后，利用R語言進行描述符優化和選擇。首先根據近零方差剔除描述符：設定頻數比率（freqCut）=25和唯一值比例（uniqueCut）=20兩個閾值，保留頻數比率小于25、唯一值比例大于20的描述符，刪除不具代表性的以及特殊的描述符，去除相關的變量，減少冗余；其次去除高共線性描述符變量：設定臨界值（cutoff）為0.95，刪除并只保留一個相似程度高于0.95的描述符，避免多個描述符描述同一特征；最后利用主成分分析（Principal Component Analysis, PCA）[16]進行描述符篩選：設定貢獻率閾值為0.5，保留貢獻率大于0.5的描述符。

1.3 苦味分子識別模型的建立和評價

建立苦味分子識別模型使用的工具是R語言[17]及其擴展包，包括：ggplot2[18]、RandomForest[19]、e1071[20]、kknn[21]、MASS[22]、sampling[23]。實驗采用了兩種算法：RF[24]和SVM[25]。RF算法是基于決策樹的分類器集成算法，其中每一棵樹都依賴于一個隨機向量，這些向量都是獨立分布的，通過生成多個分類樹，最終將分類樹結果進行匯總。通過十折交叉驗證和網格搜索，設置結點值（mtry）為9設置決策樹數目（ntree）為20構建隨機森林模型。SVM算法是在線性可分的情況下，在原空間尋找兩類樣本的最優分類超平面，在線性不可分的情況下，通過使用非線性映射將低維屬性空間的樣本映射到高維屬性空間使其變為線性情況，從而在該特征空間中尋找最優分類超平面。設置核函數為徑向基核函數（RBF核函數），gamma參數為0.1，懲罰系數（cost）為10建支持向量機模型。

所建立的苦味-非苦味數據集由139個苦味分子以及139個非苦味分子共有278個分子組成用于建立苦味分子識別模型，實驗采用無放回隨機分層抽樣，分為苦味與非苦味兩層，抽取訓練集和測試集，訓練集占總數據的3/4（208個分子），測試集占總數據的1/4（70個分子），訓練集用于模型訓練，而測試集用于模型驗證。使用準確度（Accuracy）和精確度（Precision）評估模型性能。

式中：TP-真苦，即真的苦味分子被預測為苦味；FP-假苦，即真的非苦味分子被預測成苦味，TN-真非苦，即真的非苦味被預測成非苦味，FN-假非苦，即真的苦味被預測成非苦味。

1.4 苦味分子閾值預測模型的建立和評價

建立苦味分子閾值預測模型使用的工具是R語言及其擴展包：caret、FactoMineR、factoextra、tidyverse、pls。實驗采用PLSR、RFR、kNNR、PCR四種算法進行苦味分子閾值預測。PLSR是常用的定量分析建模方法，它能有效地解決變量間的多重關聯問題。PLSR可以通過降維提取因子，設置提取的因子成分數（nt）為2，并將其作為回歸分析的目標。RFR是指通過集成學習的思想集成多棵樹的算法。它的基本單位是決策樹，基本上是一種集成學習方法，基于決策樹之間的距離進行回歸預測，設置結點值（mtry）為2，設置決策樹數目（ntree）為500。kNNR根據最接近某一未知點的k個數據點對該未知點進行回歸預測，設置鄰居個數（k）為5。PCR通過將一組高度相關的變量轉化為一組新的不相關的主成分變量來減少數據冗余來對未知點預測，設置主成分個數（ncomp）為5。實驗首先采用無放回隨機抽樣抽取訓練集和測試集，其中訓練集占總數據集的3/4，測試集占總數據的1/4，接著對實驗搜集到具有閾值的苦味分子經分子量和疏水性驗證，將驗證通過后的苦味分子用來建立苦味分子閾值預測模型。

苦味分子閾值預測模型建立之后，在苦味分子閾值預測模型中，使用五折交叉驗證（cross-validation）的方法評價模型的穩健性[26]。將訓練集劃分成5個互補相交的子集，每次選取其中一個子集做測試集，其余4個數據子集做訓練集構建模型，這個過程不斷重復，直到每個數據集的樣本都被用作測試集。通過計算模型的決定系數及誤差均方根，來評價模型的擬合優度和預測能力。當模型經過評價和驗證，證明其在統計學上具有穩健性和較好的預測能力后，該模型才可用于未知苦味分子的預測。

1.5 苦味分子及閾值驗證

對感官小組由12名評估員（5名女性和7名男性，年齡25～40歲）組成，所有參與人員均簽署感官評價知情同意書，并且沒有已知的味覺障礙病史。對于苦味的訓練和分類，MgSO4（166 mmol/L）溶液代表了一種短暫的金屬苦味品質，主要在舌頭的前部感知，水楊苷（1.4 mmol/L）則賦予了一種持久的苦味。主要在舌后部和喉嚨中感知到的苦味感和在口腔中提供持久苦味的咖啡因（8.0 mmol/L）被用作參考。感官分析在22～25 °C下進行，為了最大限度地減少任何有毒化合物的攝入，通過使用吞-吐法[27]進行感官分析，通過這種方法，測試材料不是吞咽而是吐出。苦味識別閾值由12名小組成員根據 ISO 4120[28]中詳述的方案通過三角測試確定。使用瓶裝水（pH4.5）作為溶劑和5 min的刺激間隔長度，在三個不同時間節點，將待測定呈味物質加水1:1稀釋成7個濃度梯度，每個濃度的樣品和另外兩個空白（水）組成一組，并用四位數字隨機編號呈現給受過訓練的感官小組，每人品評三次。要求小組成員將待測呈味物質在口腔中保持10 s后吐出。待7個梯度樣品品評完后，按濃度由低到高的順序整理數據，取不能感知到的濃度和上一個能感知到味覺的濃度的幾何平均值為品評人員的個人識別閾值濃度，苦識別閾值由所有個體閾值濃度的幾何平均值計算。

2 結果與分析

2.1 數據集的建立

所建立的苦味-非苦味數據集由139個苦味分子以及139個非苦味分子共278個分子組成用于建立苦味分子識別模型，苦味閾值數據集共有139個分子用于建立苦味分子閾值預測模型。苦味分子的閾值范圍在0.004～166之間，閾值分布近正態分布，說明用于建模的數據合理，適用于建立預測模型并采用顯著性分析比較模型間差異。非苦味分子數據集有139個非苦味分子，數據來源于FooDB、FlavorDB及Fenaroli的風味成分手冊（第5版）。如圖1展示了分子量與辛醇-水分配系數（Log octanol/water partition coefficient，log10P（o/w））的特征分布，log10P（o/w）為分子疏水性的物理參數，其分布反映了分子的疏水特性。苦味與非苦味分子在log10P（o/w）和分子量上的分布范圍基本一致，提高了模型的準確度。Rodgers等[4]使用的非苦味數據集是從MDL藥物數據集（MDDR）中隨機選擇的13530個假設不具有苦味的分子作為非苦數據集，與此相比，本實驗所用的非苦味數據集是經人工修正的數據集，經實驗驗證過的數據，降低了模型的噪聲，提高模型的實際應用能力；Tuwani等[29]使用甜味劑作為非苦味數據集，雖然數據集采用的是經實驗驗證過的數據，但兩種分子結構差異太大，對于結構相似的物質區分的準確性有待考證。本文搜集的數據集，其苦味、非苦味數據集貼近食品實際，增加了模型的實際應用價值。

圖1 分子量與辛醇-水分配系數特征分布圖Fig.1 Characteristics of molecular weight and octanol-water partition coefficients

2.2 分子描述符的選擇

為了降低冗余特征的影響，本文使用了近零方差、描述符共線性和主成分分析（PCA）的方法來篩選分子描述符。近零方差去除了無代表性特征和不相關特征，很多學者研究認為酚類化合物的分子量越小，苦味強度越高，對于小分子肽，其苦味與含有的疏水性氨基酸的數量和種類有關，通常苦味會隨著疏水性的增加而增強[30-31]，所以像以上兩種能夠有效描述苦味分子的特征描述符予以保留，反之剔除；描述符共線性去除并只保留一個描述同一特征的描述符；PCA線性地組合屬性，使得特征彼此正交，并捕獲數據的最大方差。以上方法都適用于訓練集分子對應的描述符。MOE 2D共生成206個描述符，有33個分子描述符被近零方差識別為與苦味分類預測無關，有89個分子描述符被共線性識別為彼此之間存在線性關系，不能區分各自對所要描述調整的具體作用，PCA優化后分子描述符集合中含有80個分子描述符。ChemoPy 2D共生成574個描述符，有230個分子描述符被近零方差識別為與苦味分類預測無關，有219個分子描述符被共線性識別為彼此之間存在線性關系，不能區分各自對所要描述調整的具體作用，PCA優化后分子描述符集合中含有97個分子描述符。Mordred 2D共生成1613個描述符，有415個分子描述符被近零方差識別為與苦味分類預測無關，有787個分子描述符被共線性識別為彼此之間存在線性關系，不能區分各自對所要描述調整的具體作用，PCA優化后分子描述符集合中含有34個分子描述符。如圖2所示，圖2a為Mordred 2D描述符篩選前熱圖，圖2b為Mordred 2D描述符篩選后的熱圖。其中紅色代表正相關，藍色代表負相關，經過描述符篩選后，紅色區域減少。結果顯示經過篩選后的描述相關性較低，表明這些描述符提供的化學信息在某種程度上是獨一無二的。

圖2 分子描述符篩選前（a）、后（b）的熱圖Fig.2 Heat map of the molecular descriptor before（a） and after（b） screening

高度的多重共線性會導致說明變量對因變量的單獨影響區分不開，變量間相互依存的強度和出現的共線性變量的重要性，會影響參數估計值[32]。PCA把多指標轉化為幾個綜合指標，降低觀測空間的維數，以獲取最主要的信息，通過少數幾個主成分最大限度地描述數據特點[16]。因此，通過以上優化可大大降低模型噪聲，提高模型準確度。

2.3 苦味分子識別模型

通過在苦味-非苦味數據集上進行分層抽樣，抽取3/4的數據（208個分子）作為訓練集，1/4的數據（70個分子）作為測試集建立苦味分子識別模型。如圖3所示，基于MOE 2D、ChemoPy 2D及Mordred 2D轉化的描述符數據建立的SVM苦味分子識別模型，以下分別簡稱MOE-SVM模型、ChemoPy-SVM模型及Mordred-SVM模型，其在測試集上的準確度（Accuracy）范圍分別為：0.900～1.000、0.886～1.000及0.929～1.000，其在測試集上的準確度（Accuracy）平均值為：0.968、0.965及0.979，其在測試集上的精確度（Precision）范圍分別為：0.857～1.000、0.800～1.000及0.943～1.000，其在測試集上的精確度（Precision）平均值分別為：0.963、0.944、0.985。基于MOE 2D、ChemoPy 2D及Mordred 2D轉化的描述符數據建立的RF苦味分子識別模型，以下分別簡稱MOE-RF模型、ChemoPy-RF模型及Mordred-RF模型，其在測試集上的準確度（Accuracy）范圍分別為：0.929～1.000、0.843～1.000及0.729～0.986，其在測試集上的準確度（Accuracy）平均值為：0.982、0.960及0.881，其在測試集上的精確度（Precision）范圍分別為：0.943～1.000、0.771～1.000及0.943～1.000，其在測試集上的精確度（Precision）平均值分別為：0.987、0.942及0.985。MOE-RF模型的準確度和精確度為0.982和0.987，均高于其他模型，說明MOERF模型可以較好地進行苦味分子識別。

圖3 苦味分子識別模型的評估結果Fig.3 Evaluation results of bitter molecular recognition model

2.4 苦味分子閾值預測模型

通過在苦味閾值數據集上進行無放回隨機抽樣，抽取3/4的數據（104個分子）作為訓練集，1/4的數據（35個分子）作為測試集建立苦味分子閾值預測模型，模型在測試集上的性能評估如圖4所示，采用MOE 2D描述符轉化的數據和基于RFR、kNNR、PCR和PLSR四種算法從而建立的苦味分子閾值預測模型，以下分別簡稱MOE-RFR模型、MOE-kNNR模型、MOE-PCR模型及MOE-PLSR模型，其決定系數（R2）分別為：0.80、0.75、0.77、0.83，其均方根誤差（RMSE）分別為0.53、0.58、0.53、0.45。采用ChemoPy 2D描述符轉化的數據和基于RFR、kNNR、PCR和PLSR四種算法從而建立的苦味分子閾值預測模型，以下分別簡稱ChemoPy-RFR模型、ChemoPykNNR模型、ChemoPy-PCR模型及ChemoPy-PLSR模型，其決定系數（R2）分別為：0.85、0.75、0.63、0.85，其均方根誤差（RMSE）分別為0.50、0.55、0.69、0.43。采用Mordred描述符轉化的數據和基于RFR、kNNR、PCR和PLSR四種算法從而建立的苦味分子閾值預測模型，以下分別簡稱Mordred-RFR模型、Mordred-kNNR模型、Mordred-PCR模型及Mordred-PLSR模型，其決定系數（R2）分別為0.81、0.72、0.68、0.83，其均方根誤差（RMSE）分別為：0.49、0.59、0.63、0.46，置信度為95%。

圖4 苦味分子閾值預測模型的評估結果Fig.4 Evaluation results of the prediction model of bitter molecule threshold

苦味分子在描述符轉化后，利用RFR、kNNR、PCR和PLSR四種算法進行模型構建。在RFR模型中，使用ChemoPy 2D描述符轉化數據建立的模型（ChemoPy-RFR）的擬合度最好，其決定系數（R2）為0.85，均方根誤差（RMSE）為0.50；在kNNR模型中，使用ChemoPy 2D描述符轉化數據建立的模型（ChemoPy-kNNR）的擬合度最好，其R2為0.75，RMSE為0.55；在PCR模型中，使用MOE 2D描述符轉化數據建立的模型（MOE-PCR）的擬合度最好，其R2為0.77，RMSE為0.53；在PLSR模型中，使用ChemoPy 2D描述符轉化數據從而建立的模型（ChemoPy-PLSR）的擬合度最高，其R2為0.85，RMSE為0.43。綜合以上結果得出：基于ChemoPy 2D描述符建立的ChemoPy-PLSR模型的擬合度最好，預測效果最準確，均方根誤差最低，模型穩定，其R2為0.85，RMSE為0.43，適合做苦味閾值預測模型。

2.5 苦味分子識別模型及苦味分子閾值預測模型應用

該研究為了預測潛在的苦味分子及其苦味閾值建立了苦味分子識別模型和苦味閾值預測模型，將上述模型應用于FooDB進行未知分子的苦味及其苦味閾值預測，并通過感官評價對對預測結果的苦味分子及閾值驗證。經本研究的最優模型Mordred-SVM模型預測，并將預測為苦味的分子使用ChemoPy-PLSR模型進一步對閾值進行預測，預測結果見“FooDB_results.csv” （https://gitee.com/wang_lab/BRATP），其中有5417個分子被預測為苦味。隨機挑選部分化合物的預測結果如表1所示：咖啡因、維生素B1、β-D-氨基葡萄糖、香橙素、葉酸。這些化合物的閾值分別為3.030、0.360、4.600、0.070、1.200 mmol/L。隨機挑選的這5種物質用于驗證模型的預測結果，雖然可能不是最具代表性的，結果如圖5所示，模型預測閾值與感官評價值擬合曲線決定系數為0.87，取得了較好的驗證結果。隨著大量研究的進行，在除這5種之外的分子上也可能取得較好的驗證結果。

圖5 苦味分子閾值預測擬合曲線Fig.5 Fitting curve of bitter molecular predict threshold and actual threshold

表1 化合物閾值預測結果Table 1 Prediction results of compound threshold

3 討論與結論

Luciana等[33]使用DRAGON plus v.5.0描述符和MobyDigs v.1.0 軟件，在遺傳算法（Genetic Algorithm, GA）方法下，使用多元線性回歸（Multiple Linear Regression, MLR）搜索最佳模型，即可變子集選擇-遺傳算法（VSS-GA）方法，得出了苦味與分子極性非直接相關且羥基和酯片段可以降低苦味的結論，由此可見，苦味強度與結構有關。但是Luciana Scotti等并未訓練苦味數據集，同時也未建立苦味分子結構與閾值的線性關系。苦味分子使用Mordred 2D轉化生成描述符從而建立的SVM模型的性能最好，其準確度和精確度分別為0.979和0.985；使用MOE 2D轉化生成描述符從而建立的RF模型的性能最好，其準確度和精確度分別為0.982和0.985。綜合以上結果得出：基于MOE 2D描述符建立的RF模型（MOE-RF模型）性能最好，其準確度和精確度分別為：0.982、0.985，最適合做苦味識別模型。與公布的最佳表現模型相比，MOE-RF苦味識別模型在準確性方面優于Zheng等[34]的e-Bitter模型和Banerjee等[7]的BitterSweet Forest模型，其中e-Bitter模型使用“Phytochemical Dictionary”數據集和1024bit-ECFP4和2048bit-ECFP6描述符，運用RF算法得到模型的準確率為85%～92%，本研究苦味識別模型的準確度較其提高了6%～13%；BitterSweet Forest模型使用苦甜數據集和Morganfeat分子描述符，運用RF算法得到模型的準確度為95%，本研究苦味識別模型的準確度較其提高了3%。由此可見，MOE-RF模型準確度高，穩定性強，更適用于苦味分子識別模型。且本研究的ChemoPy-PLSR苦味閾值預測模型對苦味分子結構及其閾值進行回歸分析，得出了較好的結果，該模型更精確地預測了苦味分子的閾值。

本研究基于MOE 2D、ChemoPy 2D和Mordred 2D三種描述符的優化及選擇，比較了由RF和SVM算法所建立的苦味分子識別模型，結果表明MOERF模型的準確度和精確度最高，其準確度和精確度分別為0.982、0.985。同時建立并比較了由RFR、PLSR、kNNR、PCR算法所建立的苦味閾值預測模型，結果表明，ChemoPy-PLSR模型的擬合度最好，均方根誤差最低，其決定系數和均方根誤差分別為0.85、0.43。由此可見，兩個模型具有良好的預測能力，可以用于苦味分子的分類預測及閾值預測。苦味識別模型使用了經實驗驗證過的數據，且對非苦味數據集進行了規范，較以往的分類模型更有可信度，準確度更高；本研究進一步實現了基于苦味分子結構預測苦味閾值的方法，模型擬合度較好，為苦味閾值預測提供了可行方法，具有一定的實際意義。所有數據和代碼存儲在https://gitee.com/wang_lab/BRATP，其他研究人員既可以利用本研究的代碼，繼續發掘其他苦味劑，也可以設計新的算法，獲得更為準確的預測結果。