ICA - Statistical Independence & Eigenvalues of Covariance Matrix

Question

I am currently creating different signals using Matlab, mixing them by multiplying them by a mixing matrix A, and then trying to get back the original signals using FastICA.

So far, the recovered signals are really bad when compared to the original ones, which was not what I expected.

I'm trying to see whether I'm doing anything wrong. The signals I'm generating are the following:

s1 = (-x.^2 + 100*x + 500) / 3000; % quadratic
s2 = exp(-x / 10); % -ve exponential
s3 = (sin(x)+ 1) * 0.5; % sine
s4 = 0.5 + 0.1 * randn(size(x, 2), 1); % gaussian
s5 = (sawtooth(x, 0.75)+ 1) * 0.5; % sawtooth

One condition for ICA to be successful is that at most one signal is Gaussian, and I've observed this in my signal generation.

However, another condition is that all signals are statistically independent.

All I know is that this means that, given two signals A & B, knowing one signal does not give any information with regards to the other, i.e.: P(A|B) = P(A) where P is the probability.

Now my question is this: Are my signals statistically independent? Is there any way I can determine this? Perhaps some property that must be observed?

Another thing I've noticed is that when I calculate the eigenvalues of the covariance matrix (calculated for the matrix containing the mixed signals), the eigenspectrum seems to show that there is only one (main) principal component. What does this really mean? Shouldn't there be 5, since I have 5 (supposedly) independent signals?

For example, when using the following mixing matrix:

A =

0.2000    0.4267    0.2133    0.1067    0.0533
0.2909    0.2000    0.2909    0.1455    0.0727
0.1333    0.2667    0.2000    0.2667    0.1333
0.0727    0.1455    0.2909    0.2000    0.2909
0.0533    0.1067    0.2133    0.4267    0.2000

The eigenvalues are: 0.0000 0.0005 0.0022 0.0042 0.0345 (only 4!)

When using the identity matrix as the mixing matrix (i.e. the mixed signals are the same as the original ones), the eigenspectrum is: 0.0103 0.0199 0.0330 0.0811 0.1762. There still is one value much larger than the rest..

Thank you for your help.

I apologise if the answers to my questions are painfully obvious, but I'm really new to statistics, ICA and Matlab. Thanks again.

EDIT

I have 500 samples of each signal, in the range [0.2, 100], in steps of 0.2, i.e. x = 0:0.1:100.

Also, given the ICA Model: X = As + n (I'm not adding any noise at the moment), I am referring to the eigenspectrum of the transpose of X, i.e. eig(cov(X')).

UPDATE

As suggested (refer to comments), I tried FastICA on only 2 signals. The results were quite good (see pic below). The mixing matrix used was A = [0.75 0.25; 0.25 0.75]. However, the eigenspectrum 0.1657 0.7732 still showed only one main principal component.

My question therefore boils down to the following: What function/equation/property can I use to check whether a number of signal vectors are statistically independent?

Answer 1

Signals 3 and 5 appear to be quite correlated -- they share their first harmonic. If I were given two mixtures of those, I wouldn't be able to separate those, I'd be tempted to put the common harmonic as one signal and the higher harmonics as a second signal. And I would be wrong! This might explain the missing eigenvalue.

Signals 1 and 2 do not look independent too.

A quick and dirty "sanity-check" for the independence of two series is to do a (x, y) plot of one signal against the other:

plot (sig3, sig5)

and then to do the same (x, y) plot with one signal shuffled:

indices = randperm(length(sig3))
plot(sig3(indices), sig5)

If the two plots have a dissimilar look, your signals are not independent. More generally, if the (x, y) plot of the data shows "features", disymmetries etc, it's a bad omen.

Proper tests for independence (and those are the objective functions used in the ICA optimization loop) include, for example, mutual information.

ICA is recovering the most independent signals a linear mixing of which yields your input data. It will work as a signal separation methods and recover the original signals only if those were maximally independent according to the optimization criterion used in your ICA implementation.

Answer 2

I'm not an expert on ICA, but I can tell you a little bit about independence.

As some of the comments have mentioned, statistical independence between two random variables can be roughly interpreted as "the amount of information that observing one variable gives about another".

However, because we're talking about statistical independence, all of this is with respect to a distribution. In particular, if you observe two random variables $X$ and $Y$, independence of $X$ and $Y$ is equivalent to a condition on their joint distribution $p(x,y)$. In particular, $X$ and $Y$ are independent if and only if $p(x,y) = p(x)p(y)$.

This condition is basically what an earlier comment suggested that you use. Plotting a scatter plot of the two variables is a kind of visualization of the joint distribution $p(x,y)$. When two variables are independent, permuting the observations should give the same underlying joint distribution.

There are lots of ways to test for independence. In practice (say, to draw scientific conclusions), you might like a careful hypothesis test that allows you to state your assumptions and the confidence of your conclusions (like the chi-square test). I'll ignore that here, and assume you have perfect information about your distributions. Then, note that the mutual information between $X$ and $Y$ is zero exactly when $X$ and $Y$ are independent. To compute the mutual information, let your joint distribution be $p(X=i, Y=j) = p_{ij}$, and $P(X=i) = p_i$ and $P(Y=j) = p_j$ be your marginals. Then,

$I(X,Y) = \sum_i \sum_j p_{ij} \log \frac{p_{ij}}{p_i p_j}$

Here's some matlab code that will generate two independent signals from a constructed joint distribution, and two from a non-independent joint distribution, and then compute the mutual information of the joints.

The function "computeMIplugin.m" is a simple function I wrote that computes the mutual information using the summation formula above.

Ndist = 25;
xx = linspace(-pi, pi, Ndist);

P1 = abs(sin(xx)); P2 = abs(cos(xx)); 
P1 = P1/sum(P1); P2 = P2/sum(P2); % generate marginal distributions

%% Draw independent samples.
Nsamp = 1e4;
X = randsample(xx, Nsamp, 'true', P1);
Y = randsample(xx, Nsamp, 'true', P2);

Pj1 = P1'*P2;
computeMIplugin(Pj1)

% I get approx 8e-15 ... independent!

% Now Sample the joint distribution 
cnt = {}; cnt{1} = xx; cnt{2} = xx; % bin centers
Pj1_samp= hist3([X' Y'],cnt); Pj1_samp = Pj1_samp/sum(Pj1_samp(:));
computeMIplugin(Pj1_samp)
% I get approx .02; since we've estimated the distribution from
% samples, we don't know the true value of the MI. This is where
% a confidence interval would come in handy. We'd like to know 
% whether value of MI is significantly different from 0. 

% mean square difference between true and sampled?
% (this is small for these parameter settings... 
% depends on the sample size and # bins in the distribution).
mean( (Pj1_samp(:) - Pj1(:)).^2)

%% Draw samples that aren't independent. 

tx = linspace(0,30,Nsamp);
X = pi*sin(tx);
Y = pi*cos(tx);

% estimate the joint distribution
cnt = {}; cnt{1} = xx; cnt{2} = xx; % bin centers
Pj2= hist3([X' Y'],cnt); Pj2 = Pj2/sum(Pj2(:));
computeMIplugin(Pj2)

% I get 1.9281  - not independent!

%% make figure
figure(1); 
colormap gray
subplot(221)
imagesc(xx,xx,Pj1_samp)
title('sampled joint distribution 1')
subplot(222)
imagesc(xx,xx,Pj2)
title('sampled joint distribution 2')
subplot(223)
imagesc(xx,xx,Pj1)
title('true joint distribution 1')

Again, this assumes you have a good estimate of the joint distribution (along with other blithe assumptions), but it should be useful as a rule of thumb.

Answer 3

As mentioned above, both signals 3 and 5 appear to be quite correlated and have a similar period.

We can think of two signals being correlated if we can shift one of the sources left or right and increase or decrease its amplitude so that it fits on top of the other source. Note we are not changing the frequency of the source, we are just performing a phase and amplitude shift.

In the above case, we can shift source 3 so that its peaks coincide with source 5. This is the sort of thing that will mess source extraction when using ICA due to the independence assumption.

Note: A nice illustration of the above concept is to think of two sinusoidal waves. These are both completely deterministic. If they both have the same frequency (even with different phase) then they are perfectly correlated and ICA will not be able to separate them. If instead they have different frequencies (that are not integer multiples of one another), then they are independent and can be seperated.

Below is some Matlab code for you to see this for yourself

%Sine waves of equal frequency
X = 1:1000;
Y(1,:) = sin(2*pi*X*10/1000);
Y(2,:) = sin(1+2*pi*X*10/1000);

figure
subplot(3,2,1)
plot(Y(1,:))
title('Initial Source 1')
subplot(3,2,2)
plot(Y(2,:))
title('Initial Source 2')
A = [1, 2; 4, -1];
Y = A*Y;
subplot(3,2,3)
plot(Y(1,:))
title('Signal 1')
subplot(3,2,4)
plot(Y(2,:))
title('Signal 2')

Z = fastica(Y);

subplot(3,2,5)
plot(Z(1,:))
title('Source 1')
subplot(3,2,6)
plot(Z(2,:))
title('Source 2')

%Sine waves of different frequency
X = 1:1000;
Y(1,:) = sin(2*pi*X*10/1000);
Y(2,:) = sin(1+2*pi*X*8/1000);

figure
subplot(3,2,1)
plot(Y(1,:))
title('Initial Source 1')
subplot(3,2,2)
plot(Y(2,:))
title('Initial Source 2')
A = [1, 2; 4, -1];
Y = A*Y;
subplot(3,2,3)
plot(Y(1,:))
title('Signal 1')
subplot(3,2,4)
plot(Y(2,:))
title('Signal 2')

Z = fastica(Y);

subplot(3,2,5)
plot(Z(1,:))
title('Source 1')
subplot(3,2,6)
plot(Z(2,:))
title('Source 2')

Note that for the waves of the same frequency, ICA just returns the input signals, but for different frequencies it returns the original sources.

Answer 4

Rachel,

From my research I have so far been able to find something called 'Chi-Squared Test For Independence', but I am not sure how it works at the moment, but it might be worth a look.