Wavelet
Analysis

The PSD (Day 9) tells you which frequencies are present but not when they occur. Neural oscillations wax and wane over time — alpha bursts during blinks, theta surges during memory encoding, beta rebounds after movement. Wavelet analysis gives you a time-frequency map: power at every frequency at every time point. This is arguably the most important analysis tool in cognitive electrophysiology.

The standard FFT computes frequency content over the entire signal — it has perfect frequency resolution but zero time resolution. If alpha power increases for just 200 ms during a trial, the FFT smears it across the whole epoch.

You could chop the signal into short windows and FFT each one (the Short-Time Fourier Transform), but this creates a fundamental tradeoff: short windows give good time resolution but poor frequency resolution; long windows give the opposite. This is the Heisenberg–Gabor uncertainty principle for signal processing: you cannot simultaneously have perfect resolution in both time and frequency.

Wavelets offer an elegant solution: adaptive resolution. Low frequencies use wide wavelets (good frequency resolution), and high frequencies use narrow wavelets (good time resolution).

The complex Morlet wavelet is a sine wave multiplied by a Gaussian envelope:

ψ(t) = exp(2πift) · exp(−t² / 2σ²)

The sine part defines the center frequency f. The Gaussian part defines the temporal width σ. The ratio between these — traditionally parameterized as the “number of cycles” n = 2πfσ — controls the time-frequency tradeoff.

More cycles = narrower frequency response, wider in time.
Fewer cycles = wider frequency response, tighter in time.

import mne
import numpy as np

# Define frequencies of interest
freqs = np.arange(4, 40, 1)   # 4 to 39 Hz in 1-Hz steps

# Number of cycles: trade off resolution
# Option 1: fixed cycles (e.g., 7)
n_cycles = 7

# Option 2: frequency-adaptive (Cohen's recommendation)
n_cycles = np.linspace(3, 10, len(freqs))  # fewer at low f, more at high f

# Compute time-frequency representation
tfr = epochs.compute_tfr(method='morlet',
                         freqs=freqs,
                         n_cycles=n_cycles,
                         return_itc=False)
tfr.plot([0], baseline=(-0.3, 0), mode='logratio')

Computing the wavelet transform is convolution: you slide the wavelet along the signal, multiplying point-by-point and summing. At each time point, the result tells you how well the signal matches the wavelet at that frequency.

Because the Morlet wavelet is complex-valued, the output at each time-frequency point is also complex. The squared magnitude gives you power; the angle gives you phase. This is how wavelet analysis simultaneously delivers time-frequency power maps and phase information.

Cohen (2019) argued that the traditional “number of cycles” parameterization is unintuitive. Instead, he advocates specifying the wavelet width using the Full Width at Half Maximum (FWHM) — the duration in milliseconds over which the wavelet’s Gaussian envelope exceeds half its peak value.

FWHM is directly interpretable: “this wavelet integrates activity over a ~300 ms window” is more intuitive than “this wavelet uses 5 cycles at 10 Hz.” The relationship: FWHM = n_cycles / f × 2√(2 ln 2).