function B = bell(X)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%% Bell & Sejnowsky のアルゴリズムで混合データ x を分離する行列 %%%%%
%%%%%% w を求めるプログラム %%%%%
%%%%%% 前処理の部分はCardosoのプログラムから一部拝借 %%%%%
%%%%%% 樺島祥介 %%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% データ行列のサイズを調べる
[N,P]=size(X);
% 定行列の設定
Id=eye(N);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%% 前処理 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%% 各行の平均値を引く %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MX=mean(X');
X=X-mean(X')'*ones(1,P);
%%%%%%%%%%%%%%%%%% ホワイトニング %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[U,D] = eig((X*X')/P) ; % 共分散行列の固有値と固有ベクトル
% U: 固有ベクトル D: 固有値
[puiss,k] = sort(diag(D)) ; % 固有値の大きさ順に並べ替え
rangeW = 1:N ;
scales = sqrt(puiss(rangeW)); % 固有値の平方根
W = diag(1./scales) * U(1:N,k(rangeW))';
% W=固有値の平方根で規格化＋固有値基底への変換を行うホワイトニング行列
X = W*X; % ホワイトニング実行
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%% 学習 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 分離行列の初期条件
B=eye(N); %単位行列に選ぶ
L=0.1; %学習率
% Bell & Sejnowsky のアルゴリズム
% tanh(u) の部分は他の非線形関数でもよい
% for t=1:P
% u=B*X(:,t);
% B=B+L*(Id -tanh(u)*u')*B;
% end;
% EPS で収束条件を設定 shiro
EPS =0.005;
DELTA=eye(N);
while norm(DELTA,'fro')>EPS,
u =B*X;
tanu=tanh(u);
DELTA=L*(-(tanu*u'-u*tanu')/P)*B;
B=B+DELTA;
end;
% ホワイトニングを含む行列に変換
B = B*W ;
return;