(完整word版)独立成分分析ICA

独立成分分析ICA

1.

PCA用于数据降维，而且只对高斯分布的数据有效。对于非高斯分布的数据，需要采用ICA进行BSS。

2.经典的鸡尾酒会问题：

假设在party中有n个人，他们可以同时说话，我们也在房间中一些角落里共放置n麦克风用来记录声音。宴会过后，我们从n麦克风中得到了一组数据

，i表示采样的时间顺序，也就是说共得到了m组

采样，每一组采样都是n维的。我们的目标是单单从这m组采样数据中分辨出每个人说话的信号。

也就是说：有n个信号源，，每一维都是一个人的声音信

号，每个人发出的声音信号独立。A是一个未知的混合矩阵（mixing matrix），用来组合叠加信号s，那么

这里的X是一个矩阵，其由采样数据构成。其中每个列向量是，

A和s都是未知的，x是已知的，我们要想办法根据x来推出s。这个过程也称作为盲信号分离。

令，那么

将W表示成

其中，其实就是将写成行向量形式。那么得到：

3.不确定性：

由于w和s都不确定，那么在没有先验知识的情况下，无法同时确定这两个相关参数。比如上面的公式s=wx。当w扩大两倍时，s只需要同时扩大两倍即可，等式仍然满足，因此无法得到唯一的s。同时如果将人的编号打乱，变成另外一个顺序，如上图的蓝色节点的编号变为3,2,1，那么只需要调换A的列向量顺序即可，因此也无法单独确定s。这两种情况称为原信号不确定。

还有一种ICA不适用的情况，那就是信号不能是高斯分布的，或者至多只能有一个信号服从高斯分布。

4.密度概率及线性变换

假设我们的随机变量s有概率密度函数（连续值是概率密度函数，离散值是概率）。为了简单，我们再假设s是实数，还有一个随机变量x=As，A和x都是实数。令是x的概率密度，那么怎么求？

公式如下：

推导过程如下：

5.数据预处理

一般情况下，所获得的数据都具有相关性，所以通常都要求对数据进行初步的白化或球化处理，因为白化处理可去除各观测信号之间的相关性，从而简化了后续独立分量的提取过程，而且，通常情况下，数据的白化处理能大大增强算法的收敛性。

6.FastICA算法

FastICA算法以负熵最大作为一个搜寻方向。由信息论理论可知，在所有等方差的随机变量中，高斯变量的熵最大，因而我们可以利用熵来度量非高斯性，常用熵的修正形式，即负熵。根据中心极限定理，若一随机变量X由许多独立的随机变量s之和组成，只要si具有有限的均值和方差，则不论其服从何种分布，随机变量X较s更接近高斯分布。换言之，si较X的非高斯性更强。在分离过程中，可通过对分离结果的非高斯性度量来表示分离结果间的相互独立性，当非高斯性度量达到最大时，则表明已完成对各独立分量的分离。

7.FastICA算法的基本步骤

ICA函数：

function z = ICA(X)

%去均值

[M,T] = size(X);

average = mean(X')';

for i=1:M

X(i,:)=X(i,:)-average(i)*ones(1,T);

end

%白化

Cx=cov(X',1); %计算协方差矩阵

[eigvector, eigvalue] = eig(Cx); %计算特征值和特征向量

W = eigvalue^(-1/2)*eigvector'; %白化矩阵

z = W*X; %正交矩阵

%迭代

Maxcount = 10000; %最大迭代次数

Critical = 0.00001; %判断是否收敛

m = M; %需要估计的分量的个数W = rand(m);

for n = 1:m

WP = W(:,n); %初始权向量（任意）

% Y = WP'*z;

% G = Y.^3; %G为非线性函数，可取y^3等

% GG = 3*Y.^2; %G的导数

count = 0;

LastWP = zeros(m, 1);

W(:,n) = W(:,n) / norm(W(:,n));

while abs(WP - LastWP)&abs(WP+LastWP)>Critical

count = count + 1; %迭代次数

LastWP = WP; %上次迭代的值

% WP = 1/T * z * ((LastWP'*z).^3)'-3*LastWP;

for i = 1:m

WP(i)=mean(z(i,:).*(tanh((LastWP)'*z)))-(mean(1-(tanh((LastWP))'*z).^2)).*LastWP(i);

end

WPP = zeros(m, 1);

for j = 1:n-1

WPP = WPP+(WP'*W(:,j))*W(:,j);

end

WP = WP - WPP;

WP = WP / (norm(WP));

if count == Maxcount

fprint('未找到相应的信号');

return;

end

end

W(:,n) = WP;

end

z = W'*z;

主程序—信号生成及分离：

N = 200; n = 1:N; %N为采样点数

s1 = 2 * sin(0.02 * pi * n);

t = 1 : N; s2 = 2 * square(100 * t, 50); %方波信号

a = linspace(1, -1, 25); s3 = 2 * [a, a, a, a, a, a, a, a]; %锯齿信号

s4 = rand(1, N); %随机噪声

S = [s1; s2; s3; s4]; %信号组成4 * N

A = rand(4, 4);

X = A * S; %观察信号

%源信号波形图

figure(1); subplot(4, 1, 1); plot(s1); axis([0 N -5, 5]); title('源信号');

subplot(4, 1, 2); plot(s2); axis([0 N -5, 5]);

subplot(4, 1, 3); plot(s3); axis([0 N -5, 5]);

subplot(4, 1, 4); plot(s4); xlabel('Time/ms');

%混合信号波形图

figure(2); subplot(4, 1, 1); plot(X(1, :));title('混合信号');

subplot(4, 1, 2); plot(X(2, :));

subplot(4, 1, 3); plot(X(3, :)); subplot(4, 1, 4); plot(X(4, :));

z = ICA(X);

figure(3); subplot(4, 1, 1); plot(z(1, :)); title('解混后的信号');

subplot(4, 1, 2); plot(z(2, :));

subplot(4, 1, 3); plot(z(3, :));

subplot(4, 1, 4); plot(z(4, :)); xlabel('Time/ms');

运行结果：

源信号

Time/ms

混合信号

020406080100120140160180200

020406080100120140160180200

020406080100120140160180200

解混后的信号

020406080100120140160180200

020406080100120140160180200

020406080100120140160180200

Time/ms