(精品论文)神经元动作电位检测分类生物医学毕业论文

多电极记录神经元动作电位的检测与分类 多电极记录神经元动作电位的检测与分类 摘要 研究神经系统群体特征，既需要得到多个神经元同一时间的信息，又需要掌握单个神 经元的放电序列。多电极细胞外记录是对神经系统进行研究的基本手段。多电极细胞外记 录的方法经过了一定时间的发展，已经得到了广泛的应用。但是对多电极记录信号的处理 一直是个难题。研究神经系统群体特征，需要掌握单个神经元的放电序列，然而细胞外记 录到的信号一方面夹杂大量背景噪音，另一方面是电极区域多个神经元放电动作的叠加。 如何将神经电信号从原始信号中准确提取出来，得知信号记录到的是多少个神经元活动的 叠加，并且将信号中的动作电位归类于单个神经元，是一切研究解码过程的基础。目前已 经有很多种方法，但是第一步动作电位检测的方法始终不尽如人意。而如果这一步的结果 不准确，后面的工作都仿佛空中楼阁。因此本文试图找到一种比较好的动作电位检测方法。 本文提出的方法是，首先利用常用的阈值检测方法对原始数据进行初步的动作电位检 测，然后利用主成分分析方法以及减法聚类获得动作电位的平均模板，以动作电位的平均 波形作为形态学滤波器，重新对原始数据进行动作电位检测。此方法使用在模拟数据中， 在各种噪声强度下效果均比阈值法有所提高，并且在存在基线漂移的情况下效果明显较好。 最后将算法用到采集的实验数据样本中。 关键词：动作电位检测，动作电位分类，阈值检测，形态学滤波器 多电极记录神经元动作电位的检测与分类 DETECTION AND SORTING OF ACTION POTENTIALS RECORDED BY MULTI-ELECTRODE SYSTEM ABSTRACT To understand the population behaviors in nervous system, we need both the real-time information each neuron carries, and the exact firing sequence of the individual neuron. Multi- electrode system is the fundamental tool for research in nervous system. With its development after years, multi-electrode system has been widely used. But the process of the signal extracted from multi-electrode system is still a big problem. We want the exact firing sequence of individual neuron, however, the signal is corrupted with a large amount of background noise and the signal may involve the firing activities of more than one neuron. To get all the spikes, count the number of neurons contributing to the signal, and find the neuron that fires each spike is the very first step of all the research. Many methods have been developed, but the first step, spike detection is still not satisfying. Therefore, in this paper, a better method of spike detection is expected. This paper proposes a method combining threshold detection and morphological filter. Firstly, apply the threshold detection to the recorded data, and cluster the spikes with principal component analysis. Secondly, take the template of one cluster of spikes as the morphological filter. Finally, filter the raw data with the best morphological filter and redetect the filtered signal. This method works better than threshold detection with varied signal noise ratio. When the baseline shifts, the new method is not influenced while the threshold detection is apparently inferior. The method is also used to process real data. Key words: spike detection, spike sorting, threshold detection, morphological filter 多电极记录神经元动作电位的检测与分类 目录 第一章 绪论1 1.1 神经元动作电位记录系统1 1.2 神经元动作电位检测的意义2 1.3 动作电位检测的方法及发展概述3 1.3.1 阈值检测法.3 1.3.2 窗口检测.4 1.3.3 基于非线性能量算子的检测.4 1.3.4 匹配滤波方法5 1.3.5 基于概率的检测5 1.3.6 基于小波变换的检测.5 1.4 神经元动作电位的分类方法6 1.4.1 模板匹配.6 1.4.2 基于特征分析的分类方法.6 1.4.3 聚类方法7 1.5 本章小结8 第二章 材料与方法10 2.1 模拟数据方法.10 2.2 神经元动作电位信号多电极记录系统11 2.2.1 多电极阵列11 2.2.2 视网膜标本.12 2.2.3 灌流系统.12 2.2.4 刺激和记录系统12 2.3 基于形态学滤波器的方法.13 2.3.1 动作电位检测方法.14 2.3.2 动作电位分类方法.15 2.4 本章小结.16 第三章 结果17 3.1 将所提出算法用于模拟数据17 3.1.1 阈值法检测与主成分分析.17 3.1.2 形态学滤波器的构建与滤波.20 3.2 两种方法对比结果22 3.3 模拟基线漂移数据结果23 3.4 真实实验数据结果25 3.5 本章小结28 第四章 总结与展望30 谢辞33 原文及译文34 多电极记录神经元动作电位的检测与分类 第 0 页 共 64 页 第一章 绪论 神经系统主导着人类的感知、思维，是体内重要的调节系统。长久以来，神经系统的 工作机理一直吸引着人类，进行不断的求知探索。然而由于神经系统的复杂性，许多问题 仍然没有答案。人们已经知道，神经元是神经系统信息处理的基本单位。神经元之间通过 丰富的突触联系构成复杂的功能性网络。对神经系统功能机制研究的一个有效途径在于对 神经元构成的信号通路以及通过突触传递所实现的细胞间通讯进行研究。神经系统内的信 息交流是以电信号和化学信号的形式进行的。动作电位是处于静息电位状态的细胞膜受到 适当刺激而产生的膜电位变化，是实现神经传导的生理基础，神经元之间往往通过一系列 动作电位信号的发放和传递来进行信息的传递、交流和处理。放电活动产生的电信号就像 电报密码一样，要想揭开神经系统的神秘面纱，就要破译这些电密码的生理意义。而准确 不失真地得到这些电信号是一切工作的第一步。神经学家不断地尝试用实验手段获得神经 电信号，对其进行分析处理，破译神经电信号的编码。 1.1 神经元动作电位记录系统 为了获得神经元的电信号信息，人们采用的方法是电极记录。 电极的发展从单电极开始。神经微电极能够测量单个神经元的电位。经历了结构和工 艺的不断完善，现在的植入微电极可以测得神经元膜上电位变化，并通过后继的电路系统 实现对这个微弱电信号的放大、滤波和变换等一系列处理。 虽然单电极记录技术日趋完善，但随着对神经系统的深入探索，单电极记录并不能满 足研究的需要。神经元信息处理中体现出的群体性越来越被人们关注。人们还发现中枢神 经系统具有非常复杂的神经网络，神经系统对于外界刺激的感受和信息处理的每个过程都 涉及多个神经元的共同作用。因此，要深入了解神经系统的工作机制，需要同时检测不同 神经元的电位变化，研究它们的相互关系。为了深入研究神经元的编码机制中的协同机制， 多电极记录技术新近得到很大发展。多电极记录使得同时记录多个神经元的同步活动成为 可能，为研究神经元之间的群体活动提供了实验基础。 单个细胞植入式微电极常采用尖端暴露的绝缘细金属丝或毛细玻璃管，但这些技术不 适合制作多个记录点的微电极。经过十多年的开发和试验，目前常用的多电极多由美国犹 他大学开发的针形硅电极阵列1和美国密歇根大学开发的线性硅电极阵列2演变而来，与 之配合的多通道记录设备也得到了良好发展。图 1-1A 所示为多电极阵列照片，B 为电极阵 列排布图以及平面电极表面放大图。 多电极记录神经元动作电位的检测与分类 第 1 页 共 64 页 A B 图图 1-1 多电极阵列实物图（多电极阵列实物图（MCS GmbH, Germany） 1.2 神经元动作电位检测的意义 多电极记录为神经科学研究提供了更好的资料，给研究带来更多更有用的信息。可同 时记录多点神经电位，还可以对各神经元之间的活动关系等进行研究。这都是以往单电极 所不能比拟的。但是，多电极采取胞外记录方式，电极可能记录到其覆盖区域周围相邻多 个神经元的电信号。同时，同一个神经元的信号，也可能被多个电极记录到。并且记录到 的信号是神经元电信号与大量噪声的加和。从多电极记录到的信号中，真实准确地还原出 每个神经元放电的信息，是对于放电编码研究的第一步。为了从夹杂噪声的信号中获得最 有用的电信号信息，对放电序列的检测甄别十分关键。 图图 1-2 多电极中一个电极记录到的放电活动多电极中一个电极记录到的放电活动 多电极记录神经元动作电位的检测与分类 第 2 页 共 64 页 图图 1-3 多电极记录数据局部放大多电极记录数据局部放大 图 1-2 中信号的采样率为 40KHz。将图 1-2 中信号的前 1000 个采样点放大，如图 1-3 所示。可以看到，一方面记录到的信号有 50Hz 工频噪声以及大量的背景噪声存在，另一 方面动作电位有多种形态，可能是记录到了不同神经元的放电活动，动作电位波形叠加产 生，也可能是动作电位被噪声干扰产生形状畸变。多电极胞外记录的背景噪声较大，来源 很多，例如：(1)来自电极本身的阻抗，(2)来自参考电极，(3)远处神经元场电位的高频分量， (4)由于电线等突然移动产生的电流3。面对这样的数据，首先需要把具有生理意义的动作 电位信息从噪声中抽取出来。这一步称为动作电位的检测。动作电位的检测是后续数据处 理过程的基础。要力求做到既不漏掉实际动作电位，又不将噪声误判为动作电位。如果第 一步动作电位检测的方法不科学，结果准确度不够高的话，之后的动作电位分类以及一切 计算及其结论都仿佛空中楼阁，难以站得住脚。 对于这个关键问题，神经科学家们提出了基于各种原理的解决办法。关于各种方法的 优劣很难一概而论。可以说，对于不同电极及数据，有不同的最优方法。 1.3 动作电位检测的方法及发展概述 1.3.1 阈值检测法 动作电位检测最初和最简单的方法是阈值检测法。对大多数动作电位来说，其波形最 明显的特征就是幅值。可以预先对数据本身设定阈值，判定动作电位波形的出现。 常用算法的原理是，先计算出信号的标准差，一般取三到五倍的标准差作为阈值，对 每一个数据点挨个检验，如果大于等于阈值，则认为检测到动作电位。 阈值检测的优点是原理清晰简单，实现方便易行，缺点在于阈值的估计有难度，设置 过低有可能将噪声算作动作电位，设置过高又有可能漏掉真实动作电位4。尤其不适用于 背景比较复杂的信号。信号被噪声影响严重畸变时，很可能漏检信号。这体现在，信号被 某种干扰叠加而“漂移”时，阈值检测难以正确提取动作电位。噪声功率与信号功率相当时， 也难以将信号检出。 多电极记录神经元动作电位的检测与分类 第 3 页 共 64 页 图图 1-3 有漂移情况的多电极记录信号有漂移情况的多电极记录信号4 图 1-3 是信号被某种干扰叠加而“漂移”(B、C 两处锋电位)的情况。若阈值设置较大 (一 3， 为数据标准差)，将把“漂移”到阈值以下的噪声判断成信号；若设置较小(一 4)， 将漏掉部分信号。 1.3.2 窗口检测 顾名思义，窗口检测的方法，借鉴了小波分析的思想，用一个移动窗沿采集到的神经 元电信号移动，给定各种判定条件，比如使用窗口内信号的幅值比、波宽或者是波形等特 征，若特征值与预先设定的标准相似，就认为检测到动作电位信号5。 这种方法的判定不是基于笼统的全局条件而是对于信号的每一个局部(窗口)加以验证, 因此精度高于普通的阈值设定。尤其是对于信号基线产生漂移的情况，效果明显优于全局 算法。但是每次移动窗口后重新计算判定条件，计算量较大，耗费时间长。并且窗口的大 小以及窗口内标准的设定有一定难度，一定程度上根据研究人员个人经验，因此不够准确， 难以找到足够理论依据。 1.3.3 基于非线性能量算子的检测 针对阈值检测法在信噪比较差时检测的缺陷，有学者提出基于非线性能量算子的检测 方法。非线性能量算子能够突出信号中的瞬时分量，是一种理想的波形峰值检测方法。非 线性能量算子方法利用即时频率和幅度信息区分动作电位，大大提高了输出结果的信噪比 6。通过找到非线性能量算子对原始信号的输出中的峰值所在，来决定原始信号中动作电 位的发放时刻。由于动作电位的发放通常被认为是“全”或“无”的，在研究中就不考虑动作 电位的幅度，只按动作电位发放的先后顺序进行排列，从而获得每个电极上的动作电位串。 如果检测得到的动作电位被看成完全相同的，只考虑放电序列的精确时间，从而可得到每 个神经节细胞的放电序列。图1.4显示了多电极阵列上记录得到的视网膜神经节细胞的信号 与其经过非线性处理过的结果，结果表明信号的信噪比得到提高，有利于继续使用阈值法 检测其结果（虚线表示阈值）7。 图图 1-4 (a) 电极记录系统的原始数据电极记录系统的原始数据 (b) 非线性能量算子作用后的输出非线性能量算子作用后的输出7 多电极记录神经元动作电位的检测与分类 第 4 页 共 64 页 但非线性能量算子本身存在着一些不足之处，当用非线性能量算子去处理一个比较复 杂的信号时，其结果输出中有一项是对应于交叉项的时变部分，它是由于非线性能量算子 的非线性操作引起的，其对于检测来说是有害的。另外，非线性能量算子对噪声也存在着 一定的敏感性。如果信号是在加性高斯白噪声的环境中，则其输出将会增加一个噪声项。 还有，非线性能量算子假设条件为背景信号是平稳的，将其用于非平稳的背景时，其效果 有时不是很理想。 1.3.4 匹配滤波方法 匹配滤波方法主要通过极值化信噪比来判断动作电位信号的存在。此算法与基于非线 性能量算子的检测方法类似，期望改善信噪比。视动作电位为待检信号，视非动作电位信 号为噪声，寻找输出信噪比最大的线性相位滤波器，在噪声背景中检测待检信号，寻求最 大信噪比。 匹配滤波的实质是要设计一个和信号中感兴趣信号片段或目标信号片段相似的滤波器, 经过匹配滤波可以获得希望得到的信号信息。匹配滤波器有以下特点: (1)最大信噪比只与 信号的能量和噪声强度成正比，而与信号的波形无关；(2)对信号幅度和时延具有适应性,但 对频移不具有适应性8。 由于匹配滤波器输出信噪比只与输入信号的能量和白噪声的功率谱密度有关, 而与输 入信号的形状和噪声的分布无关，目标信号的能量越大越好。当动作电位波形的能量远远 大于噪声能量，动作电位的主频段相对稳定时，可以采用匹配滤波。噪声滤波器对信号的 各频率分量起到幅度同相则相加的作用,对噪声的各分量起到功率相加的作用,综合而言,信 噪比得到提高。当动作电位的形态和宽度一致性好时,此方法是最佳线性滤波器。但当信号 成分复杂，动作电位波形态变化显著,目标波形不明确时，检测性能不可避免会下降。当信 号能量与噪声能量相仿时，检测性能也受到影响。 1.3.5 基于概率的检测 基于概率的检测方法中， “一个信号片段是动作电位”这个事件的概率，由该段信号的形 状、重复出现的次数等决定9。因此通过研究各信号片段的统计规律即可找出源信号中有 用的动作电位。 在原始数据中将可能的动作电位波形检出，作为一个信号片段，方法是阈值检测噪声 峰值，取阈值以上峰值及其左右两边一定数量的数据点。 “一个信号片段是动作电位”的后验 概率，是由它的形状和重复性决定的。根据这两个条件的矩阵方程，由一定的概率运算， 可以得出一个信号片段是动作电位的概率。 这种方法最大的特点是，省去了对待测信号中动作电位形态等参数的估计，不需使用 模板、滤波器等。因此克服了这些方法检测动作电位的缺陷和局限性。 对于基于概率的检测，其较大的问题是先验概率如何获得，以及判定事件为真的概率 阈值如何设定。 1.3.6 基于小波变换的检测 小波变换是近年来得到迅猛发展的一种数学方法，其作为一种信号的时间尺度分析方 法，具有多分辨率分析的特点，从而可以在时、频两域获得表征信号局部特征的能力。因 为其可变的时间窗和频率窗，使得它对于信号具有很高的适应性。 小波变换的一个重要应用就是用于信号特征点的检测。一个偶对称小波会把信号峰值 变成极大值，一个奇对称小波则会把信号峰值变成一个过零点，因此小波变换能够用于信 号特征点的检测。通过把信号进行小波变换后，对波形进行假设检验，从而检测到动作电 位信号的存在10。这种方法可以不需要模板的构造或阈值的设定，在数据条件较差时效果 仍很好。其算法较简单，可以做到实时检测。 基于小波变换的检测也有多种方法。其基本步骤为： 多电极记录神经元动作电位的检测与分类 第 5 页 共 64 页 （1）使用适当的小波基对信号进行多尺度分解。 （2）在每个尺度上分离信号和噪声。 （3）在（1）和（2）的基础上，在每个尺度上进行贝叶斯假设检验，找出动作电位。 （4）将各个尺度上的结果结合。 （5）估计动作电位出现的时间。 如果得到噪声和信号的先验信息，这种方法的效果更好。用 FIR 滤波器组来实现的小 波变换兼顾处理时间，抗干扰能力强，准确率高。但信号的干扰较大时，还有误检现象。 1.4 神经元动作电位的分类方法 将神经元信号从背景噪声中检测出来只是信号处理的初始的一步，多电极记录的一个 重要问题是单个电极往往会记录到多个神经元的活动。因此如何将检测出来的动作电位归 类于单个神经元也是一个有趣的复杂问题。 动作电位分类方法分为模板匹配法和特征聚类法两大类。 1.4.1 模板匹配 模板匹配方法顾名思义，是事先选取典型的动作电位波形作为模板，用检测到的可能 的动作电位与其进行匹配，并采用一定的评估办法根据其与某一模板的匹配程度判断是否 归于一类，再将每种模板匹配失败的事件用多种模板移动窗口去匹配，得到他们的叠加成 分。 这种方法的关键步骤，一是模板的设定，二是评价实际动作电位和模板匹配程度的分 类标准。最初的模板根据实验或者经验得到，现在多采用事先聚类，找出一类动作电位的 平均波形，将其设为模板的方法。模板匹配采用 2-test，贝叶斯分类器等方法，通过计算 动作电位信号与模板之间的距离来确定是否属于该模板11。 这种方法的缺点是第一步模板的设置有难度，需要科学依据。模板匹配的标准也因数 据而异。 目前这种方法的最新改进与完善，就在于对模板匹配筛选标准以及匹配过程的优化， 以使算法更高效。 1.4.2 基于特征分析的分类方法 由于神经元的动作电位具有“全”或“无”的特征，它可能因刺激过弱而不出现，但在刺 激到达阈值以后，就始终保持一定的大小和波形。不同神经元发放的动作电位，其特征有 显著不同。这些不同体现在很多数据量上，如幅值、波宽、发放频率等。早期的分类方法 简单的使用一些显著的特征量，如幅值、波宽、发放频率，作为分类依据。这些方法比较 直观和方便，但依据特征量的选择不可避免带来主观因素，而且容易漏掉很多有用的但比 较隐蔽的信息。因此，现在多采用主元分析方法(Principal Components Analysis, PCA)进行 动作电位分类。PCA过程实质上是把每个动作电位的各个特征量排列成一组n维空间向量， 对原坐标系进行线性重组，使得新坐标系的第一坐标轴（第一主成分）对应于数据向量变 异最大的方向，第二坐标轴（第二主成分）对应于数据变异次最大的方向，依此类推，实 现将P维向量降维到M（PM）维空间。 多电极记录神经元动作电位的检测与分类 第 6 页 共 64 页 点簇2 x Y 主成分次成分 点簇1 图图 1-5 主成分分析法降维示意图主成分分析法降维示意图 如图 1-5，原坐标空间中我们无法单独通过各个数据点的 x 或 y 坐标来区分两组数据， 但经过变换后，我们仅依据重组产生的坐标系中的主成分分量，即可轻松区分两组数据点。 主成分分析法将动作电位投射到多维空间，产生新的隐式变量，再将多个变量（指标） 化为少数几个主成分，而这几个主成分可以反映原来多个变量的大部分信息，即用最少的 变量(主成分)大致反映出所研究的波形特征。主元分析的动作电位信号的分类方法能够节 省计算资源，极大简化分类所需要的运算量，同时提高分类结果的可靠性，易于实现。并 且比一般的特征分类方法的精度高。现在的特征提取主要是利于聚类分析，提高分类的效 果。 1.4.3 聚类方法 得到特征分析的散点图后，往往要进行聚类。 聚类问题可以定义如下：给定 d 维空间的 n 个数据点，把这 n 个点分成 k 个组，即满 足最大的组内相似性和最小的组间相似性，使得不同聚类中的数据尽可能地不同，而同一 聚类中的数据尽可能地相似。 K 均值聚类 K 均值聚类方法是一种经常使用的聚类方法。它基本的思想是两类估计：首先需要估 计每个类的中心位置，然后根据元素并入哪个类来分别估计所有元素的划分结果，并且使 用其中一个估计的结果来调整另一个估计12。 假定现有一组包含 M 个类的数据，这 M 个类的总和平方差为 J。K 均值聚类的目的是 通过重复再分配类成员使各聚类本身尽可能紧凑，而各聚类之间尽可能分开。也就是要选 取 M 个类以及相关的类中心，使得 J 的值为最小。 具体思路如下： （1）将每个动作电位向量随机划分到 M 个类中任意一个； （2）计算 M 个类中每个类的均值； （3）计算每个动作电位向量到每个类均值的距离； （4）将这个动作电位向量并入最近的类中； （5）重新计算这个类的均值，作为类的中心； 重复计算步骤（1）-（5） ，直到各个动作电位向量对其类中心的距离总和最小，从而 确定各个类及其中心。 其缺陷在于运算量大，该算法本质上属于局部搜索的爬山法，容易陷入局部极小点。 最大的问题是对聚类中心的初值选择非常敏感，这样往往会因为初始值的位置或个数选择 多电极记录神经元动作电位的检测与分类 第 7 页 共 64 页 的不好减低算法性能。 模糊 C 均值聚类 经典的聚类算法是将每一个辨识对象严格地划分为属于某一类，但在实际应用过程中 某些对象并不具有严格的属性，它们可能位于两类之间，这时采用模糊聚类可以获得更好 的效果。模糊C均值聚类（Fuzzy clustering means, FCM）是模糊目标函数法中非常有效的 一种，是目前广泛采用的一种聚类算法13。与K均值聚类的硬分类方法不同，模糊C均值聚 类并非简单将聚类数据硬性划分至某个具体类，而是采用隶属度矩阵刻画各数据的划分情 况。算法使用了最小化整个权重的均方差的思想，允许每个特征向量以一个模糊值(0和1之 间)属于每个簇。它能给出每个样本隶属于某个聚类的隶属度，即使对于很难明显分类的变 量，模糊C均值聚类也能得到较为满意的效果。 模糊 C 均值聚类方法也属于局部寻优技术，其聚类质量依赖于初始聚类中心(即初始模 糊聚类矩阵)的选取。即初始聚类中心将决定模糊 C 均值聚类算法的聚类收敛速度和聚类精 度。而且模糊 C 均值聚类算法采用距离的倒数来计算模糊权重。当一个特征向量到两个簇 中心的距离相等时，它在这两个簇上的权重相同，而不管簇是如何分布的。因此，它不能 够区分特征向量分布不同的两簇。所以，模糊 C 均值聚类算法更加适合于那些在簇中心周 围或多或少呈均匀分布的数据。模糊 C 均值聚类算法会将靠近边界的具有固有形状的两个 簇合并成为一个大的簇。因此很难将非常接近的类聚类到一起。 基于神经网络的聚类方法 人工神经网络是由简单处理元构成的规模宏大的并行分布式非线性处理器。天然具有 存储经验知识和使之可用的特性。人工神经网络和人脑相似，可以通过对外界环境的学习 得到知识，并存储在突触权值里。神经网络的计算能力很明显有以下两点：（1）大规模并 行分布式结构。 （2）神经网络的学习能力以及由此而来的泛化能力。近年来，人工神经网 络日渐成为一种重要的分类工具。因此，非常适合于动作电位分类这样的复杂问题。 径向基函数网络（Radial Basis Function Network, RBF）是一种特殊的三层前馈网络， 它通过非线性基函数的线性组合实现映射关系。径向基函数网络一般包括输入层、隐含层 和输出层。从输入层到隐含层的变换是非线性的，而隐含层到输出层的变换是线性的。隐 含层采用径向基函数作为传输函数，通常采用高斯函数。径向基函数具有网络结构简单、 非线性逼近能力强、收敛速度快以及全局收敛等优点。径向基函数网络用于胞外记录信号 分类时，当信噪比降低，其分类正确率降低比较平缓14。但是，径向基函数网络仅仅只能 将严重重叠变形的波形单独列出来，仍无法完全解决重叠问题。 自组织特征映射(Self-Organizing Feature Map, SOFM)神经网络是较为广泛应用于聚类 的神经网络。它是一种无监督学习的神经元网络模型，主要功能是将输入的 m 维空间数据 映射到一个较低的维度，通常是一维或者二维输出，同时保持数据原有的拓扑逻辑关系15。 但该方法也存在一些问题，比如对聚类中心初值敏感，只能实现粗略聚类、对非球形族的 聚类效果不理想等。单个 SOFM 聚类对输入参数，尤其对输出层神经元个数非常敏感。如 果个数太小，而真实聚类数目较多，则无法实现很好的聚类；反之如果个数太多，而真实 聚类数目较少，聚成很多小聚类，也无法达到理想的聚类效果。而真实聚类数目又很难估 计。 1.5 本章小结 动作电位携带着神经元相互交流，协同工作的信息。对神经电信号的分析与研究，如 果第一步的检测不正确，研究对象的确定有误，那研究中得到的结论的可信度便无从谈起。 本章阐述了多电极记录神经元发放动作电位的意义，以及动作电位的检测及分类的意 多电极记录神经元动作电位的检测与分类 第 8 页 共 64 页 义。并分别总结了动作电位检测的常用方法，阈值检测法、窗口检测法、非线性能量算子、 匹配滤波等；以及分类的常用方法，模板匹配法和特征聚类法，并分析了这些方法的优劣 和使用条件。 可以看到，基于动作电位检测和分类对于神经方面研究的重要性，国内外的学者在提 高神经元检测分类正确率上做了大量的工作。但目前没有一种方法能够达到百分百的准确 率将电极记录到的数据中的动作电位事件毫无差错的检出，也没有一种接近完美的动作电 位分类方法。在动作电位的分类方面，目前有一些非常复杂智能的方法，比如神经网络的 运用等，能够通过尽可能完备的方法提高分类的正确率。然而，动作电位的检测仍是一大 难题。 多电极记录神经元动作电位的检测与分类 第 9 页 共 64 页 第二章 材料与方法 上一章中概述了动作电位检测与分类的一些常用方法。本章将介绍模拟数据的方法， 得到本文中真实数据的实验方法，并将着重分析本文立足于改进结合的原有动作电位检测 方法，详细说明本文提出的新的动作电位检测方法，以及本文使用的动作电位分类方法。 2.1 模拟数据方法 在使用真实实验数据前，为了能够测试所提出的方法，控制噪声强弱、动作电位宽度 等参数，与原本的情况进行对照，对结果优劣进行准确评价，首先将算法用在模拟数据上。 分别采用了描点法和指数函数两种方法模拟单个神经元动作电位。 图图 2-1 描点法模拟两个单个神经元动作电位波形描点法模拟两个单个神经元动作电位波形 图图 2-2 指数函数模拟指数函数模拟 A,B 两个独立神经元动作电位波形两个独立神经元动作电位波形 模拟中发现指数函数效果明显较好，因此本文以下模拟数据部分中，均采用指数函数 模拟单个神经元动作电位。 多电极记录神经元动作电位的检测与分类 第 10 页 共 64 页 模拟时假设：（1）每个神经元的放电频率不变；（2）每个神经元发放的动作电位波 形不变；（3）不同神经元动作电位波形有明显差异；（4）噪声为白噪声，服从高斯分布； （5）神经电信号与噪声线性叠加。 因此采用不同参数的泊松分布，再以不同信噪比叠加高斯噪声，模拟单一电极采集到 两个神经元放电的情况。将以上两个单个动作电位波形分别以 200、120 为参数泊松分布 40、60 个，并叠加。叠加后部分动作电位重叠，因此叠加后的动作电位事件数目为 81 个。 最后再添加不同信噪比的高斯白噪声，以测试本文提出的算法。下图为信噪比为 3 时的模 拟数据。 图图 2-3 模拟数据模拟数据 2.2 神经元动作电位信号多电极记录系统 对模拟数据使用算法后，还应当将算法用于真实数据。毕竟算法的意义在于处理真实 数据。本文中的真实数据来自以下介绍的多电极记录系统，由实验室学姐景玮采集。 2.2.1 多电极阵列 实验所用的多电极阵列（Multi-electrode array，简称 MEA，MMEP-4, CNNS UNT, USA） ，尺寸为 55cm2，底板为石英玻璃材质，内嵌有 64 条彼此绝缘的导线（氧化铟锡， ITO） ，64 个电极直径大约 8m，间距为 150 m（中心到中心） ，按照 88 的阵列排布。 多电极记录神经元动作电位的检测与分类 第 11 页 共 64 页 AB 图图 2-4 多电极阵列示意图多电极阵列示意图 2.2.2 视网膜标本 实验选用牛蛙为标本，实验前暗适应30-40分钟，然后在暗红光下，迅速摘出眼球，在 角膜和巩膜交界处切开眼球，并将角膜和晶状体与后半部分分开。在灌流液中将眼球后壁 剪成44 mm2的小块，此后小心分离出视网膜。最后，将其移至多电极阵列上 （MEA，MMEP-4, CNNS UNT, USA） （图2-4） ，使神经节细胞与电极紧密接触。 2.2.3 灌流系统 实验中为达到良好的灌流效果，在紧密贴合于电极阵列上的玻璃灌流小室内，放置了 一个有机玻璃灌流槽。灌流液在蠕动泵（4道） （Ismatec SA，USA）的作用下从灌流槽的 一侧小孔流出，经过视网膜表面，溢出时被蠕动泵抽回（图2-5） 。整个灌流系统能保证灌 流液平稳地流经视网膜表面，尽可能减少对细胞反应记录的影响，灌流区域与抽出管道容 积共约0.52毫升，管流速度0.8毫升/分，因此，灌流区域的的溶液可以在39秒内更换。灌流 液的切换由多通道逻辑开关进行控制。标准灌流液（Ringer液）组成（单位：毫摩尔/升）： 100 NaCl，2.5 KCl，1.6 MgCl2，2.0 CaCl2，25 NaHCO3，10 葡萄糖，并通过含95%O2 和 5%CO2 的混合气体，同时维持PH值7.50.2。在没有持续强光照射的情况下，这样的灌流 条件能使视网膜维持良好状态达8小时之久。 图图 2-5 灌流模式图灌流模式图 多电极记录神经元动作电位的检测与分类 第 12 页 共 64 页 2.2.4 刺激和记录系统 视网膜神经节细胞的对光反应通过 64 道电极同时记录，经过一个 64 通路的放大器， 由商业软件 MEA workstation（Plexon Inc. Texas, USA）对神经元反应与其刺激同步采样， 采样信息一方面显示于屏幕，并同时将数据同步存储于计算机，以备离线分析，采样率通 常设为 40kHz。 实验中给光刺激由计算机屏幕产生，通过光学聚焦系统，投射于离体视网膜上形成 11 cm2的光斑，通过计算机编程可以设定特定的刺激模式。 此实验采集数据的整体装置如图 2-6 所示。 图图 2-6 实验装置图实验装置图 2.3 基于形态学滤波器的方法 本设计提出的新方法在于将阈值检测、主成分分析法与形态学滤波器相结合，希望实 现动作电位的较准确检测。 其步骤如下： （1）按照常用的阈值设置方法设置阈值，在原始数据中将可能的动作电位波形检出， 作为一次事件。 （2）对所有的事件采用主成分分析法，并使用减法聚类，将归为一类的所有波形平均 后得到此类的模板，作为形态学滤波器原型。 （3）找到最佳的形态学滤波器结构元素长度，然后将形态学滤波器用于原始数据，滤 波后重新检测动作电位。 这种方法理论上有两个方面的优势。第一，弥补由于阈值设定依据不够而产生的误差。 多电极记录神经元动作电位的检测与分类 第 13 页 共 64 页 比如对动作电位的遗漏或者将噪声误判为动作电位的情况，可以在这次滤波中得到改善。 第二，对于形态学滤波器的设计比较完善，接近动作电位波形。 2.3.1 动作电位检测方法 初步检测 使用阈值检测和主成分分析法结合，首先依据常规方法设置出阈值，检测信号中的峰 值，取峰值前后若干个点，作为一次事件，在原始数据中将可能的动作电位波形检出。对 所有事件采用主成分分析法，得到二维散点图，并对于主成分分析结果使用减法聚类（见 2.3.2） ，最后将聚为一类的动作电位叠加取平均，得到每一类动作电位的模板。 这一步的目的是进行初步的检测、分类，以得出形态学滤波器的结构元素。 形态学滤波方法 数学形态学是建立在积分几何和随机集论等严格数学理论基础上的一门密切联系实际 的科学16。数学形态学方法用于数字信号处理的基本思想，是利用一个称作结构元素的 “探针”收集待处理信号的信息，探针在信号中不断移动，即可考察信号各个部分之间的相 互关系，从而提取信号全局或局部的有用特征。所有的形态学处理都是基于填放结构元素 的概念。基本形态学变换建立在 Minkowski 和差运算的基础之上,其基本运算包括腐蚀、膨 胀以及由此引出的形态开运算和形态闭运算17。根据分析信号的不同，可以分为二值形态 变换和灰值形态变换，本设计中用到的是一维离散信号的灰值形态变换。 形态学滤波器是一种非线性滤波器，它基于信号的几何特征，利用预先定义的结构元 素对信号进行匹配，提取出与其相似的信号特征。它利用形态学变换算法，可以将含有复 杂成分的信号分解为具有物理意义的各个部分，使信号与背景分离并保持其全局或局部的 主要形态特征。 与传统的数字滤波器相比，基于数学形态学的消噪滤波方法具有算法简便易行、物理 意义明确、实用有效等优点。 设待处理信号 f(n)是采样得到的一维多值信号，其定义域为 Df=0，1, 2, 3, , N； g(x)为一维结构元素序列，其定义域为 Dg=0，1, 2, 3, , P；其中 P 和 N 为整数。则腐 蚀与膨胀运算分别定义为 腐蚀： （2-1） 1,2,3, (n) = min s nm fgf nxg x 膨胀： （2-2） 1,2,3, (n) = max s nm fgf nxg x 开运算： （2-3） (n) = s fgfggn 闭运算： （2-4） (n) = s fgfggn 形态开、闭运算对信号处理的效果不同:形态开可以平滑信号中的正向脉冲(峰)，除去 毛刺及小桥结构；形态闭可以平滑信号中的负向脉冲(谷)，填平小沟结构。同时,形态滤波 效果还和结构元素的尺寸和形状有关,应根据信号的特点选取不同的结构元素。使用开、闭 运算的级联组合方式,构造出形态开闭和形态闭开滤波器18，19,可以同时去除信号中的正、 负向脉冲。 基于形态学开、闭运算可以构建 3 种滤波算法：交替滤波器、混合滤波器、交替混合 多电极记录神经元动作电位的检测与分类 第 14 页 共 64 页 滤波器。 交替滤波器即开闭和闭开滤波器。 交替滤波器: 开闭运算： （2-5） OC= f nf ng ng n 闭开运算： （2-6） CO= f nf ng ng n 混和滤波器 （2-7） HF= /2f nf ng nf ng n 交替混和滤波器 （2-8） 1 = OC f nCO 2 aw Mf nf n 对于相同宽度的结构元素而言（即结构元素的长度为 P） ，交替滤波器与混和滤波器的 计算速度相近，而交替混和滤波器的速度最慢，当 P 值很大时，这种时间差异会更大。 通常交替滤波器的输出会向下（或向上）偏移，这是因为开（闭）运算的输出信号总 是位于原信号下（上）方的缘故。为解决信号输出的单向偏移问题，一般选用混合滤波器 或交替混合滤波器。 本文方法中采用级联开-闭(OC)、闭-开(CO)平均组合的交替混合滤波器。 形态学处理方法是基于结构元素填充探测的思想，结构元素的形状是影响滤波性能的 关键因素，结构元素相当于一般数字滤波器中的滤波窗口，它的选取直接影响着滤波器的 滤波性能。 由于形态滤波器的性能受结构元素的影响较大，滤波结果对结构元素的宽度和幅值比 较敏感，所以应根据滤波后要保持的信号特征选取结构元素的尺寸。 往常检测动作电位的方法中常用的结构元素有三角波，抛物线等，以接近动作电位的 形态特征20。这些方法的缺陷在于，滤波器的选取没有足够的依据，不能在物理意义上足 够接近待检测动作电位。 因此，在本文使用的新方法中，将阈值检测、主成分分析、聚类得到的动作电位模板 作为形态学滤波器的结构元素。使得结构元素在最大程度上接近待测信号片断。寻找最佳 滤波器长度的方法是，将平均动作电位的峰值固定作为结构元素的中间点，对平均动作电 位峰值两边取不同的长度作为滤波器，寻找最佳滤波效果。 2.3.2 动作电位分类方法 主成分分析方法 高维数据在使用计算机技术实现自动聚类分析的时候，首先要进行降维处理，得到低 维数据，再以此为基础进行聚类分析。 主成分分析（PCA）是一种实用的多元统计决策方法，它能够对高维多变量数据系统 进行最佳综合简化，在力保数据信息损失最少的原则下，对高维变量空间进行降维处理。 其目的是在数据空间中找出一组向量来尽可能地解释数据的方差，用较少数量的特征对样 本进行描述来降低特征空间维数，将数据从原来 n 维降低到 m 维(m2 时，阈值法与形态学滤波器 结合的方法与单纯使用阈值法相比有持续较好的表现。而当信噪比）以毫秒为单）以毫秒为单 位。位。ISI 密度是伽马密度的神经元在第密度是伽马密度的神经元在第 2 栏。生成两种动作电位的神经元（第栏。生成两种动作电位的神经元（第 3 栏）的两栏）的两 个对数正态的参数值已给出。截断的对数正态密度的神经元位于第个对数正态的参数值已给出。截断的对数正态密度的神经元位于第 6 栏。最后一行表明了栏。最后一行表明了 每个神经元的事件数。神经元与颜色的对应关系如图每个神经元的事件数。神经元与颜色的对应关系如图 2b 和和 5A：1，蓝色；，蓝色；2，绿色；，绿色；3，红，红 色；色；4，棕色；，棕色；5，紫色；，紫色；6，黄色；，黄色；7，黑色。，黑色。 实现细节 代码用C语言写成，有Gnu Public License的话可以在我们的网站： http:/www.biomedicale.univparis5.fr/physcerv/Spike-0-Matic.html找到。免费的软件图书馆 （http:/www-rocq.inria.fr/scilab/）可以生成输出图形和图形用户接口，和我们的例程同时使 用。采用GNU科学图书馆（GSL：/gsl）操作矢量和矩阵，生成 （伪）随机数（均一，标准，对数正态，伽玛） 。采用了GSL的Matsumoto and Nishimura（1998年）MT19937生成器。这个生成器有2199371步。因为代码需要很多指数和 对数等式，指数和对数函数被统计和存储在内存中，从而使计算时间缩短了30。代码在 运行Linux操作系统的个人笔记本电脑（奔腾IV的CPU速度为1.6 GHz的处理器，256 MB内 存）上编译和运行。使用GCC（/software/gcc/gcc.htm）汇编器。 结果结果 数据特性 我们的算法的表现将在本文中由模拟数据集说明。假设这些数据来自四电极记录。数 原文与译文 第 55 页 共 64 页 据中存在 6 个神经元，以及一个用于模拟难以被检测到的远离四个电极的“噪声”神经元。 通过生成系统偏差，模拟出的数据集已对算法具有足够的挑战性了。假设这组数据来自 15 秒的记录时间，其间有 5058 事件被检测到。对用户显示为图 2A。图 2B 显示相同的数据， 用颜色对应于不同神经元。每个神经元的模拟参数见表 1。读者可以很容易地看到，有一 簇有很少个点（图 2B 中的紫色点集） ，其他簇都有很多点（棕色，绿色，蓝色点集） 。事 实上，最小的簇包括 73 个事件，而最大的簇有 1474 个事件。一些簇拉得非常长（例如， 绿色和红色） ，其中有一个神经元，即双生成神经元，甚至生成 2 个分隔的簇（红色） 。 “噪 声”神经元的簇（黑色）比其它的更分散。两个簇常常表现出重叠（黄色黑色对和绿色 红色对） 。 我们的模拟噪声与记录到的噪声性能之间的偏差如图 3 所示。 图 3A1 显示了第二个 记录点的第二个（绿色）神经元的峰值、ISI 与理论值之间的关系。从这里，读者可以更了 解公式（2）的意义。图 3A2 列出了相应的差值（实际值-理论值） 。从这些点的分布读者 可以看到，这组神经元的点簇在 4 维空间的形态极大地扭