基于決策樹算法的熱毒寧注射液金銀花青蒿醇沉過程質量控制研究

丁 鴻，徐芳芳，杜 慧，張 欣,3，徐 冰，吳 云,3，王振中,3，肖 偉,3*

1.南京中醫藥大學，江蘇 南京 210023

2.江蘇康緣藥業股份有限公司，江蘇 連云港 222001

3.中藥制藥過程新技術國家重點實驗室，江蘇 連云港 222001

4.北京中醫藥大學 中藥信息學系，北京 102400

熱毒寧注射液是由金銀花、青蒿和梔子3 種中藥提取精制而成[1-2]，其前提取精制過程包括金銀花提取濃縮、青蒿提取濃縮、梔子提取濃縮、金銀花青蒿（金青）醇沉、金青萃取、梔子萃取等工序，均在中藥數字化提取精制工廠中完成[3]。中藥數字化提取工廠配備集散控制系統（distributed control system，DCS）、車間制造執行系統（manufacturing execution system，ΜES）等，可實現對中藥前提取設備傳感器實時數據以及生產過程數據的抓取。在工業數據的挖掘與分析中，機器學習是常用的方法之一[4]。杜慧等[5]以金銀花和青蒿濃縮浸膏醇沉工序為研究對象，采用偏最小二乘法（partial least squares，PLS）進行數據挖掘，PLS 需要對預測變量和目標變量進行降維投影處理[6]，不利于直觀解釋生產規律，該研究篩選了關鍵變量但未劃分變量的控制范圍，不便于生產的實際調控。決策樹算法是機器學習中一個樹狀預測模型，適用于多階段決策、表查找和最優化等問題，能夠直觀解釋決策過程，目前已成功應用于醫療、廣告、交通和金融等多個領域[7-8]。本實驗采用決策樹及其相關算法對歷史數據進行分析，深入挖掘并尋找關鍵變量的優化控制范圍，期望提升熱毒寧注射液的質量控制水平。

1 材料

建模數據由江蘇康緣藥業股份有限公司數字化提取工廠提供，金青中間體批次為170105 至180905共205 批；數據的訓練集和測試集的劃分應用Python 的Scikit-Learn 工具實現。建模算法均在Salford Predictive Μodeler 8.3（美國Μinitab 公司）實現，單值過程控制圖在Μinitab 19（美國Μinitab公司）繪制。

2 方法

2.1 目標變量的確定

金青醇沉工段由金青浸膏分配、金青浸膏醇沉、上清液暫存和上清液回收4 個過程組成。金銀花和青蒿的水提濃縮浸膏通過稱量罐稱定質量后分配至2 個醇沉罐（醇沉罐1 和醇沉罐2），在醇沉過程中，向分配浸膏的2 個醇沉罐中分3 個階段按不同的體積流量加入一定量的乙醇并不斷攪拌，攪拌停止靜置一段時間并將醇沉上清液通過真空泵轉入濃縮罐，刮板濃縮后得金青醇沉浸膏。金青醇沉工段簡圖如圖1所示，各過程主要參數如表1所示。前期辨識出金青醇沉工段的關鍵因素為金青醇沉上清液和上清液傳料體積、金青總質量、醇沉罐1 及醇沉罐2 金銀花浸膏質量[5]。上清液傳料體積為上清液暫存罐傳送至濃縮罐的上清液體積；醇沉上清液為上清液暫存罐收集的上清液體積；金青浸膏總質量為2 個醇提罐加入金銀花和青蒿總質量；醇沉罐1和醇沉罐2 金銀花質量為分配至2 個醇沉罐的金銀花質量。所篩選的參數均為金青浸膏醇沉過程輸入物料及輸出物料的質量和體積參數，醇沉過程工藝較復雜，涉及影響因素較多，為進一步挖掘潛在規律，本研究選擇金青浸膏醇沉過程為研究目標建立過程模型。

表1 金青醇沉4 個過程的主要參數Table 1 Main parameters of four processes of LA alcohol precipitation

圖1 金青醇沉回收工段流程簡化圖Fig.1 Data map for LA alcohol precipitation process

為消除過程物料質量體積相關性對模型的影響，本研究擬以浸膏醇沉過程收率為目標變量。將金青浸膏醇沉收率和金青醇沉工段收率進行相關性分析，相關圖如圖2所示，相關系數（r）為0.653。控制浸膏醇沉過程對提高金青醇沉工段的生產效率十分重要。綜上，本實驗以金青浸膏醇沉收率作為目標變量，建立浸膏醇沉過程模型。

圖2 金青浸膏醇沉收率與金青醇沉工段收率相關圖Fig.2 Correlation diagram between yield of LA extract alcohol precipitation and yield of LA alcohol precipitation section

2.2 建模數據準備

本實驗建模數據主要來源于批記錄及數字化提取工廠數據庫集成數據。批記錄為手動錄入，主要為原料藥、中間體物料相關數據；數據庫收集的數據主要從DCS 系統中收集各個設備傳感器實時生產過程中抓取的時序數據，抓取頻率為每分鐘1 次。中藥生產環境復雜，不同時間點內時序參數相關性強，前期對DCS 數據進行數據清洗和特征提取后得到建模特征參數[6]。本實驗在建模特征參數中選擇30 個醇沉罐傳感器數據參數納入建模數據集，考慮到浸膏質量和加醇體積間存在較強的相關性，本實驗將金銀花浸膏質量、青蒿浸膏質量和加醇量轉化為金青浸膏質量比和加醇料液比。在中藥制劑的過程中，中藥材的質量差異會傳遞至處方藥味、中間體及成品，直接影響中藥制劑批間質量的穩定[9]。為考察藥材及金銀花浸膏及青蒿浸膏的各生產階段可能對醇沉過程的影響，收集藥材質量參數、藥材廠商及藥材至濃縮浸膏各過程工序的收率作為預測變量，樣本收集時間為2017年1月—2018年9月，共205 批歷史數據。以批號關聯建模參數，建模參數共49 個，最終得到205×49 的建模數據集，建模變量如表2所示。

表2 金青浸膏醇沉過程49 建模變量Table 2 Fourty-nine modeling variables in LA extract alcohol precipitation process

2.3 建模算法

決策樹（decision tree）是一種基本的分類和回歸方法，該學習方法主要包括特征選擇、決策樹生成和修剪，目前決策樹學習常用的算法有ID3、C4.5和分類與回歸樹（classification and regression tree，CART）[10]。當單一決策樹學習器不能很好地處理問題時，可嘗試集成學習方法（ensemble learning）提高模型預測性能。集成學習通過訓練多個學習器并進行組合來解決同一問題[11]。Bagging 算法和Boosting 算法是2 種經典的集成學習算法[12]。本實驗采用CART 算法以及隨機森林、TreeNet 2 種集成樹算法進行建模。

2.3.1 CART CART 模型是應用廣泛的決策樹算法[13]，它是由特征選擇、樹的生成及剪枝組成，既可以用于分類，也可以用于回歸[10]。相較于ID3，CART 算法選擇基尼不純度作為度量標準，并采用二分遞歸分割方法劃分數據集，最終生成二叉決策樹，二叉劃分適用更加廣泛，也更易處理數值型數據[6]。

2.3.2 隨機森林（random forests，RF） RF 是由多個Bagging 算法訓練的決策樹組合而成[14]，在解決回歸問題時，其以多個決策樹的平均值為輸出結果[12,15]。它是通過利用bootstrap 重抽樣方法從原始樣本中抽取多個樣本，沒有抽到的樣本被稱為帶外數據（out-of-bag，OOB）[16]。OOB 誤差估計是RF算法的一個無偏估計，可替代數據集的交叉驗證，明顯降低計算復雜性[12,17]。

2.3.3 TreeNet TreeNet 被稱為隨機梯度提升（stochastic gradient boosting）[18]。梯度提升是一種實現Boosting 的方法，其主要思想是每次建模都在前一模型損失函數（loss function）梯度下降方向進行[19]。TreeNet 在訓練的過程中首先生成1 棵小樹并計算該模型的預測殘差，之后建立第2 棵樹對第1 棵樹的殘差進行預測，第3 棵樹對前2 棵樹的殘差進行預測。上步驟遞歸，新產生的小樹對前面的模型進行不斷修正，最終得到最佳模型。

2.4 數據集劃分與訓練

數據按3∶1 的比例隨機劃分為訓練集和測試集。通過訓練集的數據來訓練數據，用測試集來評估模型。在訓練模型的過程中需要訓練不同模型參數的模型來篩選最佳模型，將保留訓練集的一部分來評估各個模型的性能，新的保留集稱為驗證集。本實驗采用10 折交叉驗證來評估CART 和TreeNet的訓練模型，RF 為OOB 誤差估計。本實驗以均方根誤差（root mean square error，RΜSE）、平均絕對偏差（mean absolute difference，ΜAD）以及決定系數（R2）作為模型評價指標。RΜSE 及ΜAD 越接近0，R2越接近1 表示模型誤差越小，性能越好。

n為校正集或驗證集的樣本數，t∈[1，n]，y為參考值，yi為預測值，為所有樣品參考值的平均值

3 結果與分析

3.1 各模型性能評價

基于訓練集數據分別建立CART、RF和TreeNet模型，并采用10 折交叉驗證來評估不同超參數的CART 和TreeNet 訓練模型性能，OOB 誤差估計評估RF 模型性能。最終得到各模型最佳的超參數設置如表3所示，最優模型的性能評價指標如表4所示，各模型的預測變量重要性如表5所示。由表5中可知，與CART 算法相比，RF 和TreeNet 2 種集成樹算法的驗證集和預測集均方根更小，決定系數更大，模型效果更好。說明集成學習算法有效地提高了單一決策樹算法模型的性能。

表3 3 種模型的超參數設置Table 3 Hyperparameter setting table of three models

表4 3 種最優模型的性能指標評價Table 4 Performance index evaluation of three models

由表5可知，3 種模型的共同重要變量為X16（罐1 料液比）、X34（罐2 料液比）和X5（金銀花提取液濃縮收率）。以這3 個變量建立模型，其模型效果如表6所示。

表5 CART、RF 和TreeNet 模型變量重要性Table 5 Importance chart of predicted table of CART,RF and TreeNet

表6 篩選后3 種模型的性能指標評價Table 6 Performance index evaluation of three models after screening

由變量篩選后模型性能觀察比較可知，變量篩選后測試集和預測集的模型誤差及R2變化較小，均保持較好的模型性能。說明醇沉罐1 料液比、醇沉罐2 料液比以及金銀花提取液濃縮收率為金青浸膏醇沉過程的關鍵變量。

3.2 關鍵變量依存性分析

為進一步描述3 個關鍵變量與目標變量的潛在規律。對TreeNet 模型的3 個關鍵變量的變量依存度圖（partial dependency plots，PDPs）進行分析，單變量依存性散點圖如圖3所示。

圖3 關鍵變量為金銀花提取液濃縮收率 (A)、罐1 料液比 (B) 及罐2 料液比 (C) 單變量依存度散點圖Fig.3 PDPs of univariate dependency with key variable were concentration yield of honeysuckle extract liquid (A) and solidliquid ratios in extraction tank 1 (B) and that in tank 2 (C)

通過對3 個關鍵變量的PDPs 分析，各變量與其部分依存度不是簡單的線性關系。變量依存度分析以響應值的正負表示對轉化率影響的好壞，將工藝優化問題轉化為計算問題。當金銀花提取液濃縮收率小于0.039 kg/L 時，目標變量預測響應值為正值。醇沉罐1 和醇沉罐2 料液比的依存度隨料液比的不斷增加的變化可分為3 個范圍：料液比較低，當料液比較小時依存度為負值，較小的料液比將降低目標變量的預測響應值；料液比居中，隨著料液比的不斷增加依存度也逐漸提高并轉為正值；料液比較高，當增加到一定范圍時，料液比的增加對目標變量預測響應值的影響不變。將料液比控制在合適區間，在提高金青浸膏醇沉的收率同時可避免加入過量乙醇，降低生產成本。

模型結果中醇沉罐1 的料液比控制范圍大約為5.54～5.91 L/kg，醇沉罐2 的范圍大約為5.58～5.91 L/kg。對2 罐的范圍進行中位數檢驗，分析2 罐間是否存在差異，以便于確定是否需要分罐控制。結果如表7所示。

表7 不同提取罐料液比中位數Table 7 Median of solid-liquid ratios in different extraction tanks

醇沉罐1 和醇沉罐2 相比，P值＞0.05，總體中位數相等，醇沉罐1 和醇沉罐2 間料液比優化范圍無顯著差異。綜上，醇沉罐1 和醇沉罐2 可統一控制，為控制金青浸膏醇沉過程的收率，需主要對金銀花提取液濃縮收率和各罐料液比進行控制。

3.3 關鍵變量過程控制

上文篩選影響金青浸膏醇沉的關鍵變量為金銀花提取液濃縮收率和料液比，金銀花醇沉上清液濃縮收率的影響因素較多，涉及上一工段各工藝參數及設備過程參數，需另外建模研究。本實驗主要對料液比的過程控制進行研究。

控制料液比主要控制金銀花浸膏、青蒿浸膏和加醇量。通過收集落入優化范圍內外的金銀花浸膏、青蒿浸膏和加醇量，進行中位數檢驗分析差異。結果如表8所示。統計分析表明，影響料液比的變化主要是由于金銀花浸膏質量和加醇量的波動引起的。企業內部控制標準中要求金銀花浸膏濃縮至1.10～1.12 g/cm3（70～80 ℃），并分別對密度為1.10、1.11、1.12 g/cm3的浸膏的出膏量給出不同的控制范圍。依存度分析優化范圍內的金銀花浸膏出膏溫度均符合要求，浸膏密度分為1.11、1.12 g/cm3。分別對密度為1.11、1.12 g/cm3的金銀花浸膏做控制圖，控制圖如圖4所示。

圖4 密度為1.11 g·cm-3 金銀花浸膏質量 (A) 與加醇量 (B) 及密度為1.12 g·cm-3 金銀花浸膏質量 (C) 與加醇量 (D) 的過程控制圖Fig.4 Process control diagram of quality of honeysuckle extract with a density of 1.11 g·cm-3 (A),amount of alcohol (B),quality of honeysuckle extract with a density of 1.12 g·cm-3 (C),and amount of alcohol (D)

表8 不同料液比的相關變量中位數Table 8 Median of related variables for different solidliquid ratios

如控制圖所示，通過篩選得到密度為1.11 的金銀花浸膏分配控制范圍為557.92～639.62 kg，加醇量的控制范圍為3.370～3.828 m3；密度為1.12 的金銀花浸膏的控制范圍為540.4～616.9 kg，加醇量的控制范圍為3.317～3.859 m3。綜上分析，優化范圍內加醇量的差異不大，金銀花浸膏質量的波動是影響金青醇沉過程的重要因素。需加強對金銀花浸膏質量的監控手段和方法，將金銀花浸膏質量控制在穩定的優化范圍內控制加醇量，從而控制醇沉過程高效穩定進行。

4 總結與展望

醇沉工序是中藥提取精制過程的關鍵工序之一，在實際生產控制過程中對金青醇沉過程控制的穩定性還有待提高。本實驗依托熱毒寧提取精制生產大數據，收集金青醇沉過程歷史數據，以金青浸膏醇沉上清液收率為目標變量，基于決策樹算法進行建模分析，綜合3 個模型建模結果，篩選出重要變量為金銀花上清液濃縮收率和醇沉罐料液比；通過變量依存度分析，優選批次并作金銀花質量及加醇量過程控制圖優化控制范圍。從熱毒寧注射液工業數據中挖掘質量傳遞規律，有利于該產品質量的追溯與控制，促進制造過程的優化及持續改進[17]。

決策樹算法在中藥工業大數據的建模應用報道較少，本實驗首次嘗試應用于熱毒寧注射液的醇沉過程，相較于PLS 算法，決策樹算法在缺失值處理、變量篩選及優化參數范圍的劃分等方面更具優勢，為熱毒寧注射液的數據挖掘和質量控制提供新路徑。本實驗建模數據為生產歷史數據，在后續的工作中，還需結合近紅外快速檢測技術收集中間體含量數據，多目標決策優選工藝參數控制范圍。此外，還需持續收集生產數據，不斷優化及驗證工藝參數控制范圍，提升金青醇沉過程的控制水平。

利益沖突所有作者均聲明不存在利益沖突