梯形公式预测校正matlab_近红外光谱快速预测天舒片包衣终点研究

摘 要：目的  应用近红外光谱(NIRS)技术快速预测天舒片薄膜包衣终点。方法首先，采用组合间隔偏最小二乘法(synergy interval partial least square，siPLS)优选建模区间，并对预处理方法进行筛选；其次，根据优选的建模条件，基于建立的参考光谱库，采用一致性检验和相似度匹配建立校正模型；随后，通过对比测试光谱与参考光谱间的相似度大小，设定阈值，确定相似度、阈值与包衣终点的关系；最后，采用验证集光谱对模型性能进行验证。结果  在包衣过程中，当连续多条测试光谱的相似度超过阈值时，指示为包衣终点。结论一致性检验和相似度匹配模型性能良好，能够敏锐地反映包衣薄膜的均匀性，并准确预测包衣终点。该研究对减小包衣的批间差异性和降低原料损耗，提高中药固体制剂的包衣生产效率和提升质量均一性具有重要意义。

天舒片由川芎和天麻2味药组成，临床用于治疗神经性偏头痛及高血压性头痛^[1]。为了使片剂免受光致分解或酸降解等外部影响，掩盖味道，或控制药效成分的定点定时释放，会对片剂进行包衣。包衣是一个多单元多变量的操作过程，由于包材喷涂、分布和干燥同时进行，许多关键工艺参数(critical process parameters，CPPs)都会影响片剂包衣的生长和均匀性。喷枪速率的变化会显著影响片剂到达一定厚度所需的时间^[2]。进出气的温度和喷枪速率是影响高分子聚合物成膜的重要因素。当达到最低成膜温度时，需确保片剂有适当的含水率^[3]。在工业生产时确保包衣厚度和包衣质量的均一性显得至关重要。包衣厚度的测量并非如想象中那般容易。研究表明，片剂中心和侧面的包衣厚度是存在显著性差异的^[4]。包衣时，片心粗糙的表面在主包衣层形成前就已完成了填充^[5]，因此，根据包衣前后的片剂平均质量计算厚度的预测值相比直接测量的预测值会更准确些^[6]。然而，这些方法均在少量样品上重复测量，且耗时耗力。因此更倾向于一些适合大量样品的快捷的研究方法，如X射线计算机断层技术^[4]、声显微技术^[7]。包衣过程中，一些聚合物溶液是直接喷雾干燥的，并非在片剂表面形成包衣层，这可能导致不均匀包衣层^[8]。当出现这种情况时，通常选择继续包衣，延迟包衣终点，最终包衣过厚，影响了药物的溶出。

一致性检验(conformity test)和相似度匹配(similarity match)属于定性研究方法。一致性检验主要用于药品的真伪鉴别^[9-10]，厂家、产地的归属鉴定^[11-12]及药品供应链全程质量监控^[13-14]及生物等效性分析^[15]等方面。相似度匹配则广泛应用于药材厂家或产地鉴别^[16-17]、真伪判断及药材质量控制^[18]研究。目前鲜有将这2种方法应用于片剂薄膜包衣领域的研究报道。本研究以天舒片为研究对象，收集参考样品(包衣终点样品)和测试样品(不同包衣时间的样品)，采集近红外漫反射光谱；并采用一致性检验和相似度匹配分析方法建立校正模型。以期该方法在包衣过程中的应用能够提升片剂包衣薄膜的均匀性，实现准确的终点判别。

1  仪器与材料

1.1  仪器

BGB-150 C型高效包衣机，宜春万申制药机械有限公司；Bruker MPA型傅里叶变换近红外光谱仪，配积分球漫反射检测器，德国布鲁克光学仪器公司。

1.2  材料

天舒片由江苏康缘药业股份有限公司提供，其片心为直径10 mm、厚度4 mm的双凸片；蓝色胃溶型包衣薄膜预混剂，批号F1810328，购自上海新菲尔物科技有限公司。

2  方法

2.1  包衣过程和取样方法

1批天舒片心分成4份(110、110、100、100 kg)进行薄膜包衣，随机取一份天舒片心投入包衣锅内。包衣开始后，每隔30min从设备预留取样口迅速取样，每次取约100片。取样后样品用自封袋密闭保存。本研究共完成9个批次的取样，制剂批号分别为180809、180901、180902、180903、180904、180905、180906、180907、180908，见表1。本次共收集了68个时间点的过程样品和9个包衣终点的样品^[19]。

2.2  近红外光谱采集

近红外光谱仪波数为10 000～4 000 cm⁻¹，分辨率为16 cm⁻¹，以空气为参比扣除背景，每条NIRS为64次扫描的平均值。取片心和每个时间点的包衣片10片进行近红外漫反射光谱采集。将手持式近红外探头贴紧药片凸面，使近红外光束集中在凸面的中心位置，每片正反面分别重复扫描3次，每片扫描后获得6条(2面×3次/面)NIRS。

2.3  参考光谱库的建立

以批次180809、180901、180902、180903、180904、180905、180906、180907、180908的包衣终点样品的NIRS建立参考光谱库。

2.4  样品集划分

按制剂批次顺序，以批次180809、180901、180902、180903、180904、180905共45个时间点的过程样品为校正集，以批次180906、180907、180908共23个时间点的过程样品为验证集。

2.5  一致性检验

一致性检验^[10]可用于快速比较测试光谱与参考光谱间的相似性，其原理为计算参考光谱集每个波长点的平均值和标准偏差(S)，以波长点的平均值加减一定倍数的标准偏差作为该波长点的控制区间；如测试光谱每个波长点的吸光度没有超出预定的控制范围，则说明测试样品与参考样品具有相似性。控制区间计算公式为参考光谱吸光度平均值±CI限度值×S，其中CI(conformity index)表示一致性指数，可通过以下步骤确定：(1)计算参考光谱每个波长点i处吸光度的平均值A_i,参和标准偏差S_i,参；(2)计算样品每个波长点处的CI_i值，计算公式：CI_i＝(A_i,样－A_i,参)/S_i,参，其中，A_i,样为样品测试光谱在波长点i处吸光度；(3)样品的CI_样值为该样品在每个波长点处的CI_i值的最大值，即CI_样＝Max(CI_i,样)，同样的，参考光谱的CI_参值为该参考样品在每个波长点处的CI_i值的最大值，即CI_参＝Max(CI_i,参)。

在建模阶段，通过比较合格光谱和不合格光谱的CI值，选择合适的CI限度值。而在应用阶段可直接将测试光谱的CI值与之前设定的CI限度进行比较，从而快速简单地判断。

2.6  相似度匹配

相似度匹配是一种光谱分类技术，通过比较测试光谱与参考光谱在波段范围内的光谱信息得到两者间的匹配程度。相似度匹配是基于Gram-Schmidt矢量正交方法整体评价相似度^[20-21]，基本原理如下。设参考图谱库是包含m条光谱且波长点数为n的矩阵X，将其正交分解为得分矩阵T和载荷矩阵P，其中r为矩阵X的秩，见式(1)。

其中p_j＝[ p_j1，p_j2，…，p_jn](j＝1，2，…，r)为一组相互正交的向量。由于p由正交向量组成，所以T＝XP^T。测试光谱数据为x_new，可利用式(2)计算得分向量t_new。

t_new＝x_newP^T＝[s₁，s₂，…，s_n]               (2)

j＝1时，根据e_i1＝|s₁－t_il|²(i＝1，2，…，m)计算得到e₁＝[e₁₁，e₂₁，e₃₁，…，e_m1]，并得到向量e₁中值最小的元素e_k1，对应的t_k1记为t₁′，接着，相应地计算其余光谱的t₂′，t₃′，…，t_r′，得到测试光谱与参考光谱间的残差光谱e，见式(3)，该残差光谱反映了测试光谱与参考光谱间的差异。

计算测试光谱与参考光谱之间的相似度匹配值(similarity match value，SMV)，见公式(4)。

SMV＝(1－||e||/||x_new||)×100               (4)

e为测试样品的残差光谱，x_new为测试光谱数据。SMV反映了测试光谱的残差平方和。SMV分布在0～100，越接近100，说明测试样品光谱与参考光谱信息越接近。

2.7  组合间隔偏最小二乘法^[22-23]

组合间隔偏最小二乘法(synergy interval partial least square，siPLS)是对间隔偏最小二乘法(intervalpartial least square，iPLS)的扩展，能够弥补iPLS一个子区间建模产生的不足。siPLS建模的一般步骤如下：(1)将全光谱等分为N个子区间；(2)将N个子区间排列组合成不同的子区间(一般为2、3或4个子区间进行组合)，并用于建立偏最小二乘模型；(3)根据交叉验证均方根误差(rootmean square error of cross validation，RMSECV)选择精度最高的模型，并选择相应的组合区间，RMSECV最小时的模型的精度最高；(4)优化所选组合区间，以步骤3中选择的区间为中心，对波长区域进行单向或双向扩展，得到一个最佳波长区间；(5)最后将具有高精度模型的区间进行建模预测，并根据最小RMSECV确定最优组合区间。

2.8  模型评价参数

基于所建的参考光谱库，以校正集样品的NIRS为测试光谱，建立定性校正模型。根据设定的阈值，以模型预测的正确率(W)^[16-17]为指标优化建模参数，考察模型性能，W的计算方法见公式(5)。其中，n为正确预测的样品数，N为样品总数。

W＝n/N                                 (5)

模型预测标准：当样品的光谱相似度小于相似度阈值时，模型预测为过程样品；当光谱相似度均大于或等于相似度阈值时，模型预测其为包衣终点的样品。

2.9  数据分析

一致性检验在OPUS化学计量学软件(version7.5，德国布鲁克公司)上完成；相似度匹配在TQ Analyst 光谱分析软件(version 9.0，美国Thermo公司)上完成；波段筛选在MATLAB软件(version 2016 b，美国MathWorks公司)上完成。

3  结果

3.1  一致性检验

3.1.1  波段选择 NIRS主要是分子倍频和合频的吸收，除了包含片心和高分子薄膜的化学信息外，还包含大量重叠信息和噪音，采用全光谱建模时，模型的预测精度未必能达到最佳值。因此，需要剔除冗余信息，提高分析模型的预测精度和建模效率。

基于所建的参考光谱库，以校正集原始光谱为测试光谱，光谱范围为10 000～4 000 cm⁻¹的建模参数建立模型，得到校正集光谱的CI值。以每个时间点样品的NIRS吸光度为自变量，对应光谱的CI值为因变量，构建45个样品光谱组成的数据集建立siPLS模型。将近红外全光谱等分为20个子区间，采用3个子区间组合的形式筛选波段。所得区间为7 290.06～6 996.91 cm⁻¹、5 785.76～5 477.19 cm⁻¹、4 597.75～4 312.32 cm⁻¹。从片心、包衣过程样品和终点样品光谱图(图1)上可以看出，包衣过程光谱和终点光谱在7 290.06～6 996.91 cm⁻¹和       4 597.75～4 312.32 cm⁻¹ 2个区间内有2个尖峰，而片心的光谱图上未出现该峰，推测尖峰可能是包衣薄膜的特征峰。7 290.06～6 996.91 cm⁻¹的差异对应C-H一级倍频伸缩振动基频与C-H变形振动基频组合频；在4 597.75～4 312.32 cm⁻¹的光谱差异属于C-H伸缩振动基频与C-H变形振动基频组合频和C-H二级倍频伸缩振动基频。而5 785.76～5 477.19 cm⁻¹的信息在光谱图上没有表现出十分显著性的差异，但与文献报道^[24]的包衣辅料5 986～5 369 cm⁻¹的特征吸收区相吻合，推断为包衣薄膜的特征峰。

3.1.2 预处理方法选择  受背景等环境因素的影响及仪器性能的制约，采集的光谱包含了噪音和失真信息，有必要对其进行预处理，剔除冗余信息，强化特征吸收。基于所建立的参考光谱库，以校正集光谱为测试光谱，光谱范围为优选的组合区间，分别在不同的预处理方法条件下建立模型。本研究考察了标准正态变换(standard normal variate，SNV)与导数法(一阶导数、二阶导数)结合对光谱预处理。SNV可用来消除固体颗粒大小，表面散射以及光程变化对漫反射光谱的影响，导数法可以消除基线偏移和背景干扰，分辨重叠峰，提高分辨度^[25]。不同预处理方法的建模结果见表2。本研究中，CI限度值设定为所有参考光谱CI值的最大值，即CI_限＝Max(CI_参)。

根据一致性检验预测结果可知，当一阶导数结合SNV处理光谱时，模型预测正确率最高。此时，CI_限＝Max(CI_参)＝3.3。

3.1.3  一致性检验模型建立  基于参考光谱库，以校正集光谱为测试光谱，在7 290.06～6 996.91 cm⁻¹、5 785.76～5 477.19 cm⁻¹、4 597.75～4 312.32 cm⁻¹ 3个波段内，以一阶导数＋SNV为预处理方法，CI_限为3.3的条件下进行一致性检验。

结果表明，该模型预测校正集光谱均为包衣过程中的光谱。由图2可知，CI_限为3.3时，校正集光谱和参考光谱区分明显，从整体看，随着包衣时间延长，校正集光谱的CI值向CI_限趋近，表明校正集光谱与参考光谱的相似度越高，越靠近包衣终点。此时，会出现一些光谱的CI值在CI_限以下，主要原因是包衣质量不均匀，仍需继续包衣，直至光谱CI值均在CI_限以下时，指示为包衣终点。

3.1.4  模型验证 采用验证集光谱对所建的校正模型进行验证，结果见图3，模型预测正确率为95.65%，表明所建的校正模型可以用于天舒片包衣过程样品与终点样品间的相似度分析，模型预测性良好，为提升包衣质量的均一性，准确预测包衣终点提供技术指导。

3.2  相似度匹配

3.2.1  波段选择 波段筛选的步骤如下。

首先，从所建的参考光谱库选出180906、180907、180908批次的参考光谱和校正集光谱作为测试光谱，以剩余批次的参考光谱作为参考光谱库，光谱范围为10 000～4 000 cm⁻¹的建模参数建立模型，得到每条NIRS的SMV。180906、180907、180908批次参考光谱的SMV用于计算模型的SMV阈值。

其次，以每个时间点样品的NIRS吸光度为自变量，光谱所对应的SMV为因变量，构建45个样品组成的数据集建立siPLS模型。将近红外全光谱等分为20个子区间，联合区间数设为3。一共得到20个回归模型。模型性能高低主要通过RMSECV判断，选取RMSECV最小的区间为优选区间。优选区间为9 396.074～9 102.93 cm⁻¹、6 680.625～   6 395.194 cm⁻¹、4 597.75～4 312.32 cm⁻¹。

其中，9 396.074～9 102.93 cm⁻¹和6 680.625～6 395.194 cm⁻¹均无法反映图 1中2个尖峰的信息。在此基础上，分别采用全光谱10 000～4 000 cm⁻¹，此次筛选区间9 396.074～9102.93 cm⁻¹、6 680.625～6 395.194 cm⁻¹、4 597.75～4 312.32 cm⁻¹(区间1)和一致性检验优选区间7 290.06～6 996.91 cm⁻¹、 5 785.76～5 477.19 cm⁻¹、4 597.75～4 312.32 cm⁻¹(区间2)这3个不同的区间建立相似度匹配模型，考察不同组合区间对模型性能的影响，结果见图4，10 000～4 000 cm⁻¹所建模型中，校正集光谱的SMV随时间变化趋势更明显，随着包衣时间的延长，与参考光谱的匹配度越高。全光谱所反映的样品信息丰富，能更好地呈现SMV的变化趋势，模型的灵敏度高。故选择全光谱10 000～4 000 cm⁻¹建模。

3.2.2  预处理方法选择  分别采用SNV、导数法(一阶导数、二阶导数)结合平滑法Savitzky-Golay(SG)滤波器平滑和NorrisDerivative(ND)滤波器平滑进行光谱预处理，结果见表3。

采用SNV预处理光谱所建的模型无法区分校正集光谱和参考光谱。导数法结合平滑法处理光谱时，从100 min起，校正集光谱和参考光谱的SMV较接近，但二者无法清晰区分开，包衣终点趋势不明显。SNV＋二阶导数预处理光谱时，通过SMV能够显著区分校正集光谱和参考光谱。180906、180907、180908批次参考光谱的SMV正态检验结果表明，P＞0.05，表明参考光谱的SMV符合正态分布，设定校正模型SMV阈值为参考光谱SMV的置信下限，计算公式为－σ/n^1/2z_α/2(σ为标准差，n为光谱数量，z为统计量，α为置信水平)。根据不同预处理方法所建模型的结果可知，预处理方法为SNV＋二阶导数，SMV阈值为98.47时，模型预测的正确率最高。

3.2.3  相似度匹配模型建立  基于所建的参考光谱库，以校正集光谱为测试光谱，在10 000～4 000 cm⁻¹波段内，以SNV＋二阶导数为预处理方法，建立校正模型。模型的阈值设定为98.47。校正模型结果表明，校正集光谱均为包衣过程光谱。如图5所示，校正光谱的SMV受时间变化影响较大，包衣刚开始时，SMV较小，随着包衣时间延长，SMV逐渐增大，且越靠近阈值，表明越接近包衣终点。由于包衣机的转速、喷枪角度、喷液流量均可能影响包衣质量的均匀性，180901、180903、180904、180905在200 min附近时均有部分样品的SMV大于阈值，这是同一批样品间和不同批次间样品包衣时均可能出现的现象^[26]。此时需要继续包衣，直到光谱的SMV连续大于阈值时，才达到包衣终点。

3.2.4  模型验证 将验证集光谱导入校正模型中，对模型性能进行验证。模型验证结果如图6所示，该模型预测验证集均为包衣过程中的样品，正确率为100%。虽然180906批次在临近包衣终点时，由于包衣质量不均一，部分光谱的SMV大于98.47，但是模型仍将其判定为包衣过程光谱，表明所建的相似度匹配模型能灵敏得反映包衣质量的均匀性，并且准确地预测包衣终点，可为实际生产中制剂包衣终点的判断提供技术支持。

4  讨论

本研究采用一致性检验和相似度匹配的分析方法建立模型预测包衣终点。结果显示，所建模型预测精度和可靠性较高，能够快速判断天舒片包衣终点和提升包衣质量均一性。本研究方法对检验仪器的精密度和准确度要求不高，操作简单，适用包衣过程的快速检测。

本研究按制剂批次顺序，将样品集划分为校正集和验证集。这种划分方式，能够在包衣过程中观察到测试光谱与参考光谱的相似度变化趋势，有助于提升包衣质量的均一性，准确预测包衣终点。

本研究以校正集光谱的吸光度为自变量，相应的相似度评价指标为因变量，基于siPLS算法优选建模区间。一致性检验的优选区间与光谱信息呈现较好的对应关系，模型性能良好。然而，在相似度匹配模型中，采用筛选区间建立的模型无法凸显测试集光谱SMV的变化，模型灵敏度低，这可能与SMV的计算原理有关。SMV由测试光谱与参考光谱间的残差光谱计算而得。在全光谱中，残差光谱充分反映了测试光谱与参考光谱间差异的信息，模型灵敏度高。故在相似度匹配分析中采用全光谱建模。建模结果表明，采用全光谱所建的模型具有良好的预测性。

包衣并不是一个单因素生产单元，该过程受诸多因素的影响，研究的复杂性由此可见。对于功能性衣膜(具有缓控释功能)而言，例如，治疗2型糖尿病的固定剂量组合(fixed-dosecombination，FDC)由盐酸二甲双胍的片心，格列美脲的外层包衣层和用于分隔片心和外层的中间层组成^[6]。盐酸二甲双胍具有缓释功能，格列美脲被制备成速释制剂，可见包衣膜的薄厚和均一性会影响药物的功效。促进对包衣工序的研究，并将包衣过程与药物的崩解和溶出、药效成分的含量等相结合，强化对包衣过程的理解，逐渐将包衣生产过程透明化，促进制剂工艺的发展。

参考文献(略) 

来  源：夏春燕，徐芳芳，张  欣，徐  冰，包乐伟，王振中，肖  伟. 近红外光谱快速预测天舒片包衣终点研究 [J]. 中草药, 2019, 50(21):5223-5230.