Page 233 - MATLAB Recipes for Earth Sciences

Recommended Reading 229

rand('state',1);
div = 0;
B = orth(rand(3, 3) - .5);
BOld = zeros(size(B));
while (1 - div) > eps
B = B * real(inv(B' * B)^(1/2));
div = min(abs(diag(B' * BOld)));
BOld = B;
B = (sPCA' * ( sPCA * B) .^ 3) / length(sPCA) - 3 * B;
sICA = sPCA * B;
end

We plot the separated components with (Fig. 9.6)
subplot(3,2,2), plot(sICA(:,1))
ylabel('s_{ICA1}'), title('Separated signals - ICA')
subplot(3,2,4), plot(sICA(:,2)), ylabel('s_{ICA2}')
subplot(3,2,6), plot(sICA(:,3)), ylabel('s_{ICA3}')
The PCA algorithm has not reliably separated the mixed signals. Especially
the saw-tooth signal was not correctly found. In contrast, the ICA has found
the source signals almost perfectly. The only remarkable differences are the
noise, which came through the observation, the wrong sign and the wrong
order of the signals. However, the sign and the order of the signals are not
really important, because we have in general not the knowledge about the
real sources nor their order. With

A_ICA = A_PCA * B;
W_ICA = B' * W_PCA;
we compute the mixing matrix A and the separation matrix W. The mix-
ing matrix A can be used in order to estimate the portion of the separated
signals on our measurements The components a of the mixing matrix A
i,j
correspond to the principal components loads as introduced in Chapter 9.2.
A FastICA package is available for MATLAB and can be found at

http://www.cis.hut.fi/projects/ica/fastica/

Recommended Reading

Aitchison J (1984) The statistical analysis of geochemical composition. Mathematical
Geology 16(6):531-564
Aitchison J. (1999) Logratios and Natural Laws in Compositional Data Analysis. Mathematical
Geology 31(5):563-580
Birks HJB, Gordon AD (1985) Numerical methods in Quaternary pollen analysis. Academic
Press, London