偏最小二乘法和THz-TDS在正品大黄鉴别中的应用

汪景荣，张卓勇*，张振伟，相玉红

1. 首都师范大学化学系，北京 100048 2. 首都师范大学物理系，太赫兹光电子学教育部重点实验室，北京 100048

偏最小二乘法和THz-TDS在正品大黄鉴别中的应用

汪景荣1，张卓勇1*，张振伟2，相玉红1

1. 首都师范大学化学系，北京 100048 2. 首都师范大学物理系，太赫兹光电子学教育部重点实验室，北京 100048

太赫兹技术的发展近年来受到广泛的关注并被应用于热点。中草药大黄的品质鉴定对于中药制剂的质量控制具有重要的意义。利用大黄的太赫兹时域光谱结合偏最小二乘法(PLS)模型对基于41个正品和非正品大黄的中草药鉴别模型进行了研究。首先采集大黄样品的太赫兹时域光谱(THz-TDS)信号，然后将化学计量学方法用于这些大黄样品太赫兹光谱的信号处理与建模，再建立基于太赫兹光谱的大黄品质鉴定的偏最小二乘模型方法。应用S-G一阶导数、去趋势、标准正态变换、自标度化、均值中心化等方法对原始时域谱预处理再与未经预处理的结果相比，偏最小二乘(PLS)模型的预测正确率从80%明显提高到90%。在模型建立和模型检验中，采用留一法(LOO)选取训练集和检验集样本。利用留一法交叉验证确定了PLS模型的最佳主因子数。结果表明，当采用均值中心化方法时，PLS模型的RMSECV和RMSEP的值均达到了最小，分别为0.076 6和0.169 0。研究结果表明，THz-TDS技术结合化学计量学方法能够快速、准确的对大黄的真伪进行鉴别，直接使用太赫兹时域光谱而不使用计算后的吸收谱有两个优点： (1)在分频测定和光谱信号处理时无需考虑样品的厚度；(2)使光谱信号处理过程得到简化。该技术也可以对其他中草药进行鉴别和质量控制。该法快速、简单、无污染、无需样品预处理，是一种有发展前景的中草药无损检测方法。

太赫兹时域光谱； 偏最小二乘法； 大黄

引 言

大黄是常用的传统中草药，为多年生草本植物，用药历史悠久。掌叶大黄RheumpalmatumL.、唐古特大黄RheumtanguticumMaxim.exBalf.或药用大黄RheumofficinaleBaill.的干燥根及根茎，为中国药典收录的三种正品大黄[1]。近年来，随着大黄用药量的不断增加，正品大黄已不能满足需求，一些伪品大黄如华北大黄RheumfranzenbachiiMunt.、河套大黄RheumhotaoenseC.Y.ChengetKao、藏边大黄RheumemodiiWall.等冒充正品大黄出现在市场上，严重影响了大黄的临床效用[2]。正品大黄一般具有利胆退黄、保肝、止血、抗菌、泻下等作用[3]，而伪品大黄的泻下作用较弱。因此，为了控制大黄药材的质量，保证大黄用药的安全性和有效性，需对大黄进行准确鉴定。随着中药指纹图谱技术的发展，各种色谱、光谱及其他仪器分析技术开始广泛用于大黄的鉴定。杨美华等[4]应用RAPD方法对正品和伪品大黄进行指纹图谱的研究，该法达到了在分子水平上鉴定正品和伪品大黄的目的。冯有龙等[5]利用HPLC对大黄及部分含大黄中成药的真伪进行了鉴别。李磊等[6]采用高效液相色谱法建立大黄药材的指纹图谱，该法可以明显区分正品与非正品大黄。但是这些指纹谱分析方法都需对样品进行分离、浸取、浓缩等预处理，属于有损分析，操作复杂且耗费时间，测试结果在一定程度上依赖于预处理技术。因此，发展一种无须预处理，就能快速、准确对大黄进行鉴别和质量控制的方法，有着重要的意义。随着数学、统计学、计算机等技术的发展，化学计量学结合近红外光谱技术开始广泛用于大黄的鉴别[7-8]。

太赫兹辐射，又称T射线，是指频率在0.1～10 THz范围内的电磁波，位于微波和红外之间，是微观电子学向宏观光子学过渡的区域，与许多材料具有独特的相互作用[9]。随着太赫兹辐射技术的不断发展，太赫兹光谱开始广泛用于中草药领域[10-11]。但是当待测样品在太赫兹波段没有明显的特征吸收峰时，这些方法将无法对样品进行准确鉴别。随着计算机技术的不断发展，化学计量学开始用于中草药的太赫兹光谱分析[12]。但是这些方法都需要对测得的原始时域谱进行傅里叶变换，从而得到样品的吸收谱和折射率谱，计算方法复杂，信息处理过程耗时。因此，发展一种快速、准确的中草药无损鉴别方法十分必要。

偏最小二乘法(PLS)[13-14]是一种多元统计分析方法，主要通过对光谱矩阵和性质矩阵同时进行分解，在分解时考虑两者之间的相互关系，从而建立最佳的校正模型。Kennard-Stone(KS)法[15]作为一种校正样本的选择方法，其主要优点是所选择的的校正集样本空间距离分布均匀，具有较好的代表性。

本研究发展了一种直接应用大黄的太赫兹时域谱结合偏最小二乘法来对大黄进行快速、准确鉴别的新方法。该法与近红外方法相比，具有以下优点： (1)大多数中草药中所含生物分子的低频振动或转动模式都处在THz波段，而近红外光谱不能达到这一区域。(2)THz辐射脉冲宽度在皮秒到飞秒量级，可以对中药进行时间分辨的瞬态光谱研究，而且通过采样测量技术，能有效抑制背景辐射的干扰，得到信噪比很高的时域谱。(3)THz-TDS技术可以对许多大分子进行指纹识别，从而检测物质结构的微小差异[16]。而近红外光谱主要是基于分子内部振动的倍频和合频吸收，谱峰较宽，重叠严重且吸收强度较弱[17]。

实验以41个大黄样品为研究对象，应用大黄的太赫兹时域谱建立PLS定性分析模型。应用留一法交互验证[18]选取PLS模型的最佳主因子数，同时比较了不同的预处理方法对模型预测能力的影响。

1 背景理论

1.1 PLS

采用NIPALS算法的PLS模型有效地克服了主成分分析的不足，在对光谱矩阵X进行主成分分解的同时也对性质矩阵Y进行主成分分解，且通过迭代时交换X与Y的迭代矢量使两个分解过程合二为一[19]。其主要原理[20-22]如下。

对光谱矩阵X和性质矩阵Y进行线性分解，如式(1)和式(2)所示。

X=TPT+E

(1)

Y=UQT+F

(2)

其中，X是一个矩阵，n是光谱的行数，p是光谱采样点数。Y是一个由1和0构成的性质矩阵，1代表正品大黄，0代表非正品大黄。T和U分别是X和Y的得分矩阵，P和Q是X和Y的载荷矩阵，E和F分别是X和Y的残差矩阵。

对T和U进行线性回归

U=TB

(3)

B=(TTT)-1TTY

(4)

由P得到预测集样品光谱矩阵Xpred的得分矩阵Tpred，然后通过式(5)得到Y的预测值。

Ypred=TpredBQ

(5)

在应用PLS确定主因子数时，使计算所得主因子数之间方差尽可能最大的同时，还使主因子数与性质最大程度相关。

1.2 KS法

KS法基于变量之间的欧式距离，从所有的样本中依次挑选合适的样本进入校正集[23-24]。首先将欧式距离最远的两个向量对选入校正集，随后在迭代过程中将拥有最大最小距离的样本选入校正集，以此类推，直到满足要求的样本数目[25]。该法能够保证校正集中的样本空间距离分布均匀，且具有代表性。

2 实验部分

2.1 仪器及参数

实验使用Z-3透射式THz时域光谱系统，如图1所示。系统的发射极为InAs晶体，探测极为ZnTe晶体。实验中使用的激光器为Spectra-Physics公司的锁模钛宝石激光器。脉冲的中心波长为800 nm，频率范围为0.1～3 THz，积分时间为100 ms，测量速度为0.002 5 mm·s-1，步长为0.1 mm。使用MVDS-400 controller软件来记录测量数据。

Fig.1 Schematic diagram of THz-TDS measurement

2.2 样品制备

用于实验的大黄样品共有41个，其中17个为正品大黄，24个为非正品大黄，均由北京同仁堂股份公司科学研究所提供。所有的大黄样品经干燥后粉碎，过60目筛，再次真空干燥2 h，留待制样。实验时，称量150 mg样品粉末于玛瑙研钵中研磨，在6.5 t压力下压成厚度在0.9～1.2 mm之间，直径为13 mm，两表面平行且没有裂纹的圆形薄片。

2.3 样品太赫兹光谱的采集

室温下，在系统内通入氮气保持相对湿度<1%，防止水汽对实验的影响。将样品从1～41编号，测量时随机取出一个置于太赫兹脉冲的中心位置，探测样品的时域波形。移动样品，在其不同位置测量3次，以保证样品具有代表性。每测完一个样品，同时测量一次以氮气为背景的参考信号。最后将所测得的样品信号和参考信号取平均，即可得到41个大黄样品的太赫兹时域谱。每条时域光谱有2 985个数据点。

2.4 模型的评价参数

应用PLS法建立大黄的太赫兹时域谱和性质矩阵的定性分析模型。模型分别以校正集样品的交互验证均方根误差(root-mean squared error of cross-validation, RMSECV)和预测集样品的预测均方根误差(root-mean squared error of prediction, RMSEP)作为评价指标。计算如式(6)和式(7)所示。

(6)

(7)

其中，yt, p为模型的预测值，yi, a为模型的实际值。在本文中，正品大黄的目标输出为1，非正品大黄的目标输出为0。为了使模型的预测精度更高，实验中阈值设为0.2，即当模型的预测值小于0.2时，判定为非正品大黄； 当模型的预测值大于0.8时，判定为正品大黄。n为校正集样品数，m为预测集样品数。

2.5 数据处理和建模

应用KS法选取31个大黄样本为校正集，其余10个样本为预测集。应用多元散射校正(multiplicative scatter correction，MSC)、Savitzky-Golay五点三次多项式平滑(S-G平滑)、Savitzky-Golay五点三次一阶导数(S-G一阶导数)、标准正态变换(standard normal variation，SNV)、去趋势(Detrending)、自标度化(Autoscaling)、均值中心化(Mean centering)等方法对原始时域谱进行预处理，建立PLS分类模型，整个模型的校正主因子数通过留一法交互验证选取，应用预测均方根误差和正确率来评价模型对未知样本的预测能力。同时比较不同的预处理方法对模型预测正确率的影响。

3 结果与讨论

3.1 41个大黄样品的太赫兹时域谱

41个大黄样品的太赫兹时域光谱如图2所示。两类大黄样品信号的振幅强度明显低于参考信号，可能因样品对太赫兹辐射的吸收和散射导致。样品信号相对于参考信号有一定的时间延迟，这可能是由于不同样品相对于氮气有不同的折射率。从图2可以看出，41个大黄的时域谱重叠在一起，仅凭肉眼难以区分出正品和非正品大黄，因此借助化学计量学方法进行鉴别。

Fig.2 Time-domain spectra of rhubarb

3.2 KS法划分校正集和预测集

正品和非正品大黄化学组成复杂，样本信息包含在相似性很强的时域光谱中。因此，在PLS建模过程中，选择代表性强的样本，不仅可以提高建模的速度，而且可以提高模型的预测精度和扩大模型的适用范围。应用KS法选取3/4的样本为校正集，1/4的样本为预测集，所选的校正集和预测集的样品编号如表1所示。

Table 1 Calibration set and!prediction set were selected by KS method

为了更清晰地反应所选校正集和预测集中样品的分布情况，我们对41个样品光谱做主成分分析，其得分图如图3所示。

Fig.3 Principal component scores of the calibration set and prediction set selected by KS method

31个校正集样品在主成分得分图上分布比较均匀，10个预测集样品被包含在校正集中。可见，预测集的光谱信息被包含在校正集的光谱信息之中，KS法可以提高模型的预测精度[25]。

3.3 PLS模型主因子数的选择

在应用PLS建立校正模型时，参与回归的最佳主因子数选取对模型的预测精度有重要影响。如果选取的主因子数太多，会将测量过程中产生的噪声拟合到模型中，造成过拟合。如果选取的主因子数太少，会损失原始光谱矩阵中的有用信息，降低模型的预测精度，造成欠拟合。因此，合理选择建模的最佳主因子数能够有效利用光谱信息和滤除光谱噪声[26]。本实验不用对原始时域谱进行预处理，应用留一法对校正集进行交互验证，根据RMSECV值来选取建模的最佳主因子数。选取1～10个主因子数，依次计算其RMSECV值，当主因子数为3时，其RMSECV值最小。此时，计算模型对未知样本的预测正确率为80%。模型的预测能力还不太理想。

3.4 不同的预处理方法对PLS模型的影响

由于在光谱测量时，原始时域谱中除包含样品化学信息外，还可能存在其他无关信息和噪声，如样品背景和杂散光等[17]。因此，在建模时，应用一定的预处理方法，可以滤除原始光谱中的无关信息和噪声，从而提高模型的预测能力和稳健性。分别应用MSC、S-G平滑、S-G一阶导数、SNV、detrending、autoscaling、mean centering等方法对原始时域谱进行预处理，应用留一法选取最佳的校正主因子数，随后建立PLS定性模型对预测集样本进行预测，所得结果如表2所示。

Table 2 Effects of different pretreatment methods on the PLS model

原始时域谱经S-G平滑处理后和未经预处理的建模主因子数、RMSECV、RMSEP值都相等，且对未知样品的预测正确率均为80%，没有显著的提高。这可能是由于应用S-G平滑对原始时域谱处理时，在去除光谱噪声的同时也去除了原始光谱中与性质矩阵相关的信息。应用MSC对原始时域谱进行预处理后，模型的预测正确率反而下降到60%。这可能与原始光谱数据和MSC算法两方面有关系。MSC算法假定散射与波长及样品的浓度变化无关[17]，而对于大黄样本的光谱数据，真伪大黄的化学成分复杂，组分性质变化较大，应用MSC来校正达不到预期的效果，反而可能会引入过多的噪声。

应用S-G一阶导数、去趋势、标准正态变换、自标度化、均值中心化等方法对原始时域谱预处理后，PLS模型的预测正确率从80%显著提高，达90%。S-G一阶导数和标准正态变换两种预处理方法，建模的主因子均为3。与未经预处理的PLS模型相比，主因子数没有发生变化，也就没有达到使模型更加简洁的目的。同时，RMSEP值却有明显的增大。而经去趋势预处理后的建模主因子数明显减少，由3降到1，达到了简化模型的目的，但同时RMSEP值也增加了。与S-G一阶导数、标准正态变换预处理方法相比，自标度化和均值中心化的建模主因子数增加了一倍，在保持正确率为90%的同时，RMSECV和RMSEP值却明显下降。其中，均值中

心化的RMSECV和RMSEP值最小，达到了最理想的建模效果。这可能是均值中心化方法将光谱的变化与大黄样品组成的变化关联起来，使样品光谱之间的差异增大，从而使模型的预测能力得到有效提高。应用均值中心化对原始时域谱进行预处理后，根据留一法选取的PLS建模主因子数和RMSECV值的关系如图4所示。

Fig.4 RMSECV obtained by preprocessing method of mean centering for the time domain spectra with respect to PLS component number

由图4可知，大黄的原始时域谱经均值中心化预处理后，RMSECV值随建模主因子数的变化呈现出理想的变化趋势。RMSECV值最小时主因子数为6，此为建立PLS模型的最佳主因子数。

4 结 论

应用THz-TDS系统在室温下测量了41个大黄样品的THz时域谱，所测太赫兹时域光谱直接用于光谱预处理和建模。应用不同的预处理方法结合PLS建立了大黄的定性分析模型。整个模型的最佳主因子数通过留一法交互验证来选取。实验结果表明，原始时域谱未做预处理时，建模主因子数为3，预测正确率为80%。应用S-G一阶导数、去趋势、标准正态变换、自标度化、均值中心化等方法对原始时域谱预处理后，预测正确率从80%明显提高到90%。其中，应用均值中心化预处理后，与其他方法所得结果相比，建模的主因子数虽然达到了6，但是RMSECV和RMSEP的值均达到了最小，分别为0.0766和0.1690，使模型的预测能力达到最好。该法简便、无污染，无需对样品预处理，就能准确鉴别正品与非正品大黄，还可以用于其他中草药生产中的质量控制。

[1] Chinese Pharmacopoeia Commission(国家药典委员会). Chinese Pharmacopoeia(2010 Edition) (中华人民共和国药典, 2010年版). Beijing: China Medical Science and Technology Press(北京： 中国医药科技出版社), 2010.

[2] LU Jing(卢 静). Strait Pharmaceutical Journal (海峡药学), 2008, 20(1): 63.

[3] XU Xiang, LI Bai-ping, ZHANG Hui-fen (徐 翔, 郦柏平, 张慧芬). Shanghai Journal of Traditional Chinese Medicine(上海中医药杂志), 2003, 37(4): 56.

[4] YANG Mei-hua, ZHANG Da-ming, LIU Jian-quan, et al(杨美华, 张大明, 刘健全, 等). Chinese Traditional and Herbal Drugs (中草药), 2003, 34(6): 557.

[5] FENG You-long, YU Bo-yang(冯有龙, 余伯阳). Drug Standards of China(中国药品标准), 2009, 10(4): 296.

[6] LI Lei, LIU Rui, YUAN Bo, et al(李 磊, 刘 瑞, 袁 波, 等). Chin. Pharm. J.(中国药学杂志), 2005, 9(17): 1302.

[7] ZHANG Xiao-fang, ZHANG Zhuo-yong, FAN Guo-qiang(张孝芳, 张卓勇, 范国强). Chin. J. Pharm. Anal.(药物分析杂志), 2006, 26(7): 914.

[8] Wang F X, Zhang Z Y, Cui X J, et al. Talanta, 2006, 70: 1170.

[9] Rainsford T J, Mickan S P, Abbott D. T-Ray Sensing Applications: Review of Global Developments. Proceedings of the International Society for Optical Engineering, 2005, 5649: 827.

[10] ZHANG Jian, HUANG Wan-xia, LUO Yi, et al(张 建, 黄婉霞, 罗 轶, 等). China Pharmacy(中国药房), 2011, 22(47): 4467.

[11] Zhao X L, Li J S. J. Phys.: Conf. Ser., 2011, 276: 1.

[12] Li H, Du S Q, Xie L, et al. Optik, 2012, 123(13): 1129.

[13] CHU Xiao-li(褚小立). Molecular Spectroscopy Analytical Technology Combined with Chemometrics and Its Applications(化学计量学方法与分子光谱分析技术). Beijing: Chemical Industry Press(北京： 化学工业出版社), 2011.

[14] Höskuldsson A. Journal of Chemometrics, 1988, 2(3): 211.

[15] Kennard R W, Stone L A. Technometrics, 1969, 11(1): 137.

[16] WANG Yi-fan, WEI Wan-cong, ZHOU Feng-juan, et al(汪一帆, 尉万聪, 周凤娟, 等). Progress in Biochemistry and Biophysics(生物化学与生物物理进展), 2010, 37(5): 484.

[17] CHU Xiao-li, YUAN Hong-fu, LU Wan-zhen(褚小立, 袁洪福, 陆婉珍). Progress in Chemistry(化学进展), 2004, 16(4): 528.

[18] Cawley G C, Talbot N L C. Neural Networks, 2004, 17(10): 1467.

[19] NI Li-jun, ZHANG Li-guo(倪力军, 张立国). Basis of Chemometrics and Its Applications(基础化学计量学及其应用). Shanghai: East China University of Science and Technology Press(上海： 华东理工大学出版社), 2011.

[20] Helland I S. Scandinavian Journal of Statistics, 1990, 17(2): 97.

[21] Wold S, Sjöström M, Eriksson L. Chemom. Intell. Lab. Syst., 2001，58(2): 109.

[22] Helland I S. Communications in Statistics-Simulation and Computation, 1988, 17(2): 581.

[23] Kocjancic R, Zupan J. Chemom. Intell. Lab. Syst., 2000, 54(1): 21.

[24] Siano G G, Goicoechea H C. Chemom. Intell. Lab. Syst., 2007, 88(2): 204.

[25] ZHU Xiang-rong, LI Na, SHI Xin-yuan, et al(朱向荣, 李 娜, 史新元, 等). Chemical Journal of Chinese Universities(高等学校化学学报), 2008, 29(5): 906.

[26] XIE Jun, PAN Tao, CHEN Jie-mei, et al(谢 军, 潘 涛, 陈洁梅, 等). Chinese Journal of Analytical Chemistry(分析化学), 2010, 38(3): 342.

*Corresponding author

Identification of Official Rhubarb Samples by Using PLS and Terahertz Time-Domain Spectroscopy

WANG Jing-rong1, ZHANG Zhuo-yong1*, ZHANG Zhen-wei2, XIANG Yu-hong1

1. Department of Chemistry, Capital Normal University, Beijing 100048, China 2. Department of Physics, Capital Normal University, Key Laboratory of Terahertz Optoelectronics, Ministry of Education, Beijing 100048, China

The development of terahertz technology is attracting broad intention in recent years. The quality identification is important for the quality control of Chinese medicine production. In the present work, terahertz time-domain spectroscopy (THz-TDS) combined with partial least squares (PLS) were used for the identification model building and studied based on 41 official and unofficial rhubarb samples. First, the THz-TDS spectra of rhubarb samples were collected and were preprocessed by using chemometrics methods rather than transformed to absorption spectra. The identification models were then established based on the processed terahertz time domain spectra. The spectral preprocessing methods include Savitzky-Golay (S-G) first derivative, detrending, standard normal transformation (SNV), autoscaling, and mean centering. The identification accuracy of 90% was accomplished by using proper pretreatment methods, which was higher than the classified accuracy of 80% without any preprocessing for the time domain spectra. The component number of the PLS model was evaluated by leave-one-out cross-validation (LOOCV). The minimum values of the root-mean squared error of cross-validation (RMSECV) and root-mean squared error of prediction (RMSEP) were 0.076 6 and 0.169 0 by using mean centering method, respectively. The results of this work showed that the combination of terahertz time domain spectroscopy technology with chemometrics methods, as well as PLS can be applied for the recognition of genuine and counterfeit Chinese herbal medicines, as well as official and unofficial rhubarbs. The advantage of using terahertz time domain spectra directly with no transformation into absorption spectra is: (1) the thickness of samples could not be considered in the model establishment, and (2) the spectral processing was simplified. The proposed method based on the combination of THz-TDS and chemometrics proved to be rapid, simple, non-pollution and solvent free, which is suitable to be developed as a promising tool for quality control of many other Chinese herbal medicines.

Terahertz time-domain spectroscopy; Partial least squares; Rhubarb

Aug. 11, 2014; accepted Dec. 20, 2014)

2014-08-11，

2014-12-20

国家重大科学仪器专项项目(2012YQ140005)和国家自然科学基金项目(21275101；11204191)资助

汪景荣，1986年生，首都师范大学化学系硕士研究生 e-mail: jingrong626@163.com *通讯联系人 e-mail: gusto2008@vip.sina.com

O657.3; O433.1

A

10.3964/j.issn.1000-0593(2016)02-0316-06