bsica.m

image separation using neural net

function [B] = bsica(X, k, varargin)
% BSICA    Use the Bell and Sejnowski ICA algorithm to learn the weights of 
% the matrix B for ICA.
%
% B = bsica(X, options)
%   X is the set of input vectors.  Each column of X is one sample.
%   k is the number of independent components to extract.
%   B is the seperating matrix, Y = BX are the independent signals
%
% options (specified by key/value pairs)
%   'rate' = 0.01       the learning rate (initial)
%   'B' = rand          the initial seperating matrix
%   'niter' = 1         the number of iterations to run
%   'g' = tanh(t)       Bell and Senjowski's independence forcing function

%
% David Gleich
% CS 152 - Neural Networks
% 12 December 2003
%

% dimensionality assessment
[m n] = size(X);
if (k > m)
    error('Cannot extra more sources than sensors.');
end;

% subtract the mean
X = X - repmat(mean(X,2), 1, n);

% options
options = struct(...
    'rate', 0.01, ...
    'B', rand(k,m), ...
    'niter', 1, ...
    'g', inline('tanh(t)'));
options = getopt(options, varargin);

B = options.B;
u = options.rate;
g = options.g;

success = 0;
mdS = 0;
iter = 1;

for niter=1:options.niter
    Bold = B;
    
    % update B, S
    for ii = 1:n
        y = B*X(:, ii);
        
        B = B + u*(eye(k,k) - 2*g(y)*y')*B;
        
        if (any(~isfinite(B)))
            warning(sprintf('Lost convergence at iterator %i; lower learning rate?', iter));
            success = 11;
            break;
        end;
        
        iter = iter + 1;
    end;
   
    
    if (success == 1)
        break;
    end;
end;