FastICA example

Import the spectrochempy API package

import spectrochempy as scp

Independent component analysis (ICA) is a computational method for separating a multivariate signal such as spectra into additive components. This is done by assuming that at most one subcomponent is Gaussian and that the components are statistically independents from each other.

# Load, prepare and plot the dataset
# ----------------------------------
# Here we use a dataset from :cite:t:`jaumot:2005`

X = scp.read("matlabdata/als2004dataset.MAT")[-1]

X.title = "absorbance"
X.units = "absorbance"
X.y = scp.Coord.arange(
    X.shape[0], dtype="float", title="elution time"
)  # floats to trigger sequential colormap
X.y.units = "min"
X.x.title = "wavelength"
_ = X.plot()

Create and fit a FastICA object

As argument of the object constructor we define log_level to "INFO" to obtain verbose output during fit, and we set the number of component to use at 4.

ica = scp.FastICA(n_components=4, log_level="INFO")
_ = ica.fit(X)

Get the mixing system and source spectral profiles

The mixing system \(A\) and the source spectral profiles \(S^T\) can be obtained as follows (the Sklearn equivalents - also valid with Scpy - are indicated as comments

A = ica.A  # or model.transform()
St = ica.St  # or model.mixing.T

Plot them

_ = A.T.plot(title="Mixing System", colormap=None)
_ = St.plot(title="Sources spectral profiles", colormap=None)

Reconstruct the dataset

The dataset can be reconstructed from these matrices and the mean:

X_hat_a = scp.dot(A, St) + X.mean(dim=0).data
_ = X_hat_a.plot(title=r"$\hat{X} = \bar{X} + A S^t$")

Or using the transform() method:

X_hat_b = ica.inverse_transform()
_ = X_hat_b.plot(title=r"$\hat{X} =$ ica.inverse_transform()")

Finally, the quality of the reconstriction can be checked by plot_merit()

_ = ica.plot_merit(nb_traces=15)

Uncomment the following line to display all figures when running the script directly with Python.

# scp.show()