Skip to content
Source Separation Lab

Step 01Question 1 · synthetic data

Building a brain scan with known answers

Real fMRI data never tells you which signals went in. So the assignment starts by inventing them: six time courses, six spatial maps, a little noise, and a mixing rule. Every later step is graded against these ground truths.

The model is the linear one at the heart of the course. Each of the V = 441 pixels records a time series of N = 240 samples, and the whole recording is a product of temporal sources D and spatial sources A plus an error term:

X=DA+EX∈R240×441,  D∈R240×6,  A∈R6×441\begin{gathered} X = DA + E \\ X \in \mathbb{R}^{240 \times 441},\; D \in \mathbb{R}^{240 \times 6},\; A \in \mathbb{R}^{6 \times 441} \end{gathered}

To generate data that follows it, noise is added to both kinds of source before they are multiplied together:

X=(TC+Γt)(SM+Γs)Γt∼N(0,0.25),Γs∼N(0,0.015)\begin{gathered} X = (\mathrm{TC} + \Gamma_t)(\mathrm{SM} + \Gamma_s) \\ \Gamma_t \sim \mathcal{N}(0, 0.25),\quad \Gamma_s \sim \mathcal{N}(0, 0.015) \end{gathered}

Change anything below and every figure on this page, and on the next three, updates. The defaults are the values from the assignment, and the noise is regenerated from the same seed, so out of the box you are looking at the exact matrices from my 2021 notebook.

Using the 2021 parameters. Results match the report.

Time course timing
Onset0
Increment30
Duration15
Noise
Temporal variance σ²t0.25
Spatial variance σ²s0.015

Noise comes from a port of numpy's seeded generator, so seed 30034 reproduces the exact Γt and Γs used in 2021.

Figure 1.1

Six temporal sources, standardised

TC1AV 0 · IV 30 · dur 15

Standardised time course 1 over 240 time points−202060120180240

TC2AV 20 · IV 45 · dur 20

Standardised time course 2 over 240 time points−202060120180240

TC3AV 0 · IV 60 · dur 25

Standardised time course 3 over 240 time points−202060120180240

TC4AV 0 · IV 40 · dur 15

Standardised time course 4 over 240 time points−202060120180240

TC5AV 0 · IV 40 · dur 20

Standardised time course 5 over 240 time points−202060120180240

TC6AV 0 · IV 40 · dur 25

Standardised time course 6 over 240 time points−202060120180240
Each time course is a train of boxcar pulses defined by its onset, increment and duration, then centred and scaled to unit variance (sample standard deviation, as in the original). Hover or focus a panel and use the arrow keys to read values.

Why standardise rather than normalise? In 2021 I argued that the time courses are 0/1 indicators, so dividing by the ℓ2 norm only rescales them and leaves their offset in place. Standardising centres each one (no intercept is needed) and gives it unit variance, so no source dominates the regression just because it is switched on more often.

Figure 1.2

How related are the sources to each other?

Time courses (TC)

TC1TC1TC2TC2TC3TC3TC4TC4TC5TC5TC6TC6TC1 vs TC1: r = 1.0001TC1 vs TC2: r = −0.0000.00TC1 vs TC3: r = 0.1690.17TC1 vs TC4: r = 0.0860.09TC1 vs TC5: r = −0.0000.00TC1 vs TC6: r = 0.0860.09TC2 vs TC1: r = −0.0000.00TC2 vs TC2: r = 1.0001TC2 vs TC3: r = −0.029−0.03TC2 vs TC4: r = 0.1310.13TC2 vs TC5: r = 0.0000.00TC2 vs TC6: r = −0.131−0.13TC3 vs TC1: r = 0.1690.17TC3 vs TC2: r = −0.029−0.03TC3 vs TC3: r = 1.0001TC3 vs TC4: r = 0.0440.04TC3 vs TC5: r = −0.0000.00TC3 vs TC6: r = 0.1310.13TC4 vs TC1: r = 0.0860.09TC4 vs TC2: r = 0.1310.13TC4 vs TC3: r = 0.0440.04TC4 vs TC4: r = 1.0001TC4 vs TC5: r = 0.7750.77TC4 vs TC6: r = 0.6000.60TC5 vs TC1: r = −0.0000.00TC5 vs TC2: r = 0.0000.00TC5 vs TC3: r = −0.0000.00TC5 vs TC4: r = 0.7750.77TC5 vs TC5: r = 1.0001TC5 vs TC6: r = 0.7750.77TC6 vs TC1: r = 0.0860.09TC6 vs TC2: r = −0.131−0.13TC6 vs TC3: r = 0.1310.13TC6 vs TC4: r = 0.6000.60TC6 vs TC5: r = 0.7750.77TC6 vs TC6: r = 1.0001

Spatial maps (SM)

SM1SM1SM2SM2SM3SM3SM4SM4SM5SM5SM6SM6SM1 vs SM1: r = 1.0001SM1 vs SM2: r = −0.060−0.06SM1 vs SM3: r = −0.066−0.07SM1 vs SM4: r = −0.066−0.07SM1 vs SM5: r = −0.060−0.06SM1 vs SM6: r = −0.060−0.06SM2 vs SM1: r = −0.060−0.06SM2 vs SM2: r = 1.0001SM2 vs SM3: r = −0.066−0.07SM2 vs SM4: r = −0.066−0.07SM2 vs SM5: r = −0.060−0.06SM2 vs SM6: r = −0.060−0.06SM3 vs SM1: r = −0.066−0.07SM3 vs SM2: r = −0.066−0.07SM3 vs SM3: r = 1.0001SM3 vs SM4: r = −0.073−0.07SM3 vs SM5: r = −0.066−0.07SM3 vs SM6: r = −0.066−0.07SM4 vs SM1: r = −0.066−0.07SM4 vs SM2: r = −0.066−0.07SM4 vs SM3: r = −0.073−0.07SM4 vs SM4: r = 1.0001SM4 vs SM5: r = −0.066−0.07SM4 vs SM6: r = −0.066−0.07SM5 vs SM1: r = −0.060−0.06SM5 vs SM2: r = −0.060−0.06SM5 vs SM3: r = −0.066−0.07SM5 vs SM4: r = −0.066−0.07SM5 vs SM5: r = 1.0001SM5 vs SM6: r = −0.060−0.06SM6 vs SM1: r = −0.060−0.06SM6 vs SM2: r = −0.060−0.06SM6 vs SM3: r = −0.066−0.07SM6 vs SM4: r = −0.066−0.07SM6 vs SM5: r = −0.060−0.06SM6 vs SM6: r = 1.0001
Time courses: correlations between the six TCs. The strongest pair is currently TC4 and TC5 (r = 0.775); with the 2021 timings TC4, TC5 and TC6 overlap heavily (0.78, 0.78 and 0.60), which is hard to see in the traces but obvious here. Spatial maps: the six patches never overlap, so their correlations are small and slightly negative (at most |r| ≈ 0.073).

Figure 1.3

Six spatial sources on a 21 × 21 grid

SM1

SM2

SM3

SM4

SM5

SM6

Each map is a block of ones on a background of zeros. Before mixing, each map is flattened column by column into a 441-long row of SM (the original's sm.transpose().flatten()).

Why don't the maps need standardising? They are binary masks of equal height. Rescaling them would not change which pixels belong to which source, and the time courses already carry the scale of the signal.

Figure 1.4

Gaussian noise, temporal and spatial

Temporal noise Γt ~ N(0, 0.25)

Histogram of the temporal noise with the target normal density00.20.40.60.8−2−1012ValueDensity

Spatial noise Γs ~ N(0, 0.015)

Histogram of the spatial noise with the target normal density0123−0.4−0.200.20.4ValueDensity

Correlation of Γt columns

TC1TC1TC2TC2TC3TC3TC4TC4TC5TC5TC6TC6TC1 vs TC1: r = 1.0001TC1 vs TC2: r = 0.0040.00TC1 vs TC3: r = −0.010−0.01TC1 vs TC4: r = 0.0460.05TC1 vs TC5: r = 0.0280.03TC1 vs TC6: r = −0.012−0.01TC2 vs TC1: r = 0.0040.00TC2 vs TC2: r = 1.0001TC2 vs TC3: r = −0.083−0.08TC2 vs TC4: r = 0.0080.01TC2 vs TC5: r = 0.0050.01TC2 vs TC6: r = 0.0210.02TC3 vs TC1: r = −0.010−0.01TC3 vs TC2: r = −0.083−0.08TC3 vs TC3: r = 1.0001TC3 vs TC4: r = 0.0750.08TC3 vs TC5: r = −0.045−0.04TC3 vs TC6: r = 0.0120.01TC4 vs TC1: r = 0.0460.05TC4 vs TC2: r = 0.0080.01TC4 vs TC3: r = 0.0750.08TC4 vs TC4: r = 1.0001TC4 vs TC5: r = 0.0910.09TC4 vs TC6: r = −0.109−0.11TC5 vs TC1: r = 0.0280.03TC5 vs TC2: r = 0.0050.01TC5 vs TC3: r = −0.045−0.04TC5 vs TC4: r = 0.0910.09TC5 vs TC5: r = 1.0001TC5 vs TC6: r = 0.0320.03TC6 vs TC1: r = −0.012−0.01TC6 vs TC2: r = 0.0210.02TC6 vs TC3: r = 0.0120.01TC6 vs TC4: r = −0.109−0.11TC6 vs TC5: r = 0.0320.03TC6 vs TC6: r = 1.0001

Correlation of Γs rows

SM1SM1SM2SM2SM3SM3SM4SM4SM5SM5SM6SM6SM1 vs SM1: r = 1.0001SM1 vs SM2: r = 0.0020.00SM1 vs SM3: r = −0.047−0.05SM1 vs SM4: r = 0.0060.01SM1 vs SM5: r = 0.0290.03SM1 vs SM6: r = −0.042−0.04SM2 vs SM1: r = 0.0020.00SM2 vs SM2: r = 1.0001SM2 vs SM3: r = 0.0940.09SM2 vs SM4: r = 0.0510.05SM2 vs SM5: r = −0.025−0.03SM2 vs SM6: r = 0.0120.01SM3 vs SM1: r = −0.047−0.05SM3 vs SM2: r = 0.0940.09SM3 vs SM3: r = 1.0001SM3 vs SM4: r = 0.0820.08SM3 vs SM5: r = −0.030−0.03SM3 vs SM6: r = 0.0100.01SM4 vs SM1: r = 0.0060.01SM4 vs SM2: r = 0.0510.05SM4 vs SM3: r = 0.0820.08SM4 vs SM4: r = 1.0001SM4 vs SM5: r = 0.0230.02SM4 vs SM6: r = 0.0560.06SM5 vs SM1: r = 0.0290.03SM5 vs SM2: r = −0.025−0.03SM5 vs SM3: r = −0.030−0.03SM5 vs SM4: r = 0.0230.02SM5 vs SM5: r = 1.0001SM5 vs SM6: r = 0.0170.02SM6 vs SM1: r = −0.042−0.04SM6 vs SM2: r = 0.0120.01SM6 vs SM3: r = 0.0100.01SM6 vs SM4: r = 0.0560.06SM6 vs SM5: r = 0.0170.02SM6 vs SM6: r = 1.0001
Histograms of all 1440 entries of Γt and 2646 entries of Γs, with the target normal density drawn on top. Sample moments: Γt mean 0.0010, variance 0.2492 (target 0.250); Γs mean −0.0026, variance 0.0151 (target 0.015). The correlation matrices above show the noise is uncorrelated across sources: off-diagonal values sit near zero.

Figure 1.5

The observed dataset X = (TC + Γt)(SM + Γs)

−3.603.47

100 sampled time series from X

Variance of each of the 441 variables

Scatter plot of the variance of each column of X01230100200300400Variable (pixel)Variance
The heatmap shows X as 240 time points (rows) by 441 pixels (columns), before column standardisation. Pixels inside a source patch carry that source's time course and everything else is noise. Standardising gives every column unit variance, so the background noise is amplified to the same scale as the signal. The sampled time series are 100 columns of the raw X, the same columns the original sampled, with two highlighted. The variance plot splits the columns into two groups, background pixels with tiny variance and 161 pixels inside a patch with variance above 0.5. That spread is why X is standardised before fitting.

Next · Step 02

Retrieve

Recover the sources with least squares, ridge, lasso and PCR.