Independent Component Analysis (ICA), Optional

by Seungchul Lee
iSystems Design Lab
http://isystems.unist.ac.kr/
UNIST

Table of Contents

1. The Cocktail Party Problem

You're at a crowded party. The music is loud, people are laughing, and a dozen different conversations are happening all around you. However, despite the hubbub, you're able to focus on the one voice you want to hear.

(No reason to shout. Though background noise can be distracting, the brain has the remarkable ability to track conversation and scale down unwanted noise. Courtesy of the National Archives)

Consider two conversations in a room that are happening simultaneously. How is it that the two different acoustic signals of converstaion one and two can be separated out?

$$ \begin{align*} x_1 &= a_{11}s_1 + a_{12}s_2\\ x_2 &= a_{21}s_1 + a_{22}s_2 \end{align*} $$
In [4]:
[S1,Fs1] = audioread('The_Phantom_Of_The_Opera.wav');
[S2,Fs2] = audioread('Barack_Obama.wav');

m = 300000;
S1 = mean(S1(1:m,:),2); % stereo -> mono
S2 = mean(S2(1:m,:),2); % stereo -> mono

SS = [S1 S2]';
A = [1/2, 2/3];     
mixed = A*SS;

sound(mixed,Fs1)
Out[4]:

2. ICA

Observed random vector $x$ is modelled by a linear variable model $$ x_i = \sum_{j=1}^{m} a_{ij}s_j, \qquad i = 1,\cdots, n$$

or in matrix form

$$ x= As $$

where

  • the 'mixing' matrix $A$ is constant
  • the $s_i$ are latent random variables called the independent components
  • Estimate both $A$ and $s_i$, observing only $x$

  • must assume:

    • $s_i$ are mutually independent
    • $s_i$ are non-gaussian
  • At first, we might simply say that the solution of this is trivial and it is given by

    $$ s = A^{-1}x$$ However, we do not know $A$ and the signal $s$.
  • PCA and ICA (pictorial explanation)
$$\text{original}\to$$ $$\text{mixed}\to$$ $$\text{whitened}$$
$$\text{whitened}\to$$ $$\quad \text{rotated (demixed)}$$

3. Signal Example

In [2]:
[sig, mixedsig] = demosig;
[icasig] = fastica(mixedsig);
Out[2]:
Number of signals: 4
Number of samples: 500
Calculating covariance...
Dimension not reduced.
Selected [ 4 ] dimensions.
Smallest remaining (non-zero) eigenvalue [ 0.0253466 ]
Largest remaining (non-zero) eigenvalue [ 5.96405 ]
Sum of removed eigenvalues [ 0 ]
[ 100 ] % of (non-zero) eigenvalues retained.
Whitening...
Check: covariance differs from identity by [ 1.25455e-14 ].
Used approach [ defl ].
Used nonlinearity [ pow3 ].
Starting ICA calculation...
IC 1 .........computed ( 9 steps ) 
IC 2 ......computed ( 6 steps ) 
IC 3 ......computed ( 6 steps ) 
IC 4 ..computed ( 2 steps ) 
Done.
Adding the mean back to the data.
In [3]:
subplot(4,1,1), plot(sig(1,:)')
subplot(4,1,2), plot(sig(2,:)')
subplot(4,1,3), plot(sig(3,:)')
subplot(4,1,4), plot(sig(4,:)')
Out[3]:
In [4]:
subplot(4,1,1), plot(mixedsig(1,:)')
subplot(4,1,2), plot(mixedsig(2,:)')
subplot(4,1,3), plot(mixedsig(3,:)')
subplot(4,1,4), plot(mixedsig(4,:)')
Out[4]:
In [5]:
subplot(4,1,1), plot(icasig(1,:)')
subplot(4,1,2), plot(icasig(2,:)')
subplot(4,1,3), plot(icasig(3,:)')
subplot(4,1,4), plot(icasig(4,:)')
Out[5]:

4. Sound Example

Download

  • The_Phantom_Of_The_Opera.wav
  • Barack_Obama.wav
In [1]:
[S1,Fs1] = audioread('The_Phantom_Of_The_Opera.wav');
[S2,Fs2] = audioread('Barack_Obama.wav');

m = 300000;
S1 = mean(S1(1:m,:),2);
S2 = mean(S2(1:m,:),2);
Out[1]:
In [2]:
SS = [S1 S2]';
A = [1/2, 2/3;
     4/5, 1/4];

mixed = A*SS;
Out[2]:
In [3]:
sound(mixed(1,:),Fs1)
Out[3]:
In [4]:
sound(mixed(2,:),Fs1)
Out[4]:
In [5]:
icasig = fastica(mixed);
Out[5]:
Number of signals: 2
Number of samples: 300000
Calculating covariance...
Dimension not reduced.
Selected [ 2 ] dimensions.
Smallest remaining (non-zero) eigenvalue [ 0.000294554 ]
Largest remaining (non-zero) eigenvalue [ 0.00426385 ]
Sum of removed eigenvalues [ 0 ]
[ 100 ] % of (non-zero) eigenvalues retained.
Whitening...
Check: covariance differs from identity by [ 7.10543e-15 ].
Used approach [ defl ].
Used nonlinearity [ pow3 ].
Starting ICA calculation...
IC 1 .......computed ( 7 steps ) 
IC 2 ..computed ( 2 steps ) 
Done.
Adding the mean back to the data.
In [6]:
sound(icasig(1,:),Fs1)
Out[6]:
In [7]:
sound(icasig(2,:),Fs1)
Out[7]:

5. Image Example

Download

  • Shine2.jpg
  • nature2.jpg
In [8]:
S1 = imread([pwd,'\data_files\Shine2.jpg']);
S2 = imread([pwd,'\data_files\nature2.jpg']);
Out[8]:
In [9]:
subplot(1,2,1), imshow(S1)
subplot(1,2,2), imshow(S2)
Out[9]:
In [10]:
A = [3/4 1/5;
    1/2 2/3];

X1 = double(A(1,1)*S1 + A(1,2)*S2);
X2 = double(A(2,1)*S1 + A(2,2)*S2);

subplot(1,2,1), imshow(uint8(X1));
subplot(1,2,2), imshow(uint8(X2));
Out[10]:
In [29]:
[m,n,k] = size(X1);
x1 = reshape(X1,m*n*k,1);
x2 = reshape(X2,m*n*k,1);
mixedimages = [x1 x2]';
icaimages = fastica(mixedimages);
icaimages = icaimages';

S1bar = reshape(icaimages(:,1),m,n,k);
S2bar = reshape(icaimages(:,2),m,n,k);

min1 = min(min(min(S1bar)));
S1bar = S1bar - min1;
max1 = max(max(max(S1bar)));
S1bar = S1bar*(255/max1);

min2 = min(min(min(S2bar)));
S2bar = S2bar - min2;
max2 = max(max(max(S2bar)));
S2bar = S2bar*(255/max2);

subplot(1,2,1), imshow(uint8(S1bar)), axis off
subplot(1,2,2), imshow(uint8(S2bar)), axis off