Is my implementation of FastICA right?

Question

I have written the following code to perform ICA on the mixed signals that were generated using a random mixing matrix. But the performance of the algorithm is not satisfactory. Did I do something wrong?

[x1, Fs1] = wavread('source2.wav');
[x2, Fs2] = wavread('source3.wav');
m = size(x1,1); % size of each signal
n = 2; % Number of sound sources
A = randn(n, n); %  Random 2 X 2 mixing matrix
x = A*[x1';x2']; % Mixed signals
mx = sum(x, 2)/m; % supposed to be the mean; E{x}
x = x - repmat(mx, 1, m); % x - E{x}
for i = 1:n                      % This step is done
    xi = x(i,:);                 % to whiten the data
    [E, D] = eig(xi*xi');        % (make variance one
    xi = E*sqrtm(inv(D))*E'*xi;  % and uncorrelate the data)
    x(i,:) = xi;
end
w = randn(n, 1); % Random weight vector
w = w/norm(w,2); % make 'w' a unit vector
w0 = randn(n, 1);
w0 = w0/norm(w0, 2);
while abs(abs(w0'*w)-1) > 0.000001
    w0 = w;
    w = x*G(w'*x)'/m - sum(DG(w'*x))*w/m; % This step is supposed to perform:
                                          % w = E{xg(w^{T}*x)} - E{g'(w^{T}*x)}w
    w = w/norm(w, 2);
end
sound(w'*x); % Supposed to be one of the original signals

The algorithm can be found here http://en.wikipedia.org/wiki/FastICA

In the above code G is the function $xe^{\frac{-x^2}{2}}$ (computes element wise) and DG is the derivative of G: $(1-x^2)e^{\frac{-x^2}{2}}$ (computes element wise).

Answer 1

Though the algorithm itself is correct, the part of the code that whitens the data is not correct. All I had to do to make the covariances of the data $1$ was the following.

c = cov(x');
sq = inv(sqrtm(c));
x = c*xx;

And I changed the function $G$ from the one in the question to $\tanh(x)$. It works like a charm now.