一种色谱谱峰识别算法的实现

一种色谱谱峰识别算法的实现陈淑珍麻红昭(浙江大学化学工程与生物工程学系 浙江 杭州 310027)摘 要色谱在化学分析中的应用越来越广泛，针对色谱信号数据处理要求，实现一套谱峰识别的算法。在数据预处理和提出的程序法确定合适阈值的基础上，该算法结合谱峰曲线及其一阶导数曲线的趋势对谱峰进行识别，算法经宽松变换后，用于区分化学谱图中谱峰的终点是否真正落回到基线上，并由此实现一种具有健壮性的重叠峰检测方法和及时准确的谱图基线校正方法。由于 该算法的阈值参数可以固定，从而保证了重叠峰检测较少地受到人为的影响，进一步提高了化学谱图数据处理的自动化水平。关键词数据处理 色谱峰识别 阈值重叠峰检测 基线校准中图分类号文献标识码DOI:10 3969 / j issn 1000-386x 2013 11 084TP301AEALISATION OF A CHOMATOGAPHIC PEAK ECOGNITION ALGOITHMChen ShuzhenMa Hongzhao(Department of Chemical Engineering and Biochemical Engineering，Zhejiang University ，Hangzhou 310027，Zhejiang，China)The application of chromatography is more and more popular in chemical analysis In this article，we realise a set of algorithmsAbstractfor chromatographic peak recognition in light of the requirement on chromatogram signal processing Based on data pre-processing andappropriate threshold determined with procedural method，the algorithms distinguish the peaks by combining the trend of chromatogram curve with the trend of its first derivative curve，after loose transformation，the algorithms are applied to differentiate whether the end points of peaks in chemical chromatogram have really fallen back to the baseline，and thus implement an overlapping peaks detection method with robustness and the timely and accurate method of chromatogram baseline correction As the threshold parameter in the algorithm can be fixed，therefore it ensures the overlapping peaks detection less impacted by the human，and further improves its automation level in chemical chromatogram data processingKeywordsData processing Chromatogram peaks recognition Threshold value Overlapping peak detection Baseline correction线之下，则认为它是落回基线的峰结束点，此峰是独立峰。显然，用这种方法确定两峰是否重叠，受用户预先设置的基线漂移 值的影响，带有较大的随意性。而本文采用谱图曲线和一阶导 数曲线相结合的算法能够判断峰点是否落回基线上从而能自动 识别重叠峰及后续的基线校准。常用的一阶导数法不能描述复 杂的曲线变化情况如肩峰，二阶导数虽然描述精度提高，但也同 时提高了噪音的影响，以至于在判峰过程中，峰特征点容易被淹 没，阈值设置弹性小。为了克服上述问题，本算法对谱图曲线和 一阶导数趋势曲线的描述方法进行改进，使得算法对噪音和信 号的区分度增大，即算法的分辨率增大，从而保证了重叠峰自动 识别和基线自动校准的可行性，提高了化学谱图数据处理的自 动化水平。引言0色谱仪输出的模拟量信号，经 A / D 采样板转换后形成离散的数字量信号，称之为化学谱图。化学谱图特定位置上的信号 往往对应于某一具体物质或特定信息，通常以各种类型的谱峰 形式表现。谱峰信号的强度或大小又与此物质的组分、物理化 学性质、结构活性等有关，因此从信号的位置、大小和形状可以 获取物质的定量定性信息。谱峰的准确检测是化学数据处理获 取有效信息的关键。因为仪器的分辨率有限及背景信号漂移， 化学谱图中一般还存在着谱峰重叠、基线漂移的现象，所以，噪 声滤除、峰检测、重叠峰检测及基线校正、重叠峰分解在化学谱 图数据处理中都必不可少1。本文采用滤波和平滑处理对谱峰数据预处理后进行谱峰检 测。目前最常用的方法是一阶导数法及其各种变形方法2 4， 这类方法现有的大多数化学谱图数据处理方法不能判别谱峰的 终点是否真正落回到基线上。从而不能直接区分落回基线的峰 结束点和未落回基重叠峰间的谷点，因而不能直接检测两峰是 否重叠，需用户设置某个参数来判断。例如设置“基线漂移”参 数的方法，从峰起点划斜率等于此参数值的直线，若峰终点落于数据预处理1滤 波从色谱仪检测器输出的模拟信号经模数转换器转换为数字1 1收稿日期:2012 06 29。陈淑珍，硕士生，主研领域:信息过程自动化。麻红昭，副教授。310计算机应用与软件2013 年信号5，数字信号经过接口电路输入到计算机。这些信号常伴有高频噪声。因此，在谱峰检测前，需要对输入信号进行滤波处 理。本文算法采用移动平均法对采集到的色谱数据进行滤噪。 移动平均滤噪是一种比较简单的平滑方法。具体做法是把连续 2m + 1 个采样点看成一个队列，每次猜到一个新数据放入 队尾，并扔掉原来队首的一个数据( 先进先出原则)，由连续的2m + 1 个数据求得的算术平均值作为新的观测值。1 2 平滑处理信号经过滤波后，信号斜率曲线仍会有抖动，为减小抖动对 判峰的影响，计算斜率时可以进一步平滑。设采样时刻 tn 3 ， tn 2 ，tn 1 ，tn ，tn + 1 ，tn + 2 ，tn + 3 ，的采样值分别为 y n 3 ，yn 2 ， y n 1 ，yn ，yn + 1 ，yn + 2 ，yn + 3 ; t n 时刻色谱曲线的斜率值 Sn 计算 如下:时不满足条件 h h h 和条件 d d d ，则判断0 1n 10 1n 1为峰终点，a = 5。如果以上三种情况都不满足，则继续移动窗口。3yn + 3 + 2yn + 2 + yn + 1 yn 1 2yn 2 3yn 3(1)Sn=6图 1 重叠峰曲线及其一阶差分曲线式(1)并不是固定不变的，可以根据实际扰动情况进行调整。调整可以增加或减少加权采样值的个数，也可以改变各采 样值的权系数。算法条件的宽松变换及实现33 1 算法条件的宽松变换以上是理想峰的判定过程，但在峰实际处理过程中，色谱采 样数据中存在高频噪音，使得算法中的单调递增或单调递减很2谱峰检测算法设置一数据窗口其宽度，即窗口内包含的顺序采样数据个数为 n，将这 n 个色谱数据值记为 h0 ，h1 ，hn 1 ，窗口内最左边 数据 h0 。并计算出相应采样数据点的斜率 d0 ，d1 ，dn 1每次处理结束时窗口内数据向后移动一个，即把窗口内最 左边的一个数据移出窗口，再从最右边移进一个数据，并计算 各数据单元的斜率。设检峰状态控制变量 a = 0。2 1 起点、前拐点检测当窗口内数据满足 h0 h1 hn 1 时，判断此时为峰的起 点，a = 1，继续移动窗口。如图 1 所示，当谱图的一阶导数第一 次到达顶点时，为峰的前拐点，即窗口满足 d0 d1 dn 1 时， 判断为峰前拐点，a = 2。难完全满足，因此需对算法条件进行宽松变换7。对 h h 0 1做如下变换: hn 1定义一上升趋势变量并初始化 pk_up = 0，初始化 h= h ，max0b = 0，从 0 到 n 1 循环比较 hi 和 hmax ，若 hi hmax ，则 hmax = hi ，b = b + 1，pk_up = pk_up + b，若 hi hmax ，则 b = 0，最后，比较 pk_ up 与设定的阈值 PT，如 pk _up = PT，则判断谱图曲线具有上 升趋势。类似地，对 h0 h1 hn 1 做如下变换: 定义一下降趋势变量并初始化 pk _down = 0，初始化 hmin =h0，b = 0，从 0 到 n 1 循环比较 hi 和 hmin ，若 hi hmin ，则 hmin = hi ，b = b + 1，pk_down = pk_down b，若 hi hmin ，b = 0。最后，比 较 pk_down 与设定的阈值 PT，如 pk_down = PT，则判断谱图 曲线具有下降趋势。上述算法保证了当窗口内数据较长地连续上升或下降时所 得的总体趋势评价参数绝对值，要比当窗口内数据杂乱无章没 有趋势时所得的总体趋势评价参数绝对值大得多，这样就压低 了噪声的影响，突出峰信号，从而具有排除噪声干扰，提高分辨 窗口内数据总体趋势的能力。由于总体趋势评价参数拉大了窗 口内数据有无趋势之间的差别，使得判断窗口内数据总体趋势 的阈值有较强的健壮性，可在一较为宽松范围内选取，设得稍高 或稍低不影响判定结果，加之化学谱图数据处理中已对高频噪 声进行了滤除，故判断窗口内数据总体趋势的阈值可以固定为 某一合适的值。由此可见，将上述重叠峰检测算法的条件经宽 松变换后，可以进一步提高化学谱图数据处理的自动化水平，使 重叠峰检测受到人为的影响较少。3 2 程序实现根据以上算法及条件变换，除了定义窗口曲线上升趋势参 数 pk_up，下降趋势参数 pk_down，还定义定义窗口一阶倒数上 升趋势参数 pk_dup，下降趋势参数 pk_ddown。这四个参数分别 描述窗口曲线和一阶倒数的上升下降趋势。图 2 是画出的程序 框图。2 2前肩峰、顶点、后拐点检测当峰检测到前拐点即 a = 1 后，可能出现两种情况:1) 当窗口中谱图曲线由上升变为下降，即 h0 h1 hn 1 ，则判断此时为最高点，a = 3。2) 当窗口中谱图曲线一阶导数 d0 d1 dn 1 ，则判断 为前肩峰，a = 2，准备找下一个峰顶。当检测到峰顶后，继续移动窗口，当窗口内谱曲数据的一阶 导数满足 d0 d1 dn 1 ，即一阶导数开始上升，判断为峰的 后拐点，a = 4。2 3 重叠峰谷点和后肩峰检测当峰检测到右拐点即 a = 4 后，可能出现下面三种情况，其 中第一种情况为谷点型重叠峰如图 1(a)所示，第二种情况为肩 峰型重叠峰如图 1(b)所示，第三种情况为独立峰。1) 如果谱图曲线已不再下降，而一阶导数呈上升趋势，即 窗口不满足条件 h0 h1 hn 1 但满足条件 d0 d1 dn 1 ， 则判断为重叠峰谷点，a = 1。2) 如果谱图曲线下降且一阶导数也下降，即同时满足条件h0 h1 hn 1 和条件 d0 d1 dn 1 ，则判断为肩峰，a = 4。3) 如果谱图曲线不在下降且一阶导数不再上升，即窗口同陈淑珍等:一种色谱谱峰识别算法的实现第 11 期311个大小为 n 的数据窗口用来读取色谱峰数据，计算其 pk_up 值，记录下后移动窗口，当窗口移到终点时循环比较 pk_up 值，取得 最大的 pk_up 值，作为阈值 PT 的参考值。阈值 dPT 的确定方法 与阈值 PT 的确定方法类似。PT 如果取得太大，会使得峰起点 推后，终点提前，并引起扁平峰漏判。在阈值越小，其精确度越 大5，但峰值若取得太小，容易引起起点误判，并找不到终点， 导致一些峰检测不出。在白酒标样谱图检测中，选取数据窗口 为 20，当阈值 PT 设为 50 以上时最后一个扁平峰即第 9 峰会漏 判，当阈值设为 30 以下时增加很多干扰峰，并且会由于有的峰 找不到终点程序不往下继续运行，没有谱峰结果。谱图显示4图 3 是程序根据上述检测方法对白酒标样进行检峰的示例，其中第 4 峰为拖尾峰，第 5 峰为骑峰，第 6、7、8 峰为重叠峰 群。算法在峰起点顶点右拐点终点都打上一条短竖线，其中第5 峰和第 6 峰的两条竖杠分别代表前一峰的终点和后一峰的起 点，第 6 峰第 7 峰和第 7 峰第 8 峰间是重叠峰谷点，只有一条竖 杠，算法能够自动判别，并对重叠峰群基线自动校准(6 峰起点 到 8 峰终点做了一条校准基线)，对骑峰进行判定并做了切线 平滑处理8(4 峰谷点到 5 峰切点的一条切线)，该谱图证实了 在数据预处理和程序法确定合适阈值的基础上，一阶导数和一 阶导数趋势变形法结合起来的算法能够对重叠峰谷点和峰终点 的区分。图 2 程序框图3 3阈值的确定在程序中需要确定两个参数，谱图曲线趋势阈值 PT 和一阶导数趋势阈值 dPT。虽然本文采用的算法使阈值的确定有一定的弹性，即可以在一定范围内变动而不影响检峰结果，但合适 阈值的选取对峰的确定很关键，用程序法处理未输入色谱峰数 据的空载情况下根据仪器噪声和基线信号能够取得合适的阈 值:空载情况下取一段色谱峰数据，数据个数一般大于 200，取一图 3 白酒标样谱图表 1 是白酒标样经过峰检测和基线校准后计算的详细 信息。表 1 白酒标样详细信息基于文中算法的重叠峰的自动识别和基线的自动校准可以很大程度提高定性算法的准确度，表 2 中对比了常用一阶导数 法，谱图曲线和二阶导数结合的算法及本文采用的改进的谱图 曲线和一阶导数趋势法对白酒标样的处理，三种方法的计算结 果如表所示。组分 4 乙醇含量最大，其绝对误差最大，其他几个组分由于组分含量少，因此绝对误差小而相对误差大。其中一阶导数不能自动识别重叠峰，第 4 第 5 分峰只能从谷点垂线切 割，两峰重叠重叠程度越大，垂线切割处理误差越大，特别是大 峰后面的的骑峰，垂线切割导致大峰含量会比实际含量减少很 (下转第 323 页)编号起始(min)结束(min)停留时间(min)峰高峰面积比例拐点10 3150 4500 3570 3402750 000122920 4520 5850 5251 20811 5810 002231630 5870 6580 6170 357850 000037140 6601 5780 87368 9715 168 2210 978652851 5802 3501 7730 83819 4980 00371 08262 4052 8772 6530 63017 4830 00331 63172 8783 4603 1230 63821 8620 00411 92383 4624 3573 7460 56820 7340 00392 30996 5807 9347 1550 44821 6570 00414 514崔得龙等:组合 NMF 和 SVD 的鲁棒图像零水印算法第 11 期323学报，2007，12(4):581-5853 崔得龙，左敬龙，彭志平 一种鲁棒音频零水印算法 J 华中师范 大学学报:自然科学版，2011，45(1):48-534 叶天语 抗 JPEG 压缩和几何攻击的鲁棒零水印算法J 光子学 报，2012，41(2):210-2175 叶天语 基于方差的奇异值分解域鲁棒零水印算法J 光子学 报，2011，40(6):961-9666 宋伟，侯建军，李赵红，等 一种基于 Logistic 混沌系统和奇异值分 解的零水印算法J 物理学报，2009，58(7):4449-44567 Lee D D，Seung H S Learning The Parts Of Objects By NonnegativeMatrix FactorizationJ Nature，1999，401(6755): 788-7918 李乐，章毓晋 非负矩阵分解算法综述J 电子学报，2008，36(4):737-7439 Ha Ja A，Mohammad A Digital audio watermarking based on the dis- crete wavelet transform and singular value decompositionJ European Journal of Scientific esearch，2010，39(1): 6-2110 叶天语，马兆丰，钮心忻，等 强鲁棒零水印技术J 北京邮电大 学学报，2010，33(3):126-129(上接第 311 页)多，而骑峰含量会比实际值大很多。6、7、8 重叠峰群处理成单 峰，也不会有有后面的基线校准，如图 3 所示，其计算面积的基 线都是起始点向后的水平线，而不是重叠峰群的起点到重叠峰 群的终点的直线即校准基线，因为把校准基线以下水平基线以 上的面积算进来，计算面积会比实际含量大。谱图曲线和二阶 导数结合的算法虽然也能识别重叠峰，但由于二阶导数是对一 阶导数的导数，其对谱图曲线的变化比一阶导数敏感得多，即对 噪音敏感，其阈值设定弹性小，适合一个峰型的阈值不一定适合 其他峰型，因此总体上其计算准确度虽比一阶导数有所提高，但 比本文算法低。本文的算法由于对谱图曲线和一阶导数曲线的 描述改进过，使得分辨率提高，计算精度较前两种方法都大为 提高。图 4 Lena 图像零水印与均匀分布随机序列间的相似度结语5为了解决现有数字水印中鲁棒性和不可感知性之间的矛盾，设计了一种组合非负矩阵分解和奇异值分解的鲁棒图像零 水印算法。算法主要特点是:(1) 利用 NMF 获取分块图像的总 体特征;(2) 利用 SVD 进行获取局部特征和信息压缩。实验结 果表明该方案在有效降低算法复杂度的同时对常见信号处理具 有很强的鲁棒性。此外，本方案还具有计算简单、容易实现、不 破坏原始图像载体信息完整性等优点，大大增强了其用于数字 图像作品版权保护的实用性，具有一定的应用价值。今后的工 作将集中在进一步探索新的特征提取算子、降低算法复杂度等 方面展开。参 考 文 献1 温泉，孙锬锋，王树勋 零水印的概念与应用J 电子学报，2003，31(2):214-2162 马建湖，何甲兴 基于小波变换的零水印算法J 中国图象图形表 2 重叠峰组分含量与分析结果比较通过表 2 分析，一阶导数由于没对重叠峰进行处理误差较大，谱图曲线与二阶导数方法由于能识别重叠峰误差比一阶导 数小。本文算法由于其一采用了结合算法，其二对结合算法做 了改进在，因此数据处理时分析结果的误差最小。实验结果说 明本文的重叠峰识别和基线校准的方法可行并且准确度高。材料与化学工程学院，20042 Li B Q，Siu S，Evans J W M