Function Repository Resource:

IndependentComponentAnalysis

Source Notebook

Decompose a matrix into Independent Component Analysis matrix factors

Contributed by: Anton Antonov

ResourceFunction["IndependentComponentAnalysis"][mat,k]

decomposes the matrix mat into k components.

ResourceFunction["IndependentComponentAnalysis"][mat,k,opts]

decomposes mat using the options opts.

Details and Options

The Independent Component Analysis (ICA) is a matrix factorization method that utilizes Gram-Schmidt orthogonalization by considering two vectors (two signals) to be orthogonal if their difference is Gaussian white noise.

The argument mat can be a square or rectangular matrix.

The columns of mat are interpreted as signals.

ResourceFunction["IndependentComponentAnalysis"] accepts the following options:

"InitialUnmixingMatrix"

Automatic

initial un-mixing matrix

"NegEntropyFactor"

negative entropy factor

"NonGaussianityFunction"

Automatic

non-Gaussianity measure function

"RowNorm"

False

whether to normalize the matrix rows at each step

MaxSteps

200

maximum steps of the ICA interation process

PrecisionGoal

precision goal

For some data, at least several executions are required in order to get components that provide satisfactory interpretation.

The negative entropy factor is specified with the option "NegEntropyFactor". The negative entropy factor controls how aggressively ICA learns from the "landscape" produced by the non-Gaussianity measure function (specified with "NonGaussianityFunction").

The function ResourceFunction["IndependentComponentAnalysis"] uses SingularValueDecomposition internally.

Examples

Basic Examples (3)

Here is a random integer matrix:

In[1]:=

SeedRandom[7]
mat = RandomInteger[10, {4, 3}];
MatrixForm[mat]

Out[3]=

Here are the Independent Component Analysis matrix factors:

In[4]:=

Out[4]=

Here is the matrix product of the obtained factors:

In[5]:=

Out[8]=

Scope (6)

Here is a random matrix with its first two columns having much larger magnitudes than the rest:

In[9]:=

SeedRandom[232];
mat = RandomReal[{10, 100}, {10, 2}];
mat2 = RandomReal[{-10, 10}, {10, 7}];
mat = ArrayPad[mat, {{0, 0}, {0, 5}}] + mat2;

Plot the values:

In[10]:=

Out[10]=

Compute the ICA matrix factors:

In[11]:=

Here is the relative error of the approximation by the obtained matrix factorization:

In[12]:=

Out[12]=

Here is the relative error for the first three columns:

In[13]:=

Out[13]=

Here are comparison plots:

In[14]:=

Block[{opts = {PlotRange -> All, ImageSize -> 250, PlotTheme -> "Detailed"}},
Row[{ListLinePlot[Transpose[mat], opts, PlotLabel -> "Orginal"], ListLinePlot[Transpose[A . S], opts, PlotLabel -> "Reconstructed"]}]
]

Out[14]=

Options (4)

InitialUnmixingMatrix (2)

Using the option "InitialUnmixingMatrix", we can influence the results by providing initial un-mixing matrix for ICA's iteration process. Here we compute ICA with different initial un-mixing matrices:

In[15]:=

k = 3;
SeedRandom[8];
res = Association[
MapThread[
Row[{#2, ":", Spacer[2], If[MatrixQ[#], MatrixForm[#1], ToString[#1]]}] -> ResourceFunction["IndependentComponentAnalysis"][mat, k, "InitialUnmixingMatrix" -> #1] &, {{Automatic, Automatic, RandomVariate[NormalDistribution[1, 10], {k, k}], SingularValueDecomposition[mat, k][[2]]}, Range[4]}]];

Plot the results:

In[16]:=

Grid[KeyValueMap[{#1, ListLinePlot[Transpose[#2[[1]]], PlotRange -> All, ImageSize -> 250, PlotTheme -> "Detailed"]} &, res], Alignment -> {Center, Left}]

Out[16]=

NonGaussianityFunction (2)

ICA is essentially a Gram-Schmidt process that considers two signals to be orthogonal if their difference is Gaussian noise.

Here we compute ICA with different non-Gaussianity measure functions:

In[17]:=

$res = Association[ # -> BlockRandom[ SeedRandom[22]; ResourceFunction["IndependentComponentAnalysis"][mat, 3, "NonGaussianityFunction" -> #, MaxSteps -> 2] ] & /@ {Log[ Cosh[#]] &, -Exp[-#^2/2] &, (Log[Cosh[1.2 #]]/1.2 &)}];$

Plot the results:

In[18]:=

Grid[KeyValueMap[{#1[x], ListLinePlot[Transpose[#2[[1]]], PlotRange -> All, ImageSize -> 250, PlotTheme -> "Detailed"]} &, res], Alignment -> {Center, Left}]

Out[18]=

Note that we used a small number of steps in order to obtain similar but different results. Using a sufficiently large number of steps (say, 20) we get almost the same results (with the matrix mat).

Applications (8)

Blind signal separation (8)

Here are signal generating functions:

In[19]:=

$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$

Here is a mixing matrix:

In[20]:=