一种基于EMD的睡眠脑电图梭形波自动识别方法

一种基于EMD的睡眠脑电图梭形波自动识别方法 Vol . 17 No . 1 第 17 卷第 1 期北 方 工 业 大 学 学 报 2005 年 3 月J . NORTH CHINA UNIV. OF TECH. Mar . 2005 一种基于 EMD 的睡眠脑电图梭形波 自动识别方法 ))11 ,2 杨志华齐东旭 () )) 1中山大学信息科技学院 ,510275 ,广州 2北方工业大学 ,100041 ,北京 摘 要 提出了一种基于 EMD 分解的睡眠脑电图梭形波自动识别方法 . 通过对 EEG 信号 作 EMD 分解 ,得到具有很好时频特征的 Hilbert 谱 ,对 Hilbert 谱进行处理 ,可以自动识别梭 形波出现的准确时间和持续时间 . 结果显示 ,识别准确率高达 95. 6 %. EMD 分解 ; Hilbert 谱 ; 自动识别 ; 梭形波关键词 分类号TN 318. 8 ( ) 式 分 解 Empirical mode decomposition , EMD的 9 Hilbert 谱分析方法. EMD 分解是一种基于局 1 引言 部时间尺度的自适应 、高效率的分解方法 ,适用 于非 线 性 、非 平 稳 数 据 . 以 此 为 基 础 得 到 的 睡眠是一种复杂的生理过程 ,睡眠生理的 研究在神经精神科的临床诊断和疗效评价中具 Hilbert 谱具有很好的时频局部性. 本文在研究 有重要的作用. 国际上一般将睡眠的进程分为 睡眠脑电图波形特征的基础上 ,使用这种方法 ( 对睡眠脑电图信号进行了处理 ,进而提出了一 两个时期 , 即‚无眼球快速运动 No rapid eye ) ( 种在 EGG 中识别梭形波的新方法 . movement , NREM?和‚眼球快速运动 Rapid eye ) movement , REM?. NREM 时 期 又 分 为 一 、二 、 1 三 、四级 ,睡眠深度依次加深. 睡眠进入二 、三 2 EMD 分解方法及 Hil bert 谱 ( ) 阶段的最重要的特征是脑电图 EEG中出现睡 ( ) 眠梭形波 sleep2spindles,它们是频率为 12 Hz, Hilbert 变换和瞬时频率2 . 1 2 ( ) 20 Hz 的梭状波形. 传统上对于睡眠梭形波的 对于任 意 的 时 间 序 列 X t , 其 Hilbert 变 ( ) 换 Y t定义为 : 识别是由神经内科专家和睡眠专家通过目测分 ? 析进行的 . 自 20 世纪 70 年代起 ,国际上开始自 (τ) 1 X τ()( ) d 1 P Y t= 3 ,8 π ?τt - - ? 动睡眠分析的研究. 由于梭形波在 EEG 波 其中 , P 为柯西主值 . 列中持续的时间很短 ,用传统的时域分析和频 ( ) ( ) 由 X t和 Y t构成的复信号域分析很难获得满意的结果. θ( ) i t( ) ( ) ( ) ( ) ( ) Z t = X t + iy t = a t e 2 最近 , Huang 等 人 提 出 了 一 种 基 于 经 验 模 ( ) 被定义为 X t的解析信号 ,其中 收稿日期 :2004 —10 —26 ( ) ( ) 基金项目 :国家‚九七三?重点基础研究发展规划项目 G1998030607、国家自然科学基金重点资助项目 60133020、国家自然科学基金 ()资助项目 036608 第一作者简介 :杨志华 ,博士研究生. 主要研究方向 :数值分析 、计算机图形学. 1 2 . 4 Hilbert 谱 2 2 2 ( ) ( ( )) a t = Xt+ Yt ( ) 通过对 x t 作 EMD 分解 , 我们可以得到 ( ) Y t θ( ) )()( 3 t = arctan 一组 IMF. 对每一个 IMF 作 Hilbert 变换 ,并根据 ( )X t () 式 4计算每个分量的瞬时频率之后 ,我们可以 ( ) 分别是 Z t的幅值和相角 ,它们都是时间的函数. 将数据表示成如下形式( ) 定义信号 X t 的瞬时频率为其解析信号 n 相角的导数 ,即 () ()( ω ( ) ( ) )6 t exp itd t Z t= aj j ? ? j = 1 θd ()ω( ) 4 t = d t 比较 Fourier 展 开 , 显 见 , IMF 表 达 是 Fourier 展 ( ) 它刻画了信号 X t 在每一时间的局部振荡特 开的一种更一般的表达 . 征 ,是构成 Hilbert 谱的要素之一 .() 根据式 6,我们可以在一幅三维图中把振 2. 2 固有模式函数幅和瞬时频率表示成时间的函数 ,振幅的这种 () 固有模式函数 intrinsic mode function , IMF(ω) 时频分布被定义为 Hilbert 谱 ,记为 H , t. 是指满足如下两个条件的函数 : ?在整个数据 范围内 ,穿越零点的点的个数与极值点的个数 3 梭形波的识别 相等或相差为 1 ; ?在任意点处 ,由局部极大值 点确定的上包络线和与由局部极小值点确定的 出于对分解速度的考虑 ,我们首先将脑电 下包络线的均值为零 . 图数据分成若干短段 ,每段取 T = 6s ,对每一段 固有模式函数是瞬时频率处处有意义的函 数据 ,按如下步骤进行识别 : 数 ,它代表了信号内部固有的振荡模式 . () 1对数据作 EMD 分解 ; 求各 IMF 的瞬时 2. 3 经验模式分解 频率 和 瞬 时 振 幅 , 并 将 瞬 时 频 率 规 范 到 0 , ( ) 经验模式分解 EMD过程如下 : 99 Hz ;将瞬时振幅规范到 0,255 ;() ( ) ( ) 1初始化 : rt= x t, i = 1 ; 0 () ( ) 2根据 1得到的瞬时频率和瞬时振幅 ,() i 个 IMF : 2提取第 (形成一幅宽度为 1 200 数据长度为 6s , 采样频 () ( ) ( ) a初始化 : ht= rt, j = 1 ; 0 i - 1 ) 率为 200 Hz, 高度为 100 的灰度图 , 即 Hilbert ( ) b求 h的局部极大点和局部极小点 ; j - 1 谱 ; () c以局部极大点为节点作三次样条插() 3在 Hilbert 谱 中 截 取 10 Hz , 20 Hz 的 频 值形成上包络 ,类似地形成下包络 ;带 ,组成一个 11 行 1 200 列的矩阵 M , 并对 M ( ) ( ) d计算上下包络的平均值 mt ; j - 1 ( ) α 的数据按 7式作点变换 ,取 = 0 . 95 ,以增加 () ( ) ( ) ( ) e计算 ht= ht- mt; j j - 1 j - 1 数据的对比度 () ( ) f若 ht 的穿越零点的个数与极值点 j D m1 ( ) f x= {1 + 个数相等或最多相差为 1 ,且 π 2 ( α)sin 2 T 2 ? h - h ? j - 1 jx 1 ()5 SD = ?() πα( 7 ) ?sin - } ;2 t = 0 hD 2 m j - 1 () ( )M 中 各 列 的 最 大 值 , 组 成 长 度 为( ) ,则 imf t 4求 小于某个给定的数 = ht ,否 i j ( ) ( ) ( ) 1 200 的序列 C t . C t 元素的大小反应了该 则转至 b,且 j = j + 1 ; 时刻信号能量的大小 ; () ( ) ( ) ( ) 3rt = rt - imf t ; i i i - 1 () () ( ) ( ) 5由于梭形波的持续时间至少在 0 . 5s 以4若 rt 至少有两个极值点 ,转至 2, i ( ) ( ) 上 ,所以我们对 C t 中的每一点计算从该点开 并令 i = i + 1 ,否则 ,分解结束 , rt为余量.i 始的 以 后 1 0 0 点 的 能 量 平 均 , 得 到 新 的 序 列 用 EMD 分解可以将一个信号分解为若干 ( ) 个 IMF , 并且由 于 它 是 基 于 局 部 时 间 尺 度 的 , S t ,其长度为 1 100 . () ( ) 6在 S t 中按顺序检查各点的值 , 当遇 所以 , 对非线性 、非平稳信号具有同样高的分 ) 1 . 46s.波持续了约 们选阈值为 40,该点被认为是梭形波的起点 , ( ) 继续往下搜索 ,直到遇到某个 S t 小于阈值 , 则从该点开始的后面的第 100 点被认为是该梭 形波的终点 . 图 1 是一段进入第二期睡眠的脑电图 ,图 中 A 、B 两点分别出现睡眠梭形波 . 图 2 是 10 Hz ,20 Hz 频带范围内的 Hilbert 谱 . 从图中我们看 到 ,在 0,1s 和 4s,5s 附近出现了两个明显的 高能量带 ,它们正是梭形波出现的位置. 图 3 是 经过点变换以后的 Hilbert 谱 ,梭形波的高能量 ( ) 带比图 2 更加明显 . 图 4 是 S t的变化曲线 , ‘o ’表 示 梭 形 波 出 现 的 起 点 ‘, 3 ’表 示 终 点 , 图 1 一段进入第 2 阶段睡眠期的 EEG 信号 ( ) ‘ ?’表示 搜 索 到 的 S t 的 值 小 于 阈 值 的 位 在 A 、B 两处分别出现梭形波置 ,第一个梭形波持续了约 1 . 17s ,第 2 个梭形 图 2 图 1 中 EEG 信号在 10Hz,20Hz 频率范围内的 Hilbert 谱 图 3 对图 2 中的 Hilbert 谱作点运算得到的结果 , 在 0,1s , 4s,5s 附近出现两个明显的高能量带.高能量带更加显现 ( ) 图 4 S t 的变化曲线 准 ,识别的结果如表 1 所示. 结果显示 ,此方法 识别梭形波的准确率达 95 . 6 %. EEG 波形是典型的非平稳信号 ,因此 ,使用 4 结论 平稳信号处理技术难以获得令人满意的结果. (EMD 分解的自适应性 , 使 它 对 于 非 线 性 、非 平我们用该方法对 100 段采样数据 内含 183 ) 个梭形波进行梭形波识别. 以目测的结果为标稳数据同样具有很高的分解效率 . 得到的 IMF 表 1 100 段 EEG 数据( 内含 183 个梭形波) 在 Hilbert 谱上具有很好的时频局部性 . 试验结 的梭形波识别结果 果表明 ,它确实是在 EEG 中识别梭形波的理想 梭 形 检 出 漏检率误检率 正确率 工具. 漏检数 误检数 % % % 波 数 个 数 183 177 6 0 . 033 2 0 . 011 95 . 6 参 考 文 献 1 Leukel . Essential of Physiological Psychology. USA : The work Model : Application to Automatic Analysis of Human CV Company , 1978 Sleep . Computer and Biomedical Research , 1993 ,26 :157 王笑中 ,焦守恕 . 神经系统疾病症候学 . 北京 : 人民 ,1712 卫生出版社 J Doman , C Detka , T Hoffman , et al . Automating the ,1979 7 J Smith , I Karacan , M Yang. Automated Analysis of The Sleep Laboratory : Implementation and Validation of Digi2 3 () Human Sleep EEG. Waking and Sleep ,1978 ,2:75,83tal Recording and Analysis. International Journal of E Stanus , B Lacroix , M Kerkhofs , et al . Automated Sleep 4 Biomedical Computing , 1995 ,38 :277,290 Coring : A Comparative Reliability Study of Two Algo2 N Pradhan , P K Sadasivan. The Nature of Dominant Lya2 8 rithms. Electroencephalography and Clinical Neurophysiol2 punov Exponent And Attractor Dimension Curve of EEG in Sleep . Comput Biol Med , 1996 ,26 :419,428 ogy , 1987 ,66 : 448,454 E Huang , Z Shen , R L Steven , et al . The Empirical 5 9 J Pricipe , S K Gala , T G Chang. Sleep Staging Automa2 Mode Decomposition and the Hilbert Spectrum for Nonlin2 tion Base on The Theory of Evidence . IEEE Transaction on () ear and Non2stationary Time Series Analysis. Proc . R. Biomedical Engineering , 1987 ,365:503,509 N Schaltenbrand , R Lengelle , J P Macher. Neural Net2 6 Soc . Lond. A , 1998 ,903,995 A Novel Approach Ba sed on EMD f or Automatically Detecting Spindles from Sleep EEGs ) )11 ,2Yang ZhihuaQi Dongxu () 1School of Information Science and Tech. Sun Yat2sen University , 510275 , Guangzhou , China ) )2North China Univ. of Tech. , 100041 , Beijing , China ()Abstract A novel approach for automatically detecting spindles from sleep EEGs electroencephalo2graph ( ) is presented in this paper . Empirical mode decomposition EMDis employed to decompose a sleep EEG , which ( ) is usually typically nonlinear and non2stationary data , into a finite number of intrinsic mode functions IMF. ω) (Based on these IMFs , the Hilbert spectrum H , t of the EEG can be calculated easily , which provides a high resolution time2frequency presentation. Consequently , an algorithm to detect spindles from a sleep EEG accurate2 ly is developed. Experiments show encouraging detection results. ( ) ( ) Key Words empirical mode decomposition EMD; Hilbert2Huang transform HHT; automatic detection ; spindle