图像的盲分离

1、目录摘要1一盲分离原理31.1 盲信号处理的基本概念31.2盲信号处理的方法和分类31.3 盲源分离法41.3.1 盲源分离技术41.3.2 盲分离算法实现4二 ica基本原理5三 fastica算法原理及实现63.1 数据的预处理63.2 fastica算法原理73.3 fastica算法的基本步骤及流程图9四图像的盲分离仿真与分析10五总结14六 参考资料15附录1 matlab源程序16摘要简述了独立成分分析的基本原理以及利用fastica算法进行信号分离的理论依据，并通过matlab仿真实验实现了混合图像的盲分离，取得了较好的分离效果。结果表明该算法收敛速度快，有良好的分离效果，是一种

2、行之有效的混合图像盲分离方法。关键词：独立成分分析 fastica 盲分离abstractindependent component analysis is a new signal processing method whichdevelop rapidly during last few years. this paper introduces the basic model of the algorithm,analyzes the math principle of frequently-used rapid fixed-point algorithm for independent

3、component analysis,and applies the algorithm in blind separation of three images which are mixed randomly.the result shows that the algorithm is effective and reliable.key words:independent component analysis fastica blind signal separation一盲分离原理1.1 盲信号处理的基本概念盲信号处理是现代数学信号处理、算智能学近年来迅速发展的重要方向。电子信息、通信、

4、生物医学、图像增强、雷达、地球物理信号处理等众多领域有广泛的应用前景。盲信号处理就是利用系统(如无线信道、通信系统等)的输出观测数据,通过某种信号处理的手段, 获得我们感兴趣的有关信息(如原来独立发射的信号等)。盲信号的研究是当前学术界的一个研究热点，而盲信号分离则是盲信号研究中的一个重要的课题。bss是指从观测到的混合信号中分离出未知的源信号。盲信号中的“盲”意味着两个方面：第一，对源信号一无所知或只有少许的先验知识。第二，混合本身是未知的。这看似是一个不可能的任务，然而理论和实际都证实了只需要相当简单的假设，就可以得到该问题的解。这一特点使得bss成为一种功能相当强大的信息处理方法。混合矩

5、阵分离矩阵噪声向量s(t)x(t)y(t)图1.1 盲处理原理框图1.2盲信号处理的方法和分类在盲信号处理中，就源信号进过传输通道的混合方式而言，其处理方法可分为线性瞬时混合信号盲处理、线性卷积混合信号盲处理和非线性混合信号盲处理三类。根据通道传输特性中是否含有噪声、噪声特性（白噪声、有色噪声等）、噪声混合形式，可分为有噪声、无噪声盲处理，含加性噪声和乘性噪声混合信号盲处理等。按源信号和观测信号数目的不同可以将混合方式分为欠定 、适定和超定情况 ;按源信号特性的不同分为 : 平稳 、非平稳 、超高斯 、亚高斯 、超高斯和亚高斯混合分离等 。盲处理的目的可分为盲辨识和盲源分离两大类。盲辨识的目的

6、是求得传输通道混合矩阵（新型混合矩阵、卷积混合矩阵、非线性混合矩阵等）。盲源分离的目的是求得源信号的最佳估计。当盲源分离的各分量相互独立时，就称为独立分量分析，即独立分量分析是盲源分离的一种特殊情况。1.3 盲源分离法1.3.1 盲源分离技术盲源信号分离(blind source separation, bss) 是20世纪90年代迅速发展起来的一个研究领域 。它具有可靠的理论基础和许多方面的应用潜力。其在生物医学工程 、医学图像 、语音增强 、遥感 、通信系统 、地震探测等领域有着广泛而诱人的前景 ,盲源分离成为信号处理和神经网络领域的研究热点 。盲源分离是针对从检测的混合信号中估计或恢复源

7、信号的问题提出的，是指源信号、传输通道特性未知的情况下，仅由观测信号和源信号的一些先验知识（如概率密度）估计出源信号各个分量的过程。例如最著名的鸡尾酒会问题，仅根据多个麦克风检测信号分离或恢复出某种语音源信号。1.3.2 盲分离算法实现在实际的ica盲分离算法应用中,一般有时是必需的对观测数据做一些预处理技术，如用主成分分析 ( pca )降维和白化,用滤波器进行滤波降噪处理等。另外,由于恢复准则的局限以及先验知识的缺乏,盲信号分离方法只能得到源信号的波形,而无法确定信号的幅值以及信号之间的顺序。这两点都需要人为的制定规则(如规定信号的方差为1来确定幅值)来确定。2 ica基本原理独立分量分析

8、旨在对独立信源产生且经过未知混合的观测信号进行盲分离，从而重现原独立信源，其应用主要集中在盲源分离和特征提取两方面。ica问题可简单描述为：设有n个未知的源信号构成一个列向量,其中，t是离散时刻，取值为0,1,2，设a是一个维矩阵，一般称为混合矩阵（mixsing matrix）。设是由m个可观察信号构成的列向量，且满足下列公式1： 公式1bss的问题是，对任意t，根据已知的在a未知的条件下求未知的。这构成一个无噪声的盲分离问题。设是由m个白色、高斯、统计独立噪声信号构成的列向量，且满足下列公式2： 公式2则由已知的在a未知是求的问题是一个有噪声盲分离问题。ica的目的是对任何t，根据已知的在

9、a未知的情况下求未知的，ica的思路是设置一个维反混合阵，经过w变换后得到n维输出列向量，即有公式3： 公式3整个过程可以表示成如图2.1所示： 观测信号未知源信号 恢复信号 分离系统w 混合系统ay(t)x(t)s(t)图2.1 ica的线性模型如果通过学习得以实现，则，从而达到了源信号分离目标。3 fastica算法原理及实现ica算法的研究可分为基于信息论准则的迭代估计方法和基于统计学的代数方法两大类，从原理上来说，它们都是利用了源信号的独立性和非高斯性。基于信息论的方法研究中，各国学者从最大熵、最小互信息、最大似然和负熵最大化等角度提出了一系列估计算法。如fastica算法, info

10、max算法，最大似然估计算法等。基于统计学的方法主要有二阶累积量、四阶累积量等高阶累积量方法。本课程设计主要讨论fastica算法。3.1 数据的预处理 一般情况下，所获得的数据都具有相关性，所以通常都要求对数据进行初步的白化或球化处理，因为白化处理可去除各观测信号之间的相关性，从而简化了后续独立分量的提取过程，而且，通常情况下，数据进行白化处理与不对数据进行白化处理相比，算法的收敛性较好。若一零均值的随机向量满足，其中：为单位矩阵，我们称这个向量为白化向量。白化的本质在于去相关，这同主分量分析的目标是一样的。在ica中，对于为零均值的独立源信号，有：，且协方差矩阵是单位阵，因此，源信号是白色

11、的。对观测信号，我们应该寻找一个线性变换，使投影到新的子空间后变成白化向量，即： 公式4其中，为白化矩阵，为白化向量。利用主分量分析，我们通过计算样本向量得到一个变换 公式5其中和分别代表协方差矩阵的特征向量矩阵和特征值矩阵。可以证明，线性变换满足白化变换的要求。通过正交变换，可以保证。因此，协方差矩阵： 公式6 再将式代入，且令,有 公式7由于线性变换连接的是两个白色随机矢量和，可以得出一定是一个正交变换。如果把上式中的看作新的观测信号，那么可以说，白化使原来的混合矩阵简化成一个新的正交矩阵。证明也是简单的： 公式8其实正交变换相当于对多维矢量所在的坐标系进行一个旋转。在多维情况下，混合矩阵

12、是的，白化后新的混合矩阵由于是正交矩阵，其自由度降为，所以说白化使得ica问题的工作量几乎减少了一半。白化这种常规的方法作为ica的预处理可以有效地降低问题的复杂度，而且算法简单，用传统的pca就可完成。用pca对观测信号进行白化的预处理使得原来所求的解混合矩阵退化成一个正交阵，减少了ica的工作量。此外，pca本身具有降维功能，当观测信号的个数大于源信号个数时，经过白化可以自动将观测信号数目降到与源信号维数相同。3.2 fastica算法原理 fastica算法，又称固定点(fixed-point)算法，是由芬兰赫尔辛基大学hyvrinen等人提出来的。是一种快速寻优迭代算法，与普通的神经网

13、络算法不同的是这种算法采用了批处理的方式，即在每一步迭代中有大量的样本数据参与运算。但是从分布式并行处理的观点看该算法仍可称之为是一种神经网络算法。fastica算法有基于峭度、基于似然最大、基于负熵最大等形式，这里，我们介绍基于负熵最大的fastica算法。它以负熵最大作为一个搜寻方向，可以实现顺序地提取独立源，充分体现了投影追踪（projection pursuit）这种传统线性变换的思想。此外，该算法采用了定点迭代的优化算法，使得收敛更加快速、稳健。因为fastica算法以负熵最大作为一个搜寻方向，因此先讨论一下负熵判决准则。由信息论理论可知：在所有等方差的随机变量中，高斯变量的熵最大，

14、因而我们可以利用熵来度量非高斯性，常用熵的修正形式，即负熵。根据中心极限定理，若一随机变量由许多相互独立的随机变量之和组成，只要具有有限的均值和方差，则不论其为何种分布，随机变量较更接近高斯分布。换言之，较的非高斯性更强。因此，在分离过程中，可通过对分离结果的非高斯性度量来表示分离结果间的相互独立性，当非高斯性度量达到最大时，则表明已完成对各独立分量的分离。负熵的定义： 公式9 式中，是一与具有相同方差的高斯随机变量，为随机变量的微分熵 公式10 根据信息理论，在具有相同方差的随机变量中，高斯分布的随机变量具有最大的微分熵。当具有高斯分布时，；的非高斯性越强，其微分熵越小，值越大，所以可以作为

15、随机变量非高斯性的测度。由于根据公式10计算微分熵需要知道的概率密度分布函数，这显然不切实际，于是用如下近似公式： 公式11其中，为均值运算；为非线性函数，可取，或或等非线性函数，这里，通常我们取。快速ica学习规则是找一个方向以便具有最大的非高斯性。这里，非高斯性用式公式11给出的负熵的近似值来度量,的方差约束为1，对于白化数据而言，这等于约束的范数为1。fastica算法的推导如下。首先，的负熵的最大近似值能通过对进行优化来获得。根据kuhn-tucker条件，在的约束下，的最优值能在满足下式的点上获得。 公式12 这里，是一个恒定值， ,是优化后的值。下面我们利用牛顿迭代法解方程公式12

16、。用表示式公式12左边的函数，可得的雅可比矩阵如下： 公式13为了简化矩阵的求逆，可以近似为公式13的第一项。由于数据被球化，,所以，。因而雅可比矩阵变成了对角阵，并且能比较容易地求逆。因而可以得到下面的近似牛顿迭代公式： 公式14这里，是的新值，规格化能提高解的稳定性。简化后就可以得到fastica算法的迭代公式： 公式15 实践中，fastica算法中用的期望必须用它们的估计值代替。当然最好的估计是相应的样本平均。理想情况下，所有的有效数据都应该参与计算，但这会降低计算速度。所以通常用一部分样本的平均来估计，样本数目的多少对最后估计的精确度有很大影响。迭代中的样本点应该分别选取，假如收敛不

17、理想的话，可以增加样本的数量。3.3 fastica算法的基本步骤及流程图1.对观测数据进行中心化，使它的均值为0；2.对数据进行白化，。3.选择需要估计的分量的个数，设迭代次数4.选择一个初始权矢量（随机的）。5.令,非线性函数的选取见前文。6.。7.令。8.假如不收敛的话，返回第5步。9.令，如果，返回第4步。对观测信号去均值对去均值后的信号白化优化分离矩阵是否收敛否是求源信号图3.1 fastica算法流程图四图像的盲分离仿真与分析盲分离在图像处理领域得以应用，用matlab编程后仿真如下：此次仿真时取的图4.1为3张大小为256*256的原始图片及其直方图。图4.1原图片及其直方图 产

18、生随机混合矩阵，将原始图像混合后得到混合图像，如图4.2所示，可见原始图像已经看不出来。图4.2 混合后的3张图片 用fastica对上面的混合图像进行盲分离，即假定在未知源图像和混合矩阵下对混合图像进行分离，得到分离后的结果，如图4.3所示。图4.3 分离出的3张图片及其直方图 在将分离出的图片与原来的图片进行比较，可以发现分离出的图片基本上能够保留原来图像的信息。对原图和分离出的图片进行一个差值运算，并绘出差值图的直方图后可以清晰的看到。结果如图4.4所示。图4.4 原图和分离图的差值图及其对应的直方图由于ica问题本身具有一些不确定因素，包括：1)分离后结果的排序与源信号会不一致；2)分

19、离后的信号可能会与源信号相差一个负号。由于仿真实验在源图像未知情况下进行，因此我们在分离后的图形中会发现图像的排序发生变化，不过这些不会影响该算法对实际问题的处理。进行差值运算后，因为分离出的图片与原来的图片的灰度值很相近，相减后会是比较小的数，由于灰度值也表征亮度，所以比较小的数显示的图片比较黑。从直方图中也可以看出相减以后的值接近于0。由于ica问题本身的不确定性，我们不能完全恢复出原始的图片。五总结 在此次课设中，我运用了fastica 的方法, 此法不受源及噪声分布的影响, 也不受源个数的影响都明显改善了ica 的结果, 且具有较好的收敛性。通过对混合图像的实验仿真, 显示该方法较好的

20、分离出了原图像信号。在分析独立分量分析(ica)的基本模型和方法的基础上，详细地探讨了fastica算法，并通过仿真试验，成功地用该算法将3幅混合图像有效地分离出来，但是在自然界中需要处理许多被污染的图像，这些图像含有众多未知的噪声，我们如何利用压缩稀疏编码进行图像去噪，使得基本的ica模型可用，是我们需要进一步探讨的内容。此次课设，感触颇多。进一步巩固了对仿真软件matlab的使用，加深了自己编写程序的能力。最重要的是提高了自学的能力，对于一个自己比较陌生的课题，如何开展学习，如何查资料，如何一步一步学习新的知识。对整体的能力有一个较大的考验。最终成功的完成了任务，也深感自己研究的课题有很大

21、的意义，因此也为自己的努力感到些许欣慰！6 参考资料1 王家耀. 空间信息系统原理m. 北京: 科学出版社, 2001.2 杨行峻,郑君里.人工神经网络与盲信号处理m.清华大学出版社,2003:327-330.3 石庆研，黄建宇，吴仁彪.盲源分离及盲信号提取的研究进展j.中国民航大学学报，2007,25(03):1-7.4 赵艳.盲源分离与盲信号提取问题研究d. 西安:西安理工大学;2004.5 杨福生，洪波独立分量分析的原理与应用北京：清华大学出版社，20036贾银洁，许腾飞.基于ica的混合图像盲分离算法论文.山西电子技术.2008年第6期学报.7曹慧荣，张宝镭，马莉.基于快速ica算法的

22、混合图像分离论文.电脑学习.2004年2月.附录1 matlab源程序%读取和显示图像i1=imread(图片1.tif);i2=imread(图片2.tif);i3=imread(图片3.tif);figure(1)subplot(321),imshow(i1),title(原图片1),subplot(322),imhist(i1),title(原图片1的直方图);subplot(323),imshow(i2),title(原图片2),subplot(324),imhist(i2),title(原图片2的直方图);subplot(325),imshow(i3),title(原图片3);sub

23、plot(326),imhist(i3),title(原图片3的直方图);%对信号进行随机混合得到仿真观测信号s1=reshape(i1,1,256*256);s2=reshape(i2,1,256*256);s3=reshape(i3,1,256*256); %将s1、s2、s3变换成一个1*65536的行矩阵s=s1;s2;s3; %s为一个3*65536的矩阵sig=double(s); %sig为一个double型的矩阵%a= 0.950129285147175,0.4859824687093,0.456467665168341;% 0.231138513574288,0.891298

24、966148902,0.018503643248224;% 0.606842583541787,0.762096833027395,0.821407164295253;a=rand(size(sig,1); %生成一个大小为3*3，取值在0-1之间的随机矩阵mixedsig=a*sig; %混合信号为一个3*65535的矩阵ms1=reshape(mixedsig(1,:),256,256);ms2=reshape(mixedsig(2,:),256,256);ms3=reshape(mixedsig(3,:),256,256);%将ms1、ms2、ms3四舍五入转换成无符号整型数mi1=ui

25、nt8 (round(ms1); mi2=uint8 (round(ms2);mi3=uint8 (round(ms3);%显示混合图像figure(2)subplot(131),imshow(mi1),title(混合图片1);subplot(132),imshow(mi2),title(混合图片2);subplot(133),imshow(mi3),title(混合图片3);mixeds_bak=mixedsig;%将混合后的数据备份，以便在恢复时直接调用mixeds_mean=zeros(3,1);%标准化for i=1:3 mixeds_mean(i)=mean(mixedsig(i,

26、:);end %计算mixedsig的均值和方差for i=1:3 for j=1:size(mixedsig,2) mixedsig(i,j)=mixedsig(i,j)-mixeds_mean(i); endend %白化mixeds_cov=cov(mixedsig);%cov为求协方差的函数e,d=eig(mixeds_cov);%对图片矩阵的协方差函数进行特征值分解q=inv(sqrt(d)*(e);%q为白化矩阵mixeds_white=q*mixedsig;%mixeds_white为白化后的图片矩阵i=cov(mixeds_white);%i应为单位阵%fastica算法%x=

27、mixeds_white;%以下对x进行操作variablenum,samplenum=size(x);numofic=variablenum;%在此应用中，独立元个数等于变量个数b=zeros(numofic,variablenum);%初始化列向量b的寄存矩阵，b=b1,b2,bdfor r=1:numofic i=1; maxiterationsnum=150;%设置最大迭代次数 b=2*(rand(numofic,1)-.5);%随机设置b的初值 b=b/norm(b);%对b标准化 while i=maxiterationsnum+1 if i=maxiterationsnum %循环结束处理 fprintf(n 第%d分量在%d次迭代内并不收敛.,r,maxiterationsnum); break; end bold=b;%初始化前一步b的寄存器 u=1; t=x*b; g=t.3; dg=3*t.2; b=(1-u)*t*g*b+u*x*g)/samplenum-mean(dg)*b;%核心公式 b=b-b*b*b;%对b正交化 b=b/norm(b); i