基于NSCT的独立分量分析方法

张军龙，刘立程，俞坤游

（广东工业大学 信息工程学院，广东 广州 510006）

基于NSCT的独立分量分析方法

张军龙，刘立程，俞坤游

（广东工业大学 信息工程学院，广东 广州 510006）

为了提高独立分量分析（ICA）算法的分离性能，提出了基于非下采样轮廓波变换（NSCT）的ICA方法。首先对混合图像进行NSCT，使混合图像的非高斯性在NSCT域得以增强，然后在变换域求解分离矩阵。理论分析和仿真实验均表明，该方法可以获取更好的分离精度和更快的收敛速度。

独立分量分析；非下采样轮廓波变换；固定点算法；峭度；性能指数

独立分量分析ICA（Independent Component Analysis）是20世纪90年代发展起来的一种新的信号处理技术，它是从多维统计数据中找出隐含因子或分量的方法。目前ICA算法的研究可分为基于信息论准则的迭代估计方法和基于统计学的代数方法两大类[1-2]，从原理上来说，它们都是利用了源信号的独立性和非高斯性。迄今为止，各国专家学者提出了一系列估计算法，如FastICA算法、Infomax算法、最大似然估计算法、基于二阶累积量或四阶累积量等高阶累积量算法[1-2]。

ICA问题实际上包括了两个部分：一是采用什么样的判据（目标函数或者代价函数）来判定一组信号是否相互独立或者接近相互独立；二是用什么样的算法来实现这个目标。实际工作中常用的目标函数有互信息最小化目标函数、传输信息最大化或者负熵最大化目标函数和ICA最大似然目标函数。可以证明[1]，在信息论框架下各种目标函数是统一的。因此，必须通过概率密度函数（PDF）的级数展开或者采用适当的非线性函数来简化和逼近。比如信息最大化目标函数（Infomax）通过引入非线性函数来代替高阶统计量的估计[2]。

理论分析表明[1-3]，ICA算法的性能与非线性函数的选择和源信号的概率分布有关。为了改善ICA算法的性能，通常的做法是根据源的概率密度估计来选择合适的非线性函数，使之尽量符合源的概率分布。本文换一种思路，从改变源的概率分布出发来改善ICA算法的性能。考虑到把图像变换到NSCT（Nonsubsampled Contourlet Transform）域，其高频方向子带具有更大的峭度，然后在变换域直接求解分离矩阵，就可以获得更高的分离精度和更快的收敛速度。实验结果证明了此方法的有效性。

1 非下采样Contourlet变换

1.1 NSCT变换

轮廓波（Contourlet）[4-5]作为一种多尺度和多方向的图像表示方法，已经在自然图像去噪[6]、图像融合[7]等领域[8]获得了广泛的研究，展示了这一多尺度几何分析方法的强大潜力。Contourlet最初的提出是为了寻找图像中分段光滑的轮廓信号的稀疏表示，DO M N和VETTERLI M认为Contourlet是一种 “真正”的二维图像表示方法。然而由于Contourlet在实现过程中对滤波后的图像进行了隔行隔列的重采样，使低频子带和高频子带产生频谱混叠。频谱混叠造成同一方向的信号会在几个不同的方向子带同时出现，从而在一定程度上削弱了其方向选择性[4]。

为了消除Contourlet变换的频谱混叠现象，增强其方向选择性和平移不变性。基于Contourlet变换和非下采样思想，CUNHA A L、ZHOU J P和 DO M N等于 2006年利用非下采样塔式分解和非下采样滤波器组构造出了非下采样轮廓波变换 （NSCT）[9]。由于没有下采样操作，NSCT具有平移不变性。

NSCT变换包括两个步骤：第一步采用拉普拉斯金字塔LP（Laplacian Pyramid）来实现多尺度分解，第二步通过方向滤波器组DFB（Directional Filter Bank）来实现方向分解。

1.2 非下采样LP分解

在Contourlet变换的LP分解中,首先对上一尺度低频图像用低通滤波器进行低通滤波,然后进行下采样得到低频图像；再对该低频图像进行上采样，然后用高通滤波器对上采样后的图像进行高通滤波，并将高通滤波后的图像与上一尺度的低频图像进行差分，得到塔式分解后的高频部分。

非采样的LP分解不同于Contourlet变换中的LP分解，在NSCT中，采用Atrous算法的思想，对低通滤波器和高通滤波器分别进行上采样，然后对上一尺度低频图像采用上采样后的低通滤波器进行低通滤波，得到低频图像；对上一尺度低频图像用上采样后的高通滤波器进行高通滤波，得到LP分解后的高频图像。图像经N级非采样塔式滤波后,可得到N+1个与源图像具有相同尺寸大小的子带图像。

1.3 非下采样DFB分解

2 基于NSCT的ICA算法

2.1 NSCT域高频子带概率分布特性

理论分析和实验[4,10]都表明，近高斯分布的自然图像（大多为亚高斯分布）其Contourlet域的高频子带为超高斯分布（接近 Laplace分布），而作为 Contourlet算法改进的NSCT，其变换域的高频方向子带也呈现非高斯分布。图2给出了图像pepper作一次NSCT分解产生的低频子带和两个高频子带的统计直方图和归一化峭度值。其中，归一化峭度定义为[3]：

图2 NSCT变换前后统计直方图与归一化峭度

由图2可以看出，经过NSCT变换后，高频子图像的峭度约为原图像峭度的15倍（图像信号一般为亚高斯信号，峭度为负值）。实验表明，对于其他的图像也有类似的结果。

2.2 算法的可行性分析

ICA最基本的模型为

其中,x=（x1,x2,x3,…,xm）T是m维的零均值随机观测信号，它是由n个未知的零均值独立源信号s=（s1,s2,…,sn）T线性混合而成的，A=[a1,a2,…,an]是 m×n阶（通常 m和 n相等）满秩的混合矩阵；aj为混合矩阵的m维列向量。考虑到NSCT是一种线性变换，对式（2）和式（3）两边进行NSCT（忽略时刻 t），有

其中,T表示NSCT变换。所以时域的混合矩阵A及分离矩阵W与NSCT域的A和W是一致的，在NSCT域得到的分离矩阵W可直接用来分离原混合信号，从而得到源信号,没有必要对NSCT变换域的独立分量求取反变换。

2.3 算法的收敛性分析

收敛性分析如下[2,12]。 估计信号y（忽略时间参数t）的任何一个分量为:

其中，giT为全局传输矩阵的第i行。由（6）式知，估计信号y的任何一个分量可以表示为相互独立的源信号的线性组合。由高阶累积量的性质[1]：kurt（s1+s2）=kurt（s1）+kurt（s2）和 kurt（βs1）=β4kurt（s1），得到以峭度作为非高斯性准则的代价函数表达式[2]:

其中，α为步长参数。由式（5）知，以峭度作为非高斯性度量的算法，其收敛速度与源信号的峭度成正比。可以推出，基于ICA各种代价函数的等价性，在使用其他目标函数的情况下也能得出相同的结论。

2.4 算法的具体实现步骤

假设混合的N幅图像已经按行堆叠，每一行代表一幅混合图像，共计N行。具体实现步骤如下：

（1）把混合信号每一行还原成图像信号，则可以获得N幅图像；

（2）对还原的图像作NSCT变换，分别获取相应的低频子图和高频子图（共有N幅低频子图，其他的全是高频子图）；

（3）选择某个分解尺度、某个方向的高频子图，将其按行重新排列成矩阵；

（4）在新矩阵上运行ICA算法（如FastICA算法）,得到分离矩阵W和相互独立的NSCT域分量；

（5）分离矩阵W直接乘以原时域混合图像信号，可以得到时域分离的独立分量（矩阵），把每一行还原成二维的图像形式，就可以得到N幅源图像。

3 实验仿真与分析

实验参数如下：NSCT分解层数为两层，Pyramidal滤波器为 maxflat，Directional滤波器为 dmaxflat7，这样就可以获取一个低频子带和两个高频方向子带（实验取第2个方向子带）。采用的ICA算法为FastICA，分离性能评价标准采用式（1）的PI性能指数，实验图像大小取为256×256的标准灰度图像，实验时混合矩阵A由系统随机生成。

图3给出了混合前的图像，图4给出了混合后的图像。

图3 混合前的灰度图像

图4 混合后的灰度图像

图5为直接运行FastICA算法得到的图像，图6为在NSCT域运行FastICA算法得到的图像。PI性能指数分别为 0.223和0.007 3。

图5 FastICA直接解混后的灰度图像

图6 变换域FastICA算法解混后的灰度图像

为进一步验证算法的性能，反复实验100次,结果如图7所示。

图7 迭代次数与PI性能指数

从图7可以看出，在NSCT域ICA算法的性能非常好，PI指数基本稳定在0.008 7左右，迭代次数基本稳定在25次左右，视觉效果也非常好。而在时域直接运行FastICA算法，视觉效果并不理想，PI性能指数基本上迭代次数超过150仍未收敛（概率约为0.08，绘制图7时已排除此情况）。

实验证明，在小波域、Contourlet域以及NSCT域，FastICA算法的收敛精度及收敛速度均有较大的提高。总之，从改变源信号的概率分布的角度出发，通过线性变换把混合信号变换到变换域,在变换域进行相关处理。理论和实验均证明，该方法行之有效，同时也是一个重要的值得研究的方向。

[1]杨福生,洪波.独立分量分析的原理与应用[M].北京:清华大学出版社,2006.

[2]HYVARINEN A,KARHUNEN J,OJA E.Independent component analysis[M].New York:Wiley,2001.

[3]阴法明,梁瑞宇.改进的FastICA算法在图像盲分离中的应用[J].电子技术应用,2013,39（1）:116-118.

[4]DO M N,VETTERLI M.Contourlet:a directional multiresolution image representation[C].Proceedings of IEEE International Conference on Image Processing,Rochester,NY,2002:357-360.

[5]DO M N,VETTERLI M.The contourlet transform:an efficient directional multiresolution image representation[J].IEEE Transactions on Image Processing,2005,14（12）:2091-2106.

[6]陈增境,汪西原,马润渊.基于ICA的有噪图像特征提取研究[J].微型机与应用，2010,29（8）:40-42.

[7]刘胜鹏,方勇.基于Contourlet变换和PCNN的融合算法及其在可见光与红外线图像融合中的应用[J].红外与毫米波学报,2007,26（3）:364-369.

[8]焦李成,侯彪,王爽,等.图像多尺度几何分析理论与应用-后小波分析理论与应用[M].西安:西安电子科技大学出版社,2008.

[9]CUNHA A L,Zhou Jianping,DO N M.The nonsubsampled contourlet transform:theory,design,and applications[J].IEEE Transactions on Image Processing,2006,15（10）:3089-3101.

[10]PO D D Y,DO M N.Directional multiscale modelings of image using the contourlet transform[J].IEEE Transactions on Image Processing,2006,15（6）:1610-1620.

[11]CROUSE M S,NOWAK R D,BARANIUK R G.Waveletbased statistical signal processing using hidden markov models[J].IEEE Transactions on Signal Processing,1998,46（4）:886-902.

[12]舒朗,舒勤,苏静.新息模型的独立分量分析方法[J].计算机应用,2011,31（2）:556-558.

Independent component analysis based on NSCT

Zhang Junlong,Liu Licheng,Yu Kunyou

（Faculty of Information Engineering,Guangdong University of Technology,Guangzhou 510006,China）

In order to improve the separation performance of the independent component analysis（ICA）algorithm,a new method based on nonsubsampled contourlet transform（NSCT）is proposed in this paper.Firstly,the mixed images are decomposed using NSCT,which increases the non-Gaussianity in the NSCT domain,then the separated matrix can be obtained in the transform domain directly.Theory analysis and simulation experiments all show that the new method can achieve superior performance in both convergence rate and separation accuracy.

ICA;NSCT;FastICA;kurtosis;performance index

TP391

A

1674-7720（2013）20-0044-04

2013-08-16）

张军龙，男，1988年生，硕士研究生，主要研究方向：图像处理，信号处理。

刘立程，男，1972年生，博士，副教授，主要研究方向：通信信号处理。

俞坤游，男，1988年生，硕士研究生，主要研究方向：仪器仪表工程。