How to do Independent Component Analysis?

Question

The new experimental function in 10.3.1 DimensionReduce[], has the following three options for Method:

"PrincipalComponentsAnalysis"
"LatentSemanticAnalysis"
"LowRankMatrixFactorization"

What about Independent Component Analysis (ICA) - is that already implemented somewhere, or am I missing it?

If it is not currently available, there are at least two packages in R for ICA which can be accessed from Mathematica. — JimB, Feb 10 '16 at 04:57
You could easily link to the FastICA algorithm with either RLink or MATLink. — dr.blochwave, Feb 10 '16 at 15:21
It seems to me that the question and accepted answer of Disentangling the data are related. To be clear, Non-Negative Matrix Factorization (NNMF) does dimension reduction, but its norm minimization process does not enforce variable independence. (It enforces non-negativity.) — Anton Antonov, Feb 10 '16 at 16:23
I think it is a good idea to add to the question an example of data for which ICA would give results that help interpreting it. E.g. the "cocktail party problem". — Anton Antonov, Feb 10 '16 at 16:29

Anton Antonov · Accepted Answer · 2016-02-11T00:07:55.363

It seems that Non-Negative Matrix Factorization (NNMF) can be applied for doing ICA. At least in some cases.

In order to demonstrate this I will make up some data in the spirit of the "cocktail party problem". Then I am going to apply an NNMF algorithm.

To be clear, NNMF does dimension reduction, but its norm minimization process does not enforce variable independence. (It enforces non-negativity.) There are at least several articles discussing modification of NNMF to do ICA. For example this one: "A new nonnegative matrix factorization for independent component analysis". (From it I took the data generation formulas.)

Data

(*Signal functions*)
Clear[s1, s2, s3]
s1[t_] := Sin[600 \[Pi] t/10000 + 6*Cos[120 \[Pi] t/10000]] + 1.2
s2[t_] := Sin[\[Pi] t/10] + 1.2
s3[t_?NumericQ] := (((QuotientRemainder[t, 23][[2]] - 11)/9)^5 + 2.8)/2 + 0.2

(*Mixing matrix*)
A = {{0.44, 0.2, 0.31}, {0.45, 0.8, 0.23}, {0.12, 0.32, 0.71}};

(*Signals matrix*)
nSize = 600;
S = Table[{s1[t], s2[t], s3[t]}, {t, 0, nSize, 0.5}];

(*Mixed signals matrix*)
M = A.Transpose[S];

(*Signals*)
Grid[{Map[
   Plot[#, {t, 0, nSize}, PerformanceGoal -> "Quality", 
     ImageSize -> 250] &, {s1[t], s2[t], s3[t]}]}]

(*Mixed signals*)
Grid[{Map[ListLinePlot[#, ImageSize -> 250] &, M]}]

Application of Non-Negative Matrix Factorization (NNMF)

Load NNMF package (from MathematicaForPrediction at GitHub):

Import["https://raw.githubusercontent.com/antononcube/MathematicaForPrediction/master/NonNegativeMatrixFactorization.m"]

After several applications of NNMF we get signals close to the originals:

{W, H} = GDCLS[M, 3];
Grid[{Map[ListLinePlot[#, ImageSize -> 250] &, Normal[H]]}]

score 8 · Answer 2 · edited Apr 13 '17 at 12:56

The package IndependentComponentAnalysis.m can be used for Independent Component Analysis (ICA).

This answer uses the generated data from my previous answer (which is about opportunistic application of general Non-Negative Matrix Factorization for ICA).

Load the package:

Import["https://raw.githubusercontent.com/antononcube/MathematicaForPrediction/master/IndependentComponentAnalysis.m"]

It is important to note that the usual ICA model interpretation for the factorized matrix $X$ is that each column is a variable (audio signal) and each row is an observation (recordings of the microphones at a given time). The matrix $M \in R^{3 \times 1201}$ was constructed with the interpretation that each row is a signal, hence we have to transpose $M$ in order to apply the ICA algorithms, $X=M^T$.

X = Transpose[M];

{S, A} = IndependentComponentAnalysis[X, 3];

Check the approximation of the obtained factorization:

Norm[X - S.A]    
(* 3.10715*10^-14 *)

Plot the found source signals:

Grid[{Map[ListLinePlot[#, PlotRange -> All, ImageSize -> 250] &, 
   Transpose[S]]}]

Because of the random initialization of the inverting matrix in the algorithm the result my vary. Here is the plot from another run:

The package also provides the function FastICA that returns an association with elements that correspond to the result of the function fastICA provided by the R package "fastICA". See

Here is an example usage:

res = FastICA[X, 3];

Keys[res]    
(* {"X", "K", "W", "A", "S"} *)

Grid[{Map[
   ListLinePlot[#, PlotRange -> All, ImageSize -> Medium] &, 
   Transpose[res["S"]]]}]

Note that (in adherence to the cited documents) the function FastICA returns the matrices S and A for the centralized matrix X. This means, for example, that in order to check the approximation proper mean has to be supplied:

Norm[X - Map[# + Mean[X] &, res["S"].res["A"]]]
(* 2.56719*10^-14 *)

How to do Independent Component Analysis?

2 Answers2

Further readings

Linked