CN101957311A - 一种素片干燥过程的定性分析方法 - Google Patents

Abstract

本发明属于药物制剂领域，具体涉及一种素片干燥过程的定性分析方法。本发明提供的方法通过结合近红外光谱、主成分分析法和支持向量机法等方法快速建立了素片含水量定性模型，实现了对素片含水量合格与否的快速评价，并且监测准确率基本达到了工业生产的要求。相对于传统方法或经验判断，本发明提供的定性分析方法能够实现在素片干燥过程中快速、无损、可靠、准确地在线监测和/或反馈其含水量是否合格的目的。

Description

一种素片干燥过程的定性分析方法

技术领域

本发明属于药物制剂领域，具体涉及一种素片干燥过程的定性分析方法。

背景技术

直接压制而成的片剂，若其药物性质不稳定、有不良气味、对胃有刺激，或易被胃液破坏等，通常还需要在其表面上包一层辅料，使片内药物与外界隔离，从而进一步保证片剂质量。在药物制剂领域中，通常将未覆盖包衣层的片剂称为素片。

素片的含水量是影响片剂质量的重要因素之一。对于糖衣片，如果素片含水量过高，则水分容易渗入到粉衣层，糖衣片由内向外逐步渗透吸湿，糖衣层潮解易造成隔离作用减弱，甚至发生吸潮霉变或胀片爆裂的现象。如果素片水分过少，则素片的硬度会较差，在滚转时易引起素片的破碎或片粉的脱落。如果素片的水分分布不均匀，则会造成糖衣层的局部隔离层潮解，进而导致出现花斑现象。

因此，在素片干燥过程中及时监测其含水量是十分必要的。目前，素片水分测定的药典方法为干燥失重法，该方法的分析时间较长，且具有破坏性，不适用于工业现场的在线分析。而近红外光谱技术能够无损快速地测定素片的含水量，可用于快速检测。曾有文献报道过，运用近红外光谱技术来控制样本水分取得了较满意的结果。

但是，现有技术并没有提供关于工业现场在线分析素片含水量的任何方法，也没有提供关于在素片干燥过程中在线监测和/或反馈其含水量是否合格的任何方法。

发明内容

发明目的

为了解决上述技术问题，并实现在素片干燥过程中在线监测和/或反馈其含水量是否合格的目的，本发明提供了一种素片干燥过程的定性分析方法。

技术方案

为了实现上述目的，本发明提供了一种素片干燥过程的定性分析方法，该方法包括以下步骤：

(1)采集光谱：提供素片含水量合格与不合格的批量样本，采集各样本的近红外光谱数据，其中所述近红外光谱数据均为全波长数据，即NIR光谱扫描范围为4,000～10,000cm^-1；

(2)变量筛选：使用主成分分析法对全部样本的近红外光谱数据进行波段的选择，筛选出与素片水分相关性大的波段；

(3)光谱预处理：对步骤(2)中所选出波段的光谱数据进行预处理，以增强有效信号并提高信噪比；

(4)选取训练集样本：使用KS法选取步骤(3)中所有近红外光谱数据中的一部分作为训练集样本，将上述数据中的另一部分作为预测集样本；

(5)建立并评价定性模型：使用步骤(4)中训练集样本的数据通过支持向量机(SVM)法建立定性模型，然后使用上述定性模型对步骤(4)中的预测集样本进行含水量监测，以此评价所建的定性模型；和

(6)使用上述定性模型对其它素片样本的含水量进行监测。

具体地，在上述变量筛选步骤(步骤(2))中，使用主成分分析法进行光谱波段的选择。主成分分析法是将原始指标重新组合成一组新的互相无关的综合指标(主成分)来代替原始指标，同时根据实际需要从中选取几个尽可能保留原始指标的信息、且互不相关的主成分的过程。确定主成分的个数是其中的关键环节。本发明提供的方法采用交叉验证准确率来考察模型性能。随着主成分个数的增加，模型的交叉验证准确率增大，当主成分个数达到一定数值时，模型的交叉验证准确率达到最大值。选择与这些主成分相关的波段作为建模波段，不仅可以剔除不相关的波段，并且可以大大简化模型，从而实现压缩变量的目的，并提高模型的学习精度与预测能力。

在该步骤中，选择波段的依据是：NIR谱图上既包括有用信息也包括干扰信息，全谱建模由于干扰信息的存在不仅使运算时间增加，同时还会使定性模型的预测精度降低。优选的波段就是那些包含大量丰富的有用信息，而噪音等干扰信息影响较小的波段。主成分分析法的原理是对全谱建立主成分模型，选出特征信息量最丰富且对两组样本的区分度最大的波段。

在上述光谱预处理步骤(步骤(3))中，对不同的预处理方法进行了优化组合，并对预测结果进行比较，最终确定预处理方法选自标准化法、均值中心化法、多元散射校正法及其组合中。例如，可以是其中的一种或者两种及两种以上的组合。优选为标准化法与多元散射校正法的组合。

在上述选取训练集样本步骤(步骤(4))中，使用KS法进行训练集样本的选择。即将所有的样本都看作训练集候选样本，用每个样本的自变量计算样本之间的欧氏距离，在选择下一个候选样本进入训练集之前要求被选入的样本与已经选入训练集的样本之间的欧氏距离最大，并基于此选取部分样本作为训练集样本，使之与作为预测集的剩余部分样本有相同的或相近的空间分布。与随机选取法相比，KS法选取的训练集更具代表性，有利于得到外推能力更强的模型。优选地，取所有样本的2/3作为训练集样本，1/3作为预测集样本。

在上述建立并评价定性模型步骤(步骤(5))中，支持向量机(SVM)法的核函数为径向基核函数；其中，使用k折交叉验证准确率来确定径向基核函数的惩罚因子C和参数γ。优选地，使用网格搜索法结合k折交叉验证准确率来确定径向基核函数的惩罚因子C和参数γ。SVM法以训练误差作为优化问题的约束条件，置信范围值最小化作为优化目标，即采用结构风险最小化准则的学习方法，提高了模型的泛化能力。

有益效果

本发明提供的方法通过结合近红外光谱、主成分分析法和支持向量机法等方法快速建立了素片含水量定性模型，实现了对素片含水量合格与否的快速评价，并且监测准确率基本达到了工业生产的要求。相对于传统方法或经验判断，本发明提供的定性分析方法能够实现在素片干燥过程中快速、无损、可靠、准确地在线监测和/或反馈其含水量是否合格的目的。

本发明提供的方法以支持向量机(SVM)法为建模基础，结合KS法以保证总体样本的广泛性与代表性，结合主成分分析法以保证各个样本的有效性与合理性，并通过对光谱的预处理进一步降低无关信息的干扰，最终获得预测能力良好的模型。

附图说明

图1为本发明实施例1的步骤(1)中各乳块消素片样本的近红外原始光谱图；

图2-1为本发明实施例1的步骤(2)中前六个主成分与吸收波长之间的相互关系的图；

图2-2为本发明实施例1的步骤(2)中后六个主成分与吸收波长之间的相互关系的图；

图3-1为本发明实施例1的步骤(3)中优选波段的近红外原始光谱图；

图3-2为图3-1中的光谱经标准化法处理后样本的近红外光谱图；

图3-3为图3-1中的光谱经均值中心化法处理后样本的近红外光谱图；

图3-4为图3-1中的光谱经多元散射校正法处理后样本的近红外光谱图；

图3-5为图3-1中的光谱经多元散射校正-标准化法处理后样本的近红外光谱图；

图3-6为图3-1中的光谱经多元散射校正-均值中心化法处理后样本的近红外光谱图；

图4为本发明实施例1的步骤(4)中经KS法选取的训练集样本和预测集样本在主成分空间的分布情况示意图，其中1为训练集样本，2为预测集样本；

图5-1为本发明实施例1的步骤(5)中经网格搜索后的计算结果示意图；

图5-2为本发明实施例1的步骤(5)中使用支持向量机法建立定性模型的参数C和γ的优化过程示意图；以及

图5-3为本发明实施例1的步骤(5)中乳块消素片在干燥过程中的含水量分布图，其中○表示未干燥样品的含水量，◇表示样品干燥65min后的含水量；□表示样品干燥150min的含水量，表示样品干燥235min的含水量

表示样品干燥360min的含水量。

具体实施方式

本发明提供的定性分析方法可以应用于各种不同的药物制剂中所包括的素片，并且不受包衣类型(如糖衣或薄膜包衣)的限制。

下文，将以乳块消素片(中成药)为例，举例说明本发明提供的定性分析方法。应该理解，下述优选实施方式只是对本发明进行的例证性说明，并不会限制本发明。

实施例1

1仪器

Antaris傅立叶变换近红外光谱仪(美国Thermo Nicolet公司制造)配有InGaAs检测器、积分球漫反射采样系统、Result操作软件和TQ Analyst V6光谱分析软件。

2NIR扫描条件

采用积分球漫反射采样系统，NIR光谱扫描范围4,000cm^-1～10,000cm^-1；扫描次数32；分辨率8cm^-1；以内置背景为参照。

3数据处理及软件

本实施例使用主成分分析法选择波段，对光谱进行适当的预处理，使用KS法选取训练集样本和预测集样本，并使用SVM法建立近红外光谱的定性模型，使用训练集k折交叉验证(k-fold cross-validation)准确率来确定径向基核函数的参数γ和惩罚因子C，模型的学习精度和预测能力分别用训练集交叉验证准确率和预测集交叉验证准确率进行考察。其中，所述训练集交叉验证准确率是这样获得的：将训练集分为k份，其中k-1份作为训练数据集，训练得到模型，而另外的1份作检验，计算准确率；这样重复n次，选取准确率最高的模型，该准确率即为训练集交叉验证准确率，此值表征的是模型拟合的好坏。

本实施例使用的是台湾大学林智仁教授编制的支持向量机软件libsvm-2.89(参见，Chih-Chung Chang，Chih-Jen Lin，LIBSVM：a library forsupport vector machines，2001，http://www.csie.ntu.edu.tw/～cjlin/libsvm)，PCA工具包使用的是由

等人提供的网络共享(http://www.models.kvl.dk/source/iToolbox/)，SPXY法由

提供(参见，Roberto Kawakami Harrop

Mário César Ugulino Araujo，Gledson EmidioJosé，et al.A method for calibration and validation subset partitioning[J].Talanta，2005，67：736-740)，采用MATLAB 7.1软件工具(Mathwork Inc.)计算。

步骤(1)：采集光谱

由“北京中医药大学药厂”提供210个乳块消素片样本。

使用Antaris傅立叶变换近红外光谱仪采集上述各样本的近红外光谱数据，全部样本的数据如图1所示。如图1所示，各样本的原始光谱图严重重叠，无法直接区分出合格和不合格样本，需要进一步的数据分析。

步骤(2)：变量筛选

本发明使用主成分分析法，以原始光谱数据为输入数据，以交叉验证准确率为指标，筛选出与素片含水量相关的变量。

下表1给出了主成分分析的结果，前两个主成分的方差贡献率达到95％以上，可以很好地表征原始光谱的整体信息，但其建模效果却不理想。究其原因，样本之间的差异主要是由于水分和水分相关的性质引起的，但也存在部分差异来自于其它性质，如化学成分的含量。从图2-1中可以看出，第一个主成分与所有波长点关系均较为密切，故可能是因为包含了很多干扰信息，造成交叉验证准确率较低。

进行主成分分析时，仅是针对光谱数据进行的分析，不考虑因变量的影响，这增加了出现较大偏差的风险。采用主成分分析法处理数据，减少了测量噪声和共线性的信息，但没能把专属于目标性质的信息有效地提取出来。为了能够获得仅与水分或含水量有关的光谱信息，对12个主成分与吸收波长之间的关系做了进一步的分析。水分是影响样本光谱信号的首要因素，波长与主成分的相关系数越大，其与水分有关的可能性越大。从图2-1和2-2中可以看出，与12个主成分相关系数最大的波长是4004～4057cm^-1，4123～4126cm^-1，4250cm^-1，4319～4339cm^-1，4405～4493cm^-1，4921～4925cm^-1，4933cm^-1，4936cm^-1，4952cm^-1，5064～5238cm^-1，5762～5789cm^-1，6977～7035cm^-1，这些波长是影响综合指标的主要因素，用上述波段建模，交叉验证准确率达到94.29％。与文献报道(《近红外光谱解析实用指南》，Practical Guide toInterpretive Near-Infrared Spectroscopy，Jerry Workman，Ji.和Lois Weyer编著，褚小立、许育鹏、田高友译，化学工业出版社，2009年10月，57～62页)相比，5064～5238cm^-1和6977～7035cm^-1可能属于水分的响应区，仅以5064～5238cm^-1和6977～7035cm^-1建模，交叉验证准确率88.10％，其它的波段有助于提高分类的准确性，故推测可能反映的是受水分影响的性质的信号。

表1：各样本的近红外原始光谱主成分分析的结果

步骤(3)：光谱预处理

本实施例比较了标准化法(Autoscaling)、均值中心化法(Mean centering)、多元散射校正法(MSC)及其组合对优选波段所建模型性能的影响，结果如表2所示，并得到处理后的光谱数据图3-1至图3-6。

表2：不同数据预处理方法对结果的影响

为改善模型性能，需要对光谱数据进行预处理。从上表2可以看出，采用多元散射校正法与标准化法结合的方式取得的效果最好。

图3-1为本发明实施例1的步骤(3)中优选波段(即与12个主成分相关系数最大的波长)的近红外原始光谱图。从图3-1可以看出，即使是与目标性质相关的波段，其原始光谱之间差异也不明显，谱带重叠较为严重。

图3-2为图3-1中的光谱经标准化法处理后样本的近红外光谱图。从图3-2可以看出，经标准化法处理后的近红外光谱图的y轴取值，较原始光谱的y值在更大的范围内变化。标准化法属于数据增强算法，对整体的光谱信号起到增强作用，放大样本之间的细小差异。与图3-1所示的原始光谱图相比，图3-2中经标准化法处理后的光谱图，各个样本之间的距离增大了，但两组样本(合格品和不合格品)仍存在谱带重叠的情况。从建模效果上看，经标准化法处理后的光谱优于原始光谱，略差于经多元散射校正-标准化法处理后的光谱。究其原因，造成样本光谱之间差异的因素除了水分，还有颜色、密度、光散射、仪器噪音等，经标准化法处理后这些信号均被不同程度的放大了，影响了光谱的分布和模型性能。

图3-3和图3-6分别为图3-1中的光谱分别经均值中心化法和多元散射校正-均值中心化法处理后样本的近红外光谱图。均值中心化法是在原始光谱上扣除平均谱，使其变量均值为零，并使谱图在x轴上下分布更加均匀。这种处理方法简化了模型，但如表2所示，这种处理方法没有改善模型的效果。

图3-4为图3-1中的光谱经多元散射校正法处理后样本的近红外光谱图。多元散射校正可以修正光谱中因样本不均匀产生的基线平移和偏移。经多元散射校正法处理后的光谱仅消除了由光散射产生的噪音，其它吸收信号几乎不受影响，提高了信噪比，模型性能有所提高(见表2)。

图3-5为图3-1中的光谱经多元散射校正-标准化法处理后样本的近红外光谱图。

如表2所示，在所有预处理方法中，采用多元散射校正法和标准化法二者组合的建模效果最好，优于单独使用多元散射校正法或标准化法处理的模型。这可能是因为先经多元散射校正法处理，消除了样本内部不均匀性对光谱的影响，减少了背景干扰，起到了减少噪音的作用；再采用标准化法处理，使得所有数据具有相同的基点和变化范围，增强了样本之间的可比性，起到了放大信号的作用。两种预处理方法的结合大大增加了信噪比。

步骤(4)：选取训练集样本

本实施例使用KS法进行训练集样本的选择。KS法将所有的样本都看作训练集候选样本，用每个样本的自变量计算样本之间的欧氏距离，在选择下一个候选样本进入训练集之前要求被选入的样本与已经选入训练集的样本之间的欧氏距离最大，并基于此选取部分样本作为训练集样本，使之与作为预测集的剩余部分样本有相同的或相近的空间分布。实施过程如下：

首先，根据下述样本距离公式(1)选择欧氏距离最远的两个向量进入训练集：

$d_{i,j}=\sqrt{{(x_i-x_j)}^2}=|x_i-x_j|;i,j\∈(1,n)---\left(1\right);$

其中，x为样本向量；d_i，j为第i个样本向量到第j个样本向量的欧氏距离。

对每个候选样本向量i₀，计算它与已选入训练集中每个样本向量欧氏距离的最小值

然后将具有最大的最小距离max_i0(min_i(d_i，i0))的样本向量选入训练集。重复该过程直到满足训练集对样本数的要求。

本实施例共提供了210个样本，选取其中的2/3(即140个样本)作为训练集样本，另1/3(即70个样本)作为预测集样本，其计算结果如下表3所示。图4为经KS法选取的训练集样本和预测集样本在主成分空间的分布情况示意图。从图4可以看出，前两个主成分的累积方差贡献率大于90％，可以很好的表征原谱信息。由该图还可以看出，所有预测集样本均落在训练集样本范围内，实现了训练集样本在整体样本空间均匀分布。

表3：KS法计算结果

(5)建立并评价定性模型

SVM法是Vapnik等人基于统计学习理论(Statistical Learning Theory，SLT)提出的一种新的机器学习算法(参见，Vapnik V.Statistical Learning Theory，John Wiley，New York，1998)。此前算法的重要基础是传统的统计学，前提是有足够多样本，当样本有限时难以取得理想的效果。而基于统计学习理论的支持向量机(SVM)，以训练集误差作为优化问题的约束条件，以置信范围值最小化作为优化目标。预期得到的是一个对所有样本预测准确率都很高的模型，但如果仅考虑训练集误差，有可能出现如下情况：训练集误差很小(即训练集样本的准确率很高)，保证了分类精度很高，但由于学习机器过于复杂，置信范围增大，推广性很差，对非训练集样本的预测能力很差。所以，要同时保证训练集样本和预测集样本的准确率，即采用结构风险最小化准则的学习方法，提高了模型的泛化能力。SVM的基本思想是把原始特征数据从输入空间非线性地映射到一个高维特征空间(Hilbert空间)，然后在此空间中求解凸优化问题(典型二次规化问题)，可以得到唯一的全局最优解。

给定一训练集{(x_i，y_i)，i＝1，2，3，...，n}，其中y_i∈{-1，1}表示任一样本x_i的分类标识。如果训练集是线性可分的，SVM就是寻求超平面：

f(x)＝ω·x+b    (2)

使正样本(y_i＝+1)和负样本(y_i＝-1)可分，且使其边界上的点到该超曲面的距离最大。

由于很多两类情形并非线性可分，为此，SVM将样本点x(各样本的近红外光谱数据)通过函数φ(x)投影到高维特征空间以使其线性可分。但SVM并不是直接引入φ(x)，而是通过核函数K(x_i，x)方法间接引入的：

K(x_i，x)＝φ(x_i)·φ(x)(3)

核函数包括线性、径向基(RBF)、多项式和Sigmoid等多种形式，可根据实际的应用数据或环境不同而选用适合的核函数，大多数文献使用径向基核函数，即：

K(x_i·x_j)＝exp(-γ||x_i-x_j||²)，γ＞0(4)

在公式(4)中，x和y分别表示不同样本的测量数据(x为样本的光谱值，y为样本属性值)，γ为径向基核函数的可调参数，其数值需要在模型优化的过程中确定。

系数ω和b由极小化风险得到

$R_{{\mathrm{SVM}}_s}\left(C\right)=C\frac{1}{n}{\mathrm{\Σ}}_{i=1}^n{L}_{\mathrm{\ϵ}}(d_i,y_j+\frac{1}{2}{\left|\mathrm{\ω}\right|}^2)---\left(5\right)$

其中ε为不敏感损失函数，C为惩罚因子，表示对ε的样本的惩罚，C值越小，惩罚越小，训练误差越大；C值太大，训练误差虽然小，但是学习机器过于复杂，置信区间将随之增大，模型外推能力下降，故C值的大小需要在模型优化过程中确定。

使用k折交叉验证(k-fold cross-validation)准确率来确定径向基核函数的参数γ和惩罚因子C。k折交叉验证(k-fold cross-validation)准确率(即训练集交叉验证准确率)是这样获得的：将训练集分为k份，其中k-1份作为训练数据集，训练得到模型，而另外的1份作检验，计算准确率；这样重复n次，选取准确率最高的模型，该次准确率即为训练集交叉验证准确率，此值表征的是模型拟合的好坏。

本实施例采用网格搜索法结合5折交叉验证法优化参数C和γ。网格搜索法是将以C和γ为横纵坐标的平面划分成网格，将(C，γ)看做是网格上的节点，按照设定的步长来进行计算各点交叉验证准确率而找出最大值。C值的范围设为为-5～15，γ值的范围设为3～-15，步长均为2，平面上共有110种C-γ网格点，分别计算并比较不同C-γ下模型的5折交叉验证准确率，其计算结果如图5-1所示。从图5-1可以看出，γ值对建模效果的影响较为显著。γ＝8时，无论C取何值，建模结果均很差。C值代表惩罚程度的大小，C值越大，训练误差越小，但C值过大，会使模型过分学习，导致预测能力的下降，所以大部分情况曲线呈现以下趋势：准确率先随着C值的增加而提高，C到达一定值后准确率基本不再变化。γ值不同，其达到最大准确率的C值也不同。图5-2示出了5折交叉验证准确率在95.00％以上时的C-γ值，最终获得最佳的参数组合为C＝32768，γ＝0.00048828125。

在上述优化参数下，定性模型训练集5折交叉验证准确率92.14％(n＝140)，训练集样本的预测准确率为100％(n＝140)，预测集样本的预测准确率为98.57％(n＝70)。“训练集样本预测准确率-预测集样本预测准确率”较小，说明此法选取的训练集样本具有代表性，所建的模型同样适用于预测集样本。预测集样本的预测准确率高，说明模型的预测能力好。

为了进一步分析建模结果，对样品含水量进行了测量，共11个样品，置于40℃烘箱内，干燥不同时间，并间隔取样。实验结果如图5-3所示，合格品和低水分样品含水量较接近，可能是导致错判的原因。

步骤(6)：使用定性模型监测含水量

使用上述定性模型，通过预测集样本的近红外光谱数据，监测70个预测集样本。

在上述优化参数下，对70个预测集样本(包括43个合格品和27个不合格品)进行监测，错判个数为1，准确率为98.57％。

由上述实验结果看见，在各项优化参数的条件下模型对样本的监测准确率基本达到工业生产的要求。相对于传统方法或经验判断，本发明提供的定性分析方法能够实现在素片干燥过程中快速、无损、可靠、准确地在线监测和/或反馈其含水量是否合格的目的。

Claims (7)

1.一种素片干燥过程的定性分析方法，该方法包括以下步骤：

(1)采集光谱：提供素片含水量合格与不合格的批量样本，采集各样本的近红外光谱数据，其中所述近红外光谱数据均为全波长数据；

(2)变量筛选：使用主成分分析法对全部样本的近红外光谱数据进行波段的选择，筛选出与素片水分相关性大的波段；

(3)光谱预处理：对步骤(2)中所选出波段的光谱数据进行预处理，以增强有效信号并提高信噪比；

(4)选取训练集样本：使用KS法选取步骤(3)中所有近红外光谱数据中的一部分作为训练集样本，将上述数据中的另一部分作为预测集样本；

(5)建立并评价定性模型：使用步骤(4)中训练集样本的数据通过支持向量机法建立定性模型，然后使用上述定性模型对步骤(4)中的预测集样本进行含水量监测，以此评价所建的定性模型；和

(6)使用上述定性模型对其它素片样本的含水量进行监测。

2.根据权利要求1所述的方法，其中，在步骤(3)中，所述预处理的方法选自标准化法、均值中心化法、多元散射校正法及其组合中。

3.根据权利要求2所述的方法，其中，在步骤(3)中，所述预处理的方法为标准化法与多元散射校正法的组合。

4.根据权利要求1所述的方法，其中，在步骤(4)中，使用KS法选取步骤(3)中所有近红外光谱数据中的2/3作为训练集样本，上述数据中的另1/3作为预测集样本。

5.根据权利要求1所述的方法，其中，在步骤(5)中，所述支持向量机法的核函数为径向基核函数。

6.根据权利要求5所述的方法，其中，使用k折交叉验证准确率来确定径向基核函数的惩罚因子C和参数γ。

7.根据权利要求6所述的方法，其中，使用网格搜索法结合k折交叉验证准确率来确定径向基核函数的惩罚因子C和参数γ。

Images