✦ Day 06 · Signal Separation & Geometry

Independent Component
Analysis (ICA)

Solving the cocktail party problem from first principles: Stone's geometric mixing model, the Central Limit Theorem, state-space projection, and the two-step whitening-rotation pipeline.

🍸 Cocktail Party Problem 📐 Stone's Geometric Model ⚡ Quasi-Static Maxwell 🔔 Central Limit Theorem 🔀 3 Indeterminacies 🎯 PICARD Algorithm 📏 Makeig $k \cdot N^2$ Rule 🐍 mne.preprocessing.ICA
1

The Cocktail Party Intuition
Physical Independence & Blind Source Separation (Stone Ch. 1)

Imagine standing in the center of a crowded room where two people are speaking simultaneously. You place two microphones at different locations. Each microphone records a mixture of both voices, but with different proportions because of the physical distance from each speaker:

$$\begin{aligned} \text{Microphone 1: } \quad x_1(t) &= a \cdot s_1(t) + b \cdot s_2(t) \\ \text{Microphone 2: } \quad x_2(t) &= c \cdot s_1(t) + d \cdot s_2(t) \end{aligned}$$

Here, $s_1(t)$ and $s_2(t)$ are the pure source signals (the individual voices), and $x_1(t), x_2(t)$ are the measured mixtures. In matrix notation: $$\mathbf{x}(t) = \mathbf{A}\,\mathbf{s}(t)$$ where $\mathbf{A} = \begin{pmatrix} a & b \\ c & d \end{pmatrix}$ is the unknown mixing matrix.

Why are the original voices independent?
The amplitude of Speaker 1 at any millisecond is physically unrelated to Speaker 2's vocal chords. They are governed by two separate brains and physical processes. Therefore, knowing the instantaneous amplitude of $s_1(t)$ tells you nothing about $s_2(t)$. This is **Statistical Independence**.

🧠 Biophysical Grounding: Why Does Linear Mixing $\mathbf{x}(t) = \mathbf{A}\,\mathbf{s}(t)$ Hold in the Human Brain?

In an actual recording setup, microphones are scalp electrodes, and speakers are cortical neuronal populations and physiological artifact generators. The linear instantaneous model is physically validated by three biological realities:

  • ✦ 1. Pyramidal Columnar Dipole Generators: Layer V pyramidal neurons in the cerebral cortex are arranged in strict parallel columns perpendicular to the cortical sheet. Synchronous postsynaptic potentials (PSPs) across thousands of adjacent neurons summate constructively into an Equivalent Current Dipole (ECD) with vector moment $\mathbf{d} = I \cdot \mathbf{l}$.
  • ✦ 2. Functional Autonomy of Cortical Areas: Distant functional areas (e.g. primary visual V1, auditory A1, motor M1, and frontal eye fields) have independent local microcircuitry and incoming synaptic drive. Their instantaneous local field potential (LFP) fluctuations are statistically independent across time.
  • ✦ 3. Quasi-Static Electromagnetic Linearity: In biological tissues at EEG frequencies ($0.1\text{–}100\text{ Hz}$), tissue displacement currents are negligible ($\epsilon \partial \mathbf{E}/\partial t \approx 0$). Maxwell's equations reduce to Laplace's equation $\nabla \cdot (\sigma \nabla V) = 0$. Volume conduction across the cerebrospinal fluid, skull, and scalp to electrodes is strictly instantaneous, purely linear, and free of frequency-dependent phase delays.

Why PCA Fails Here

Principal Component Analysis (PCA) looks for orthogonal axes of maximum variance. PCA decorrelates signals (removes 2nd-order dependencies), but orthogonal mixtures are still mixtures! If the physical mixing matrix $\mathbf{A}$ sheared the coordinate axes non-orthogonally, PCA will rotate the data into an artificial orthogonal coordinate frame where the sources remain hopelessly blended.

Why ICA Succeeds

Independent Component Analysis (ICA) does not restrict unmixing axes to be orthogonal. Instead, it searches across all possible coordinate shears and rotations to find directions where the recovered signals are statistically independent across all higher-order statistical moments.

Interactive Cocktail Party Simulator

Adjust the microphone mixing coefficients ($a, b, c, d$) and see how pure voices blend into mixtures and are reconstructed via unmixing matrix $\mathbf{W} = \mathbf{A}^{-1}$.

Live 2×2 BSS

Pure Sources

s₁(t) Oscillation s₂(t) Pulse

Microphone Mixtures

x₁(t) = a·s₁ + b·s₂ x₂(t) = c·s₁ + d·s₂

Recovered Signals

ŝ₁(t) Recovered 1 ŝ₂(t) Recovered 2
2

The Visual Intuition & Geometry of Mixing
Three Hallmarks of Signal Mixtures (Stone Ch. 2 & 3)

How does ICA separate mixed brain signals without knowing what the original brainwaves looked like? In Chapters 2 & 3 of Independent Component Analysis, James V. Stone reveals that mixing signals produces three simple, visual geometric effects. Once you see them visually, the entire algorithm clicks:

1

Hallmark 1: Independence Looks Like a Cross — Mixing Tilts It into a Diamond

When two physical processes are completely independent, knowing the value of one tells you nothing about the other. If you plot them on a 2D graph, the data points naturally form a crisp, symmetrical plus sign (+) or box aligned with the horizontal and vertical axes.

👁️ The Tale of Two Signals: Occipital Alpha vs. Eye Blink

Imagine two distinct things happening in a subject's head at the same time:

  • The Brainwave (Alpha): A gentle, continuous oscillation coming from the visual cortex at the back of the head. It hums back and forth like a pendulum.
  • The Artifact (Eye Blink): A giant electrical surge from the front of the head. 95% of the time it is completely silent (zero voltage), but once every few seconds it spikes with massive amplitude.
✨ In Pure Source Space: When eyes are still, alpha traces a horizontal line. When a blink fires, it shoots straight up. Together, they form an unmistakable orthogonal cross (+).
🌀 On the Scalp Electrodes: Because head tissues conduct electricity everywhere, both sensors record a blend of both. The crisp cross gets sheared and squished into a tilted diamond.
2

Hallmark 2: The "Smoothie" Effect — Mixtures Always Look Like a Bell Curve

This is the famous Central Limit Theorem, and it is ICA's golden compass. If you mix several distinct non-Gaussian signals together, their combined sum always becomes smoother, rounder, and more bell-shaped (more Gaussian) than any individual ingredient!

🥤 The Fruit Smoothie Analogy

Think of independent signals as fresh fruits with sharp, unmistakable flavors — super-sour lemon (spiky eye blink) and sweet strawberry (steady alpha).

  • In the Blender: Once you puree them together, the distinct sharp flavors blend into a uniform, generic smoothie (a smooth Gaussian bell curve).
  • The Unmixing Mission: At any scalp electrode, you are recording a blended smoothie of 10+ neural and artifact sources.
  • The ICA Trick: Instead of guessing what the brain is doing, ICA simply asks: "In which direction does the signal taste least like a smoothie and most like a pure single fruit?" Maximizing non-Gaussianity automatically recovers the pure source!
3

Hallmark 3: Unmixing by Looking "Edge-On" — Annihilating Unwanted Sources

How do you remove a huge artifact without knowing its shape? You look at it from an angle where it has zero width! In geometry, this is called an orthogonal projection.

📏 The Sheet of Paper Analogy

Hold a large sheet of cardboard in front of you. Viewed face-on, it completely blocks your view of whatever is behind it (like an eye blink drowning out tiny microvolt brainwaves).

  • If you turn your head 90 degrees perpendicular to the cardboard, you are now looking directly down its thin edge. The cardboard becomes virtually invisible!
  • ICA does the exact same thing mathematically: it finds a digital viewing angle that is perpendicular to the eye blink's direction vector. Along that line of sight, the blink is multiplied by zero — completely vanishing from the recording while leaving the underlying brainwaves untouched.
🎭

The Three Signal Personalities (Stone Ch. 2 & 3)

ICA classifies signals not by their frequency, but by their statistical personality — how spiky or flat their distribution is compared to a standard bell curve:

🎵

The Steady Hum
Sub-Gaussian (Flat / Negative Kurtosis)

Personality: Never rests. It continuously swings smoothly between minimum and maximum values without hanging around zero.

Everyday Examples:
• 10 Hz Alpha Rhythm: Smooth rhythmic sine wave.
• 50/60 Hz Mains Hum: Constant electrical line buzz.
• Synthesizer Tones: Clock ticks & square waves.
2D Shape: Flat Square / Solid Box
⚡

The Silent Lightning Bolt
Super-Gaussian (Spiky / Positive Kurtosis)

Personality: Quiet 90% of the time, punctuated by sudden, explosive bursts of high amplitude.

Everyday Examples:
• Eye Blinks (EOG): Silent baseline + massive pulses.
• Muscle Twitches (EMG): High-voltage tension bursts.
• Human Speech: Words separated by quiet pauses.
2D Shape: Sharp Cross / 4-Point Star
🌫️

The Featureless Fog
Gaussian (Zero Kurtosis — The Blind Spot)

Personality: Completely uniform in all directions, like a perfectly round disk. No edges or sharp corners.

Everyday Examples:
• Thermal Hardware Noise: Random wire electron jitter.
• Blended Background Brain Noise: Sum of thousands of tiny distant synaptic potentials.
2D Shape: Round Circular Disk (No Angle Preference)

Why ICA Fails on Pure Gaussians: Because a circular disk looks identical no matter how much you spin it, ICA cannot find a preferred angle. But the moment even one source is spiky (super-Gaussian) or box-like (sub-Gaussian), the disk develops corners — giving ICA a clear compass needle to point towards!

Signal / Artifact What It Looks Like Personality Can ICA Clean It?
Eye Blink (EOG) Huge 150–300 µV frontal upward peaks every few seconds ⚡ Spiky (Super-Gaussian) ✔ Cleaned effortlessly
Jaw / Muscle Tension (EMG) Fast, scratchy high-frequency electrical bursts ⚡ Spiky (Super-Gaussian) ✔ High isolation
Heartbeat (ECG) Rhythmic periodic sharp pulses across electrodes ⚡ Spiky (Super-Gaussian) ✔ High isolation
Occipital Alpha (10 Hz) Smooth sinusoidal waves when eyes are closed 🎵 Steady (Sub-Gaussian) ✔ Separated cleanly
Mains Hum (50/60 Hz) Continuous, persistent line-noise oscillation 🎵 Steady (Sub-Gaussian) ✔ Separated cleanly
Thermal Sensor Noise Random microvolt background floor 🌫️ Fog (Gaussian) ✘ ICA Blind Spot

2D State-Space Geometry & Non-Gaussianity Sandbox

Experience Stone's geometric model (Chapters 2 & 3) and the Central Limit Theorem. Click & drag on the canvas to rotate projection vector $\mathbf{w}$, adjust mixing shear, and watch the projected histogram transform from a blended Gaussian bell curve into pure independent source distributions!

Stone Ch. 2–3 Interactive

2D Joint State Space & Vector $\mathbf{w}$

Points $(x_1, x_2)$ Vector $\mathbf{w}$ (Drag on canvas) Source Axes

Projected Signal $y = \mathbf{w}^T\mathbf{x}$ vs. Gaussian Reference

Empirical Histogram Gaussian Fit $\mathcal{N}(\mu, \sigma^2)$

0% = Pure Unmixed Sources · 100% = Maximum Non-Orthogonal Shear

Excess Kurtosis: 0.00 (Gaussian Mixture)
CLT in action: Arbitrary projection blends sources into a Gaussian bell curve. Rotate $\mathbf{w}$ to isolate maximum non-Gaussianity!

Intuition Summary (Stone Ch. 2 & 3): When $\mathbf{w}$ points at an arbitrary angle, it captures a mix of both sources $\to$ the histogram matches the theoretical Gaussian bell curve (low non-Gaussianity). When $\mathbf{w}$ aligns perpendicular to an unwanted source axis, the projection isolates a pure source and maximally deviates from Gaussianity!

⚠️ The Three Inescapable Ambiguities (Indeterminacies) of Blind Source Separation

Because ICA reconstructs sources without knowing the mixing matrix $\mathbf{A}$ or the original sources $\mathbf{s}(t)$, there are three mathematical indeterminacies inherent to all BSS algorithms:

1. Permutation (Order) Indeterminacy

ICA has no natural concept of "Component 1" vs "Component 2". Any permutation matrix $\mathbf{P}$ yields equally independent components ($\mathbf{P}\mathbf{W}\mathbf{x} = \mathbf{P}\mathbf{s}$). Software packages (MNE, EEGLAB) sort components post-hoc by percent variance accounted for ($PVAF$).

2. Polarity (Sign) Indeterminacy

Since $\mathbf{a}_i \cdot s_i(t) = (-\mathbf{a}_i) \cdot (-s_i(t))$, flipping the sign of both the scalp topography column and the activation waveform leaves the measured sensor data unchanged. The sign of an IC waveform is physically arbitrary until projected back to sensor space!

3. Scaling (Energy) Indeterminacy

Multiplying a source by a constant $c$ while dividing the mixing column by $c$ preserves the exact same scalp voltages: $(c^{-1}\mathbf{a}_i)(c\,s_i(t)) = \mathbf{a}_i s_i(t)$. By convention, source activations are normalized to unit variance ($\sigma^2 = 1$), and all physical microvolt units ($\mu\text{V}$) are absorbed into the scalp map $\mathbf{A}$.

3

The Two-Step Geometric Separation
Whitening (Sphering) & FastICA / Infomax

How do modern ICA algorithms actually find the unmixing matrix without checking every possible angle? They divide the problem into two clean, decoupled geometric stages:

Step 1: Centering & Whitening (PCA)

First, subtract the mean ($E[\mathbf{x}] = 0$). Then, apply a whitening matrix $\mathbf{V} = \mathbf{D}^{-1/2}\mathbf{E}^T$ computed via eigendecomposition of the covariance matrix. Geometric effect: The skewed elliptical data cloud is transformed into a perfect unit sphere where all directions have equal variance $\mathbf{I}$.

Step 2: Pure Rotation (FastICA)

Because the data is now spherical, any further transformation must be a rigid orthogonal rotation $\mathbf{R}$ (no more stretching or shearing). This reduces the search space from $N^2$ unknown parameters to just $N(N-1)/2$ rotational degrees of freedom!

Algorithm Core Optimization Principle Convergence Speed EEG Application & Recommendation
PICARD
Ablin et al., 2018		 Preconditioned Infomax using quasi-Newton L-BFGS optimization Super-Linear / Quadratic (10x–50x faster than runica) ★ Recommended Modern Standard in MNE-Python (`method='picard'`)
Infomax / Ext. Infomax
Bell & Sejnowski, 1995		 Maximizes information transfer (joint entropy) through sigmoidal transfer functions Linear (natural gradient ascent) Worldwide gold standard in EEGLAB (`runica`, Extended Infomax)
FastICA
Hyvärinen, 1999		 Fixed-point iteration maximizing Negentropy / Non-Gaussianity approximation Cubic (ultra-fast, few iterations) Classic default in MNE-Python (`method='fastica'`)
JADE
Cardoso & Souloumiac, 1993		 Joint Approximate Diagonalization of 4th-order cumulant tensors Exact algebraic batch solution Suitable for small montages ($\le 32$ channels); $O(N^4)$ scaling

📏 The Makeig / Onton $k \cdot N^2$ Data Length Rule of Thumb

A common experimental pitfall is fitting ICA on insufficient data. Because ICA estimates 4th-order statistical moments across $N$ channels, the required number of time samples scales with the square of the channel count:

$$\text{Minimum Data Points: } \quad N_{\text{samples}} \ge k \cdot N_{\text{channels}}^2 \quad (k \ge 25\text{–}30)$$
Channel Count ($N$) Min. Data Samples ($25 \cdot N^2$) Min. Recording Time (@ 250 Hz) Min. Recording Time (@ 500 Hz)
32 Channels $25 \times 32^2 = 25,600$ ~1.7 minutes ~51 seconds
64 Channels $25 \times 64^2 = 102,400$ ~6.8 minutes ~3.4 minutes
128 Channels $25 \times 128^2 = 409,600$ ~27.3 minutes ~13.6 minutes

Laboratory Takeaway: If your recording is shorter than the minimum sample requirement, reduce the input dimensionality using PCA pre-whitening (`n_components = 30`) prior to fitting ICA to prevent overfitting random noise.

Interactive Whitening & Rotation Studio

Watch how raw correlated EEG channels morph into independent neural components through 4 discrete linear transformations. Step through manually or click Auto-Play Pipeline to watch the continuous geometric morphing and time-series decorrelation.

4-Stage BSS Studio

2D State-Space Geometry: Skewed Ellipse

Data Cloud Boundary / Axes

Channel Time-Series (x₁, x₂)

Channel 1 Channel 2
Transformation: x(t) = A s(t) + μ
Covariance: E[x₁x₂] ≠ 0 (Correlated) · Var ≠ 1

Raw recorded EEG channels have non-zero means (baseline offset), unequal variances, and high cross-channel correlation due to volume conduction, tracing a tilted, off-center ellipse.

4

End-to-End Pipeline & Code Reference
Python (MNE-Python) & MATLAB (EEGLAB)

Below are the standard production workflows for performing Blind Source Separation, component inspection, and artifact subtraction across modern EEG pipelines: 

import mne
from mne.preprocessing import ICA

# 1. Load continuous raw EEG dataset and set reference
raw = mne.io.read_raw_edf('sub-01_task-rest_eeg.edf', preload=True)
raw.set_eeg_reference('average')

# 2. Filter a temporary copy for ICA (1.0 Hz highpass removes non-stationary baseline drift)
raw_filt = raw.copy().filter(l_freq=1.0, h_freq=100.0)

# 3. Compute true data rank (Accounts for average referencing: N - 1)
rank = mne.compute_rank(raw_filt, rank='info')['eeg']

# 4. Fit ICA using PICARD (exact Infomax loss with quasi-Newton L-BFGS speed)
ica = ICA(n_components=rank - 1, method='picard', max_iter='auto', random_state=42)
ica.fit(raw_filt, picks='eeg')

# 5. Automatically detect blink (EOG) and heartbeat (ECG) artifact components
eog_indices, eog_scores = ica.find_bads_eog(raw_filt, ch_name='Fp1', threshold=3.0)
ecg_indices, ecg_scores = ica.find_bads_ecg(raw_filt, method='correlation', threshold=0.25)
ica.exclude = list(set(eog_indices + ecg_indices))
print(f"Marked {len(ica.exclude)} artifact components for rejection: {ica.exclude}")

# 6. Visualize components before rejection
ica.plot_components(picks='eeg')           # 2D interpolated scalp maps
ica.plot_sources(raw_filt)                 # Activation time-courses s(t)
ica.plot_overlay(raw, exclude=ica.exclude) # Sensor-level subtraction verification

# 7. Apply unmixing & reconstruction to the ORIGINAL unfiltered raw dataset
raw_clean = ica.apply(raw.copy(), exclude=ica.exclude)
% 1. Load dataset into EEGLAB
EEG = pop_loadset('filename', 'sub-01_task-rest.set');

% 2. High-pass filter a temporary copy at 1.0 Hz to stabilize Infomax gradient
EEG_filt = pop_eegfiltnew(EEG, 'locutoff', 1.0);

% 3. Compute true rank (Accounts for average referencing: N - 1)
dataRank = sum(eig(cov(double(EEG_filt.data'))) > 1e-7);

% 4. Run Extended Infomax ICA (extracts both sub- and super-Gaussian sources)
EEG_filt = pop_runica(EEG_filt, 'extended', 1, 'pca', dataRank - 1);

% 5. Transfer unmixing weights back onto the original unfiltered dataset
EEG.icaweights = EEG_filt.icaweights;
EEG.icasphere   = EEG_filt.icasphere;
EEG.icawinv     = EEG_filt.icawinv;
EEG = eeg_checkset(EEG);

% 6. Plot component topographies and activations
pop_topoplot(EEG, 0, [1:dataRank-1], 'Independent Components', [4 8]);
pop_eegplot(EEG, 0, 1, 1);

% 7. Subtract marked artifact components (e.g. ICs 1 and 3 for ocular blinks)
EEG_clean = pop_subcomp(EEG, [1 3], 0);

The Golden Pipeline Sequence: Mark bad channels $\rightarrow$ Filter a temporary copy at $1.0$ Hz $\rightarrow$ Compute matrix rank $\rightarrow$ Fit ICA with $\text{rank}-1$ components $\rightarrow$ Reject artifacts $\rightarrow$ Apply unmixing back onto original raw data $\rightarrow$ Interpolate bad channels last!

📚

Further Reading & Resources

  1. Independent Component Analysis: A Tutorial Introduction (Stone, 2004) 📖
    The definitive intuition-first book on ICA covering the cocktail party problem, state-space mixing geometry, Central Limit Theorem, and FastICA algorithms from geometric first principles.
  2. An Information-Maximization Approach to Blind Separation and Blind Deconvolution (Bell & Sejnowski, 1995) 📑
    Landmark Neural Computation paper introducing the Infomax algorithm that established Blind Source Separation as the worldwide gold standard for biomedical EEG processing.
  3. Independent Component Analysis of Electroencephalographic Data (Makeig, Bell, Jung, & Sejnowski, 1996) 🔬
    The historic NeurIPS publication that first applied Infomax ICA to electroencephalography, showing that scalp EEG signals separate into distinct cortical ERP sources and artifact streams.
  4. Faster Independent Component Analysis by Preconditioning with PICARD (Ablin, Cardoso, & Gramfort, 2018) ⚡
    IEEE Transactions on Signal Processing paper detailing the state-of-the-art PICARD algorithm that accelerates Infomax and Extended Infomax by 10x–50x using quasi-Newton L-BFGS preconditioning.
  5. EEGLAB Workshop: Biophysics and Geometry of Independent Component Analysis 🎥
    In-depth lecture by Scott Makeig (UCSD SCCN) on why cortical local field potentials form spatially stationary, statistically independent dipolar generators in human EEG.
  6. MNE-Python Guide: Repairing Artifacts with Independent Component Analysis (ICA) 🐍
    Official MNE tutorial covering highpass filtering prerequisites, fitting FastICA/PICARD, inspecting spatial topographies and time-courses, and artifact component subtraction.