转自:https://blog.csdn.net/lg1259156776/article/details/51132160

ICA算法的研究可分为基于信息论准则的迭代估计方法和基于统计学的代数方法两大类,从原理上来说,它们都是利用了源信号的独立性和非高斯性。基于信息论的方法研究中,各国学者从最大熵、最小互信息、最大似然和负熵最大化等角度提出了一系列估计算法。如FastICA算法, Infomax算法,最大似然估计算法等。基于统计学的方法主要有二阶累积量、四阶累积量等高阶累积量方法。

一:最大似然估计算法

1.1 目标函数部分

我们假设信号Si有概率密度函数Ps(t),由于我们假定信号源是相互独立的,其实经过白化处理后就变成独立的了;那么在给定时刻的联合分布函数为:

知道了信号源的联合分布Ps(t),再由分解矩阵S=WX,可以得出信号x的联合分布函数。

其中|W|为W的行列式。

由于没有先验知识,只知道原信号之间特征独立,且最多有一个是高斯分布,所以我们没有办法确定Ps(t)的分布,所以我们选取一个概率密度函数Ps^(t)来近似估计Ps(t)。

在概率论中我们知道概率密度函数由累积分布函数F(x)(cumulative distribution function,cdf)求导得到。

F(x)要满足两个性质:1,单调递增;2,值域在[0 1]范围

我们发现sigmoid函数的定义域是负无穷到正无穷,值域为0到1,缓慢递增的性质。基于sigmoid函数良好的性质,我们用sigmoid函数来近似估计F(x),通过求导得到Ps^(t)。

            

两函数图像如图:

怎么感觉pdf 函数g(s)这么像高斯分布函数呢??这个有点不理解(以后在研究吧)。

如果我们预先知道Ps(t)的分布函数,那就不用假设了;但是在缺失的情况下,sigmoid函数大多数情况下能够起到不错的效果。由于Ps(t)s是个对称函数,所以均值E[s]=0,那么E[x]=E[AS]=0,x的均值也是0。

知道了Ps(t),就剩下W了,在给定训练样本{Xi(Xi1,Xi2,........Xin),i=1,2....m个样本,样本的对数似然估计如下:

 (   W=H’)

T=m为独立同分布观测数据的样本数。最大化此似然函数就可获得关于参数W 的最佳估计。

1.2 优化部分(梯度下降算法)

接下来就是对W求导了,这里牵涉一个问题是对行列式|W|进行求导的方法,属于矩阵微积分

最终得到的求导后公式如下,logg'(s)的导数为1-2g(s)(可以自己验证):

当迭代求出W后,便可得到S=WX来还原出原始信号.

注意:我们计算最大似然估计时,假设了Xi与Xj之间是独立的,然而对于语音信号或者其他具有时间连续依赖特性(比如温度)上,这个假设不能成立。但是在数据足够多时,假设独立对效果影响不大,同时如果事先打乱样例,并运行随机梯度上升算法,那么能够加快收敛速度。

二:负熵最大的FastICA算法

2.1目标函数部分

2.1.1负熵判别准则

由极大熵原理可知,在方差相同的条件下,所有概率分布中,高斯分布的熵最大;因而我们可以利用熵来度量分布的非高斯性。

因此通过度量分离结果的非高斯性,作为分离结果独立性的度量;当非高斯性达到最大时,表明已完成对各个分量的分离。

因为FastICA算法以负熵最大作为一个搜寻方向,因此先讨论一下负熵判决准则。由信息论理论可知:在所有等方差的随机变量中,高斯变量的熵最大,因而我们可以利用熵来度量非高斯性,常用熵的修正形式,即负熵。

负熵的定义:   其中XG是和X具有相同协方差的随机变量,H()为变量的微分熵

微分熵定义:

联系极大熵原理,XG为高斯分布,所以J(X)>=0;当且仅当X本身也为高斯分布时=0;所以J(x)的值越大,证明X的非高斯性越强,

2.1.2负熵与独立性关系

假设n维随机变量X=[X1,X2……Xn],其互信息为I(X):

互信息即为:独立分布乘积分布与联合分布之间的负熵J(X),当Xi相互独立时,互信息为0;

由于计算J(X)需要联合分布函数和各个分量的分布函数,,这个显然不切实际;所以采用非线性变换g(x)后的均值期望来近似替代。

其中g(x)可以为:g(x)=tanh(a*x),g(x)=x^3;g(x)=x*exp(-x^2/2)  (为啥我也不清楚,都是对称的奇函数)

由于Xi即为观测数据X分离后的独立变量Si,再由中心极限定理可知,若随机变量X有许多相互独立的随机变量信号源Si相互组合而成,则不论Si为何种分布,观测变量数据X比Si具有更强的高斯性,换言之Xi的非高斯性更强。

所以,负熵J(X)的值越小,即此时的互信息I(X)越小,此时分离的变量Si独立性越好。

2.2 优化部分

快速ICA算法是找一个方向以便WX具有最大的非高斯性,也即最大的相互独立性;这里的独立性通过负熵来给出,通过均值近似估计来计算。这里通过白化处理,使W的范数为1,即使WX的方差估计为1;

优化过程推导比较复杂,公式太多就通过截图给出吧!

实践中,FastICA算法中用的期望必须用它们的估计值代替。当然最好的估计是相应的样本平均。理想情况下,所有的有效数据都应该参与计算,但这会降低计算速度。所以通常用一部分样本的平均来估计,样本数目的多少对最后估计的精确度有很大影响。迭代中的样本点应该分别选取,假如收敛不理想的话,可以增加样本的数量。

三:MATLAB源程序及说明

%下程序为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);    %计算协方差矩阵Cx

[eigvector,eigvalue]= eig(Cx); %计算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));

whileabs(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

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

return;

end

end

W(:,n)=WP;

end

Z=W'*Z;

%以下为主程序,主要为原始信号的产生,观察信号和解混信号的作图

clear all;clc;

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([0N -5,5]);title('源信号');

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

subplot(4,1,3);plot(s3);axis([0N -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');

参考文献:

1:http://wenku.baidu.com/view/941f6782e53a580216fcfe03.html

2:http://wenku.baidu.com/view/8fa7f71b6bd97f192279e9aa.html

3:http://wenku.baidu.com/view/d362d8330b4c2e3f57276387.html

其他的一些关于ICA的介绍:

https://blog.csdn.net/sinat_27763257/article/details/80568792

【机器学习】ICA算法简介相关推荐

  1. AI - 常见算法简介(Common Algorithms)

    机器学习常见算法简介 - 原文链接:http://usblogs.pwc.com/emerging-technology/machine-learning-methods-infographic/ 应 ...

  2. AI:人工智能领域算法思维导图集合之有监督学习/无监督学习/强化学习类型的具体算法简介(预测函数/优化目标/求解算法)、分类/回归/聚类/降维算法模型选择思路、11类机器学习算法详细分类之详细攻略

    AI:人工智能领域算法思维导图集合之有监督学习/无监督学习/强化学习类型的具体算法简介(预测函数/优化目标/求解算法).分类/回归/聚类/降维算法模型选择思路.11类机器学习算法详细分类(决策树/贝叶 ...

  3. ML与Optimality:最优化理论(GD随机梯度下降/QN拟牛顿法/CG共轭梯度法/L-BFGS/TR置信域/GA遗传算法/SA模拟退火算法)在机器学习中的简介、常用方法、案例应用之详细攻略

    ML与Optimality:最优化理论(GD随机梯度下降/QN拟牛顿法/CG共轭梯度法/L-BFGS/TR置信域/GA遗传算法/SA模拟退火算法)在机器学习中的简介.常用方法.案例应用之详细攻略 目录 ...

  4. 机器学习与数据挖掘,机器学习算法简介

    ​什么是数据挖掘 数据挖掘就是从大量的数据中去发现有用的信息,然后根据这些信息来辅助决策.听起来是不是跟传统的数据分析很像呢?实际上,数据挖掘就是智能化的数据分析,它们的目标都是一样的.但是,又有很大 ...

  5. 推荐系统算法_机器学习和推荐系统(二)推荐算法简介

    推荐算法简介 一. 基于人口统计学的推荐算法 二.基于内容的推荐算法 三. 基于协同过滤的推荐算法 协同过滤(Collaborative Filtering , CF) 基于近邻的系统过滤 基于用户( ...

  6. python线性回归算法简介_Python机器学习(二):线性回归算法

    机器学习研究的问题分为分类问题和回归问题.分类问题很好理解,而回归问题就是找到一条曲线,可以最大程度地拟合样本特征和样本输出标记之间的关系.当给算法一个输入时,这条曲线可以计算出相应可能的输出.回归算 ...

  7. 流数据机器学习平台SAMOA简介

    流数据机器学习平台SAMOA简介 SAMOA (Scalable Advanced Massive Online Analysis)是由雅虎开源的流数据在线机器学习的框架.和其它绝大多数大数据处理框架 ...

  8. 开根号的笔算算法图解_机器学习KNN算法之手写数字识别

    1.算法简介 手写数字识别是KNN算法一个特别经典的实例,其数据源获取方式有两种,一种是来自MNIST数据集,另一种是从UCI欧文大学机器学习存储库中下载,本文基于后者讲解该例. 基本思想就是利用KN ...

  9. ML之Clustering之K-means:K-means算法简介、应用、经典案例之详细攻略

    ML之Clustering之K-means:K-means算法简介.应用.经典案例之详细攻略 目录 K-means算法简介 1.K-means算法适用的数据类型​ 2.K-Means算法的全局最优解和 ...

最新文章

  1. 最新版动手学习深度学习和GAN电子书免费下载!
  2. mount /mnt/cdrom 提示错误信息“mount: can't find /mnt/cdrom in /etc/fstab or /etc/mtab”
  3. pcDuino上如何安装wordpress
  4. StringTokenizer类的使用
  5. 正面迎战英特尔,英伟达要和博世一起卖自动驾驶系统
  6. Serverless 实战——使用 Rendertron 搭建 Headless Chrome 渲染解决方案
  7. 五款常用邮件管理系统评测
  8. os引导程序boot 在根目录区寻找os加载程序文件loader 对应的根目录条目
  9. 【shell资源限制】RLIMIT_MEMLOCK too small
  10. java kotlin lateinit_15. Kotlin 究竟该不该用 lateinit?
  11. zlib-1.2.11手册
  12. 自然语言处理从零到入门 成分句法分析
  13. 蓝懿iOS 技术内容交流和学习心得 11.11
  14. 构建CA证书详解过程步骤
  15. BZOJ1778 [Usaco2010 Hol]Dotp 驱逐猪猡
  16. 自动化测试框架[各自动化测试框架大比拼]
  17. 服务器型号类型,cpu类型和型号服务器
  18. Jsp(二)——简单的Jsp程序
  19. iOS 客户端基于 WebP 图片格式的流量优化(下)
  20. 什么是大数据4v 指的是哪四个

热门文章

  1. 转战物联网·基础篇03-从JSON数据到短指令谈思维的转变
  2. 7621:硬币面值组合
  3. 欢迎大家关注我创办的北京圣思园科技有限公司IT高端就业培训
  4. padStart 和 padEnd的使用
  5. Shell脚本读取mysql结果集各数据项的值
  6. SolidWorks_画螺杆
  7. Native Bayes
  8. 四种常用聚类及代码(三):birch(一种层次聚类)
  9. 电脑开机提示 EFI PXE 0 for IPv6
  10. OpenCV中的Shi-Tomasi角点检测器