独立分量生物医学运用探究

1、独立分量生物医学运用探究1独立分量分析；独立分量分析（13）是统计信号处理近年 来的一项发展。顾名思义，这是一种分解技术，其特点是把信 号分解成若干相互独立的成分。主分量分析（ps和奇异值 分解（svd）是人们较熟悉的分解信号的线性代数方法,ica与 它们的主要不同之处表现在:；（1）后者只要求分解出来的各 分量互相正交（不相关），但并不要求它们互相独立。用统计 信号处理的语言来表达，即:后者只考虑二阶统计特性，而前 者则要更全面考虑其概率密度函数的统计独立性。；（2）后 者按能量大小排序来考虑被分解分量的重要性。这样的分解 虽然在数据压缩和去除弱噪声方面有其优点,但分解结果往 往缺乏明确的生

2、理意义。前者虽然分解出的分量其能量大小 存在不确定性，但当测量值确实是由若干独立信源混合而成 时，分解结果往往具有更好的生理解释。由于测得的生理信 号往往是若干独立成分的加权迭加（例如，诱发脑电总是被 自发脑电所淹没,而且常伴随有心电、眼动、头皮肌电等干 扰），此13是一项值得注意的分解方法。；此外，神经生理 研究认为，人类对认知、感知信息的前期处理有"去冗余” 的特点。13在这方面也表现出类似特性，因为互相独立的分 量之间互信息是最少的。13是伴随着盲信号处理，特别是盲 信源分离发展起来。其研究热潮方兴未艾，也正在引起生物 医学工程界的注意,ieeetransbs正在组织出版以它为

3、重点 的专辑。就国际范围看，以下几个研究单位目前工作比较领 先:(1)美国加州大学生物系计算神经生物学实验室,(2)日 本ri ken脑科学研究所脑信息研究室,(3)芬兰赫尔辛基工业大学计算机及信息科学实验室，目前发表有关文献较多的刊 物有 ieeetrans 的 sp 和 nsi 以及 isfeural co叩ut at i on 等。本 文目的是对13的原理、算法及应用作一简述，以引起国内 同行对它的关注。将侧重于概念说明，而不追求数学上的严 谨性。；2原理；2.1问题的提法,s-(n)是一组互相独立的 信源,a是混合矩阵,x- (n)是观察记录，即x- (n) =as-(n)。问 题的任

4、务是:在a阵未知且对s- (n)除独立性外无其它先验知 识的情况下，求解混矩阵b,使得处理结果y- (n) =bx- (n)中各 分量尽可能互相独立，且逼近s(n)o容易理解,解答不是唯一 的，它至少受以下条件的限制:(1)比例不定性:s-(n)中某一 分量大k倍时,只要使相应的a阵系数减小k倍,x- (n)便保 持不变。;因此,求解时往往把s-(n)假设成具有单位协方差 阵，即s-中各分量均值为零,方差为1,且互相独立。(2)排序 不定性:y-与s-中各分量排序可以不同。因为只要对调b阵 中任意两行,y-中相应元素的位置也便对调。(3)s-(n)中至 多只能有一个高斯型信源:这是因为高斯信源

5、的线性组合仍 是高斯型的，因此混合后便无法再区别。(4)信源数目n只能 小于或等于观测通道数m nsgt;mt青况目前尚未解决。以下讨论设因此,y- (n)只是在上述条件下对s- (n)的逼近。换名话说，任务的实质是优化问题，它包括两个主要方面:优 化判据(目标函数)和寻优算法。；2.2目标函数；这一领域 的研究者已经从不同角度提出了多种判据。其中以互信息极 小判据(m ni mzati onof l4itual infornati onz 简记 w)和信 息或%商极大判据(inforrrax或晞xi niizati onof entropy,简记 m)应用最广。由于最基本的独立性判据应由概率

6、密度函数 (probabi i i tydensi tyf uncti on,简记 pdf)引出，而工作时 pdf 般是未知的，估计它又比较困难，因此通常采用一些途 径绕过这一困难。；常用的方法有两类:把pdf作级数展 开,从而把对pdf的估计转化为对高阶统计量的估计;在图 1的输出端引入非线性环节来建立优化判据。后一作法实际上隐含地引入了高阶统计量。(1)互信息极小判据:统计独立 性的最基本判据如下：令p(y)是y的联合概率密度函 数,pi(yi)是y-中各分量的边际概率密度函数。当且仅当y- 中各分量独立时有：p(y-) = -|n=1pi(yi)因此用p(y-)与口 i =1 pi (y

7、i)间的kul i back- lei bl er散度作为独立程度的定量 度量:i(y-)=klp(y-), n n =1pi (yi) = j p(y-)i ogp(y-) r|n=1pi(yi)dy-(1)显然fi(y-)of当且仅当各分量独立时i(y-)=0o因此,互信息极小判据的直接形式是:在y二bx条件下寻找b使(1)式的i (y)极小为了使判据实际可用，需要 把1( y-)中有关的pdf展成级数。；由于在协方差相等的概 率分布中高斯分布的燔值最大，因此展开时常用同协方差的 高斯分布作为参考标准。例如，采用gramcharlier展开时w : p(yi)pqyi) =1+13! k2

8、yi h3( y-i)+u! k4yi h4( yi) + 式中pqyi)是与p(yi)具有同样方(o 2=1)和均值(p=0)的高斯分布。k3yi、k4yi是yi的三、四阶累计量(cunwl ant), hn(yi)是n阶herrnt多项式。此外还有许多其他展开办法，如edgewrt h展开，利用负熾hbgentropy)等。不论采用何种展开方式，经推导后总可把式(1)近似改成k3sk4的函:i(y)=f(k3y-f k4y-f b)(1)f f()的具体形式多种 多样，视推导时的假设而异。；这样就得到互信息判据的实 用近似形式:在y=bx条件下寻找b,使式(1)的i(y)极小 (2)inf

9、onax判据:这一判据的特点是在输出端逐分量地引入 个合适的非线性环节把yi转成ri (如图2)。可以证明，如 果gi(-)取为对应信源的累积分布函数cdf(它也就是概率 密度函数的积分)，则使r-=(r1.rmt的燔极大等效于使 i(y-)极小,因此也可达使y-中各分量独立的要求。从而得到 inforax判据:在选定适当gi( )后，寻找b使® hr-)极大 需要指出的是，虽然理论上gi( )应取为各信源的cdf,但实践证明此要求并不很严格，有些取值在01之间的单调升函 数也可以被采用，如si gm)i d函数、tanh()等。估计hr-)固然也涉及pdf,但由于其作用已通过刃()

10、引入，所以可以 不必再作级数展开而直接用自适应选代寻优步骤求解。文献 中还提出了一些其他判据,如极大似然、非线性p3等，但它 们本质上都可统一在信息论的框架下，所以不再一一列举1o ； 3处理算法优化算法；可大致分为两类，即批处理与 自适应处理。；3.1批处理批处理比较成熟的方法有两类。较早提出的是成对旋转法2,其特点是把优化过程分解成 两步。先把x-(n)经wp车加以 球化”得z-(n),使z-(n)t=in 即:各分量不相关且方差为%然后再寻找合适的正交归一阵u达到使y各分量独立的目的。前一步类似于pca后一步则 可利用g vens旋转,根据目标函数，将z-中各分量两两成对反复旋转直到收敛

11、。这种方法计算量较大。1999年,gadoso 提出几种方法对它作了进一步改进3,其中包括：sxkurt 法、j ade法、sh bbs法等，限于篇幅，本文不再叙述。近年 来，提出的另一类方法是所谓固定点"法(fi xedpoi nt陇t hod) 4f 5,其思路虽来源于自适应处理，但 最终算法属于批处理。；简单地说,通过随机梯度法调节b 阵来达到优化目标时，有:b(k+1)=b(k) +a b(k)a b(k)=-p £ kb(k)式中k是选代序号疋k是瞬时目标函数。当到达稳 态时必有e是总集均值算子：ea b( k) =0( 2)如果a b(k) 与b( k)有关，就

12、可由(2)式解出b的稳态值。不过由于(2)式 总是非线性方程，因此求解时仍需要采用数值方法(如牛顿 法、共辄梯度法等)迭代求解。实践证明，不论是收敛速度还 是计算量，此法均优于前一种方法，而且它还可以根据需要 逐次提取最关心的yi,因此是一类值得注意的方法。；3. 2 结合神经网络的自适应处理结合神经网络的自适应处理算 法的框图。1994年ci chocki提出的调节算法 是：b(k+1)=b( k)+ a b(k) a b(k)= p ki-屮(y-k) 0) t(y-k)b(k)式中屮、0都是n维矢量，其各元素都是单调升 的非线性函数：屮(yk) =sgnyk y2k, e ty- k=3

13、tanh( 10yk)所 得结果虽令人鼓舞，但是方法是经验性的。其后学者们从理 论上沿着这一方向作了更深入的讨论,并发展出多种算法。概括地说，主要发展有以下几点:；(1)引入自然梯度(或相对 梯度)。按照最陡下降的随机梯度法推导出的系数调节公式 往往具有如下一般形式:a b(k)=p kb-t(k)-屮(y-k)x-tk 式中的屮(y-k)视具体算法而异。infonax法中屮()由所选 用的g( )决定;w法中则与yk的三、四阶矩有关。b-t(k) 是矩阵求逆再转置，它的计算量很大。arar i 7在1998年提 出将最陡下降梯度改为"自然梯度"，两者间关系是:自然 梯度珂

14、最陡下降梯度 bt( k) b( k)于是有：a b(k) = p k&t(k) 屮 (yk)xtkbt(k)b(k)二 p ki- 屮 (y- k)y-tk b( k)由于此式避免了矩阵求逆，因此计算量明显 降低且收敛加快。目前,这一作法已被普遍接受。；(2)引入 自然梯度后，采用不同的优化判据得出的调节公式虽各有千 秋，但大致都可表示为如下的"串行更新形 式:b(k+1)=b(k)+a b(k)=i+hy-k)b(k)只是 ky-k)的具体形式各不相同。串行矩阵更新的算法还具有一些理论上值 得注意的性质，如均匀特性(uni f or npr operty)和等变性 (eq

15、ui vari ant)等8, 9。； (3)四阶累计量 k48gt; 0 的超高 斯信号和k4& t; 0的欠高斯信号，其处理过程应当予以区别。 采用同一算法效果往往不好。目前的办法多是在调节公式中 引入一个开关。根据估计得k4的符号来切换不同算法，如扩 展的infonax法就是一例10o此法的系数调节公式是:b(k) =p ki-ktanh( y- k) - y-tk-y- ky-tk b( k)其中 k 是对角 阵，其对角元素之值为+1或视该信号分量k48gt; 0或 &t;0而定。为了实时应用，估计k4也可采用递归算法。总之，自适应算法是目前采用较广的方法。；4应用举例

16、；4.1 仿真计算为检验经ica算法分解信源的能力，左图是一组源 信号，它们对系统来说是未知的。这一组信号经混合后的观 察信号作为(中图所示)13算法的输入，分解后的结果如右 图所示。可以看到，除了波形的次序、极性和波幅发生变化 之外,源信号的波形被很好地分解出来。一般情况下,临床脑 电信号中既有超高斯成分（如诱发电位），也有亚高斯成分 （如肌电和工频干扰）o为了检验扩展i nforrax算法处理这类 情况的能力，我们又用此法进行了如图6所示仿真实验。左 图第一行是一段自发脑电信号,第二行是仿真的视觉诱发电 位,第三行是肌电干扰。混合后的信号（图中第二列所示）经 ica分解得到如右图所示的结果

17、。这一结果表明扩展13算法在同时存在超高斯和亚高斯信号的情况下,仍然能够很好 地实现盲分解。但应指出：这一仿真结果并不说明通过ica 分解就能直接得到视觉诱发电位，因为还没有涉及头皮上的 多导数据。；4.2实验已卩分析（1）多导脑电观察中£卩的增 强:需要强调，把多导脑电作13分解后直接取出其中与vep 有关的成分,得到的并不是头皮电极处的vep分量，因为它们 只是分解出来的信源，而这些信源的位置并不在头皮上，为 了得到电极处测量值中的vep成分，需按下述步骤处理:用训 练得的曲车直接对头皮上取得的多导脑电数据进行ica分解, 得到各独立分量组成的矩耻y=bx（见图7a）;再根据各分

18、量的 波形特征及产生时段，选择与vep有关的一部分分量（例如在 前300ns中具有较大幅度的分量），并将其余分量置0,得到新 的独立分量矩阵y'再反变换回头皮各电极处得/ =1-y 这样才能得到去除噪声和干扰后各电极处的vepo ；采用这 样的方法可显着地减少提取vep所需要的累加次数。左图是 经3次累加所得vep,中图是经50次累加所得结果，右图则是 用左图经图7中13处理后提取的vep。比较中、右两图， 两者波形趋势基本相同，但后者比前者其主要峰、谷显然更 清楚,而累加次数由50减到3o (2)13分量的空间模式:把 某一个13分量的瞬时值经& 1逆推回头皮各电极处得x-' 后,就可以按断层图的插补方法得到该时该分量在头皮上的 空间分布模式。这个空间分布模式也可以用更简单办法得到: 只要把逆矩阵中相应于某ica分量的列中各元素的值赋 与头皮各电极处，再作断层图插值,就可以表现该13分量在 任意时刻的空间分布模式。也就是:xf i(t)=bf ijy1 j(t),i=1n式中h ij是81的第i行第j列元素。；可见 13分量/ j(t)在头皮各电极处的对应值等于用逆阵b-1第j列各元素来对/ j(t)加权。因