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

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

独立成分分析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

相关文档
最新文档