EEGLAB runica.m

EEGLABの独立成分分析の御本尊はrunica.mのようです。 command lineから直接runicaを操作してみようと思います。
path(path,'.../eeglab4.512');
path(path,'.../eeglab4.512/functions');

help runica
とすると説明が表示されます。


・・・
引数・返り値が豊富に用意されているようです。
下記のdataで独立成分分析を行ってみます。

tic;[weights,sphere]=runica(data,'pca',8);toc

・・・

これも結構速そうです。 sphereはMatlabの球面座標を作成する関数ですが、 runicaの返り値の変数名がdefaultでsphereとなっていますので、そのまま使っています。
topographyと波形をみてみます。
load('.../vv_lout');
W=120;x=vv_mag(:,2);y=vv_mag(:,3);
xlin=linspace(min(x),max(x),W);
ylin=linspace(min(y),max(y),W);
[X,Y]=meshgrid(xlin,ylin);
B=weights*sphere;
IC=B*data;
AA=pinv(B);
A=sum(AA.^2,1);
P=diag(A)*IC;
V=AA./(ones(size(AA,1),1)*A);
figure;set(gcf,'color',[1,1,1]);
set(gcf,'renderer','zbuffer');
nn=8;
for n=1:nn;...
subplot('position',[0,1-n/nn,0.2,1/nn]);...
z=(V(1:2:204,n).^2+V(2:2:204,n).^2).^0.5;...
Z=griddata(x,y,z,X,Y,'linear');...
h=surf(X,Y,Z);axis off;set(h,'edgecolor','none');...
axis tight;daspect([1,1,0.01]);view([0,90]);...
subplot('position',[0.25,1-n/nn+0.02,0.7,1/nn-0.06]);...
plot(time,P(n,:));axis tight;grid on;end;

いろいろparameterをいじれそうですが、よくわかりません。

infomaxのウリは、チャンネル数が多い場合、jadeRより、速いことだそうです。
tic;[weights,sphere]=runica(data,'pca',20,'verbose','off');toc

B=weights*sphere;
IC=B*data;
AA=pinv(B);
A=sum(AA.^2,1);
P=diag(A)*IC;
V=AA./(ones(size(AA,1),1)*A);
figure;set(gcf,'color',[1,1,1]);
set(gcf,'renderer','zbuffer');
nn=20;
for n=1:nn;...
z=(V(1:2:204,n).^2+V(2:2:204,n).^2).^0.5;...
Z=griddata(x,y,z,X,Y,'linear');...
ha=subplot(nn,2,n*2-1);...
position=get(ha,'position');...
position=[0,position(2),1/nn,position(4)];...
set(ha,'position',position);...
hh=surf(X,Y,Z);axis off;axis tight;...
set(hh,'edgecolor','none');...
daspect([1,1,1]);view([0,90]);...
hb=subplot(nn,2,n*2);...
plot(time,P(n,:));axis tight;grid on;...
position=get(hb,'position');...
position=[1/nn+0.05,position(2),0.9-1/nn,position(4)];...
set(hb,'position',position);...
end;

各独立成分が生理学的に意味があるのかどうかよくわかりません。

Pentium4 3.2GHzでの計算時間(5回平均)
1) [W,S]=runica(data,'pca',n,'verbose','off');
2) B=jadeR(data,n);
3) [E,D]=pcamat(data,1,n);[nv,wm,dwm]=whitenv(data,E,D);[A,W]=fpica(nv,wm,dwm);
固有値数
n %
infomax
[sec]
jadeR
[sec]
FastICA
[sec]
10 (74.4)1.49360.36884.3250
20 (81.8)3.43125.215620.8968
30 (85.2)5.806250.887646.8129
40 (87.5)10.8000623.337670.2158
50 (89.2)13.7094
60 (90.7)18.5624
70 (92.0)18.6718
80 (93.1)23.8000
90 (94.1)28.8780
100 (95.0)36.6500
チャンネル数が多いときは、Infomaxがぶっちぎりで速そうです。尚、%は累積寄与率です。