PCA on an image time series: from stack to parametric maps

You captured a hundred frames of the same scene, and the answer is not in any single frame — it is in how each pixel's value moves over time. That is the situation behind every cooling curve filmed through a fixed camera, every curing adhesive, every slow drift or batch process you want to characterize rather than just detect. This guide shows you how to turn the stack into a pixels-by-time matrix, eigendecompose its covariance, and read out ranked temporal behaviors plus one parametric map per pixel for each. It also gives the memory arithmetic that decides whether your recording rig and computer can carry a series, and the checks that keep the maps honest.

The decision: decompose time, not space

Read the stack as a table. Every pixel is a variable, every frame is one observation of all of them, and the behaviors that dominate the recording are the eigenvectors of that table's covariance. The first eigenvector is the strongest common behavior over time, the second the strongest behavior left after removing the first, and so on. An eigenvalue's share of their sum is the variance that mode carries.

One choice follows before any code: pixels always outnumber frames in machine vision — millions against tens or hundreds — so form the frame covariance of size T × T and never the pixel covariance of size N × N. Once the stack is centered, the two share every nonzero eigenvalue, and the memory ratio between them is (N/T)². Everything else in this guide is bookkeeping around that choice.

Reshape the stack: one row per pixel

Store the series as a matrix X with N = H × W rows and T columns: column t is frame t flattened, and row n is one pixel's value across all frames — its temporal signature. A volumetric series in (x, y, z, t) reshapes the same way to x·y·z rows by t columns: every spatial position, whatever the volume, becomes one row. The rearrangement moves no data and resamples nothing; it is bookkeeping.

Plot a few rows before any math: single pixel traces, or one image row at one instant, as in the guide to plotting an image row intensity profile. Glitches, saturated stretches and dead pixels show up there immediately, and any of them poisons a decomposition.

Fix the flattening order now — row by row or column by column, both work — and write it down. You need the same order to reshape scores back into maps, and mixing the two directions transposes every map silently.

Center first, then eigendecompose the frame covariance

Subtract each pixel's mean over time before anything else. Skip the centering and the first component becomes the mean image — a brightness map, not a behavior. Centering also caps the rank at T − 1, which is why 100 frames carry at most 99 nonzero modes.

Now build C, the T × T covariance of the centered frames, and eigendecompose it. Each eigenvalue λk over the sum of all eigenvalues is mode k's share of the total variance, and the cumulative sum of those shares tells you where to stop keeping modes. In image processing this decomposition is the Karhunen–Loève transform, also called the Hotelling transform; it is the same construction that gives the optimal linear approximation of a random vector, as the Karhunen–Loève theorem overview records. A factor-analysis routine runs this same machinery plus one extra idea — an explicit independent noise term per pixel, fitted by iteration — which is why its loadings drift slightly from plain PCA scores on noisy stacks.

The shortcut is exact, not approximate. With more pixels than frames, the centered stack has rank at most T − 1, so the frame covariance and the pixel covariance carry identical nonzero eigenvalues, and each eigenvector of one maps to an eigenvector of the other. Decompose the small matrix and project back; you lose nothing.

Worked example: a 256 × 256 stack of 100 frames

Inputs: 256 × 256 pixels per frame, T = 100 frames covering one full process cycle, stored as float32 (4 bytes per value). The stack sits in memory as 65,536 rows by 100 columns.

Rule: center each row, form the frame covariance C = (1/99)·XcᵀXc, eigendecompose the 100 × 100 result, and score every pixel against the eigenvectors.

• Stack size: 65,536 × 100 × 4 B = 26,214,400 B ≈ 26 MB, or 52 MB in float64.
• Forming C: 65,536 × 100 × 100 ≈ 6.6 × 10⁸ multiply-adds ≈ 1.3 × 10⁹ FLOP, under a second on a modern multi-core CPU.
• C itself: 100 × 100 = 10,000 entries, 80 KB in float64; its eigendecomposition is effectively instant.
• Scores S = Xc·V: 65,536 × 100 again, another 26 MB and another fraction of a second.
• The pixel covariance you avoided: 65,536² ≈ 4.3 × 10⁹ entries ≈ 34 GB in float64, which is why this example needs the shortcut to run on a bench machine at all.

Result: the whole decomposition runs in seconds with a working set under about 110 MB, and every kept mode returns one 256 × 256 parametric map. Decision: keep modes until the cumulative variance share passes 90%. On a stack carrying one slow drift, one periodic behavior and noise, that cut lands at mode 2 or 3. If the cumulative share is still under 90% at mode 5, the variation is not low-rank — chase acquisition problems (illumination stability, trigger jitter, exposure) before trusting any map.

Score every pixel: from loadings to parametric maps

Scores give every pixel a coordinate on every mode. Column k of S, reshaped to your H × W in the flattening order you fixed, is the k-th parametric map: bright where a pixel follows mode k, dark where it opposes it, near zero where it ignores it.

Read the loading before the map. Plot the eigenvector's T values first: a near-constant positive shape is illumination drift, a smooth ramp is cooling or drying, a clean oscillation is a periodic mechanism. The map then says where in the scene that behavior lives. A loading that spikes at one frame index is a trigger glitch, not physics.

Fix eigenvector signs with one convention — for example, force each loading's mean to be positive — so maps stay comparable across runs. If a map looks like pure noise, first check the images carry enough structure at all: the spatial-regularity pattern-strength measure gives a quick number before you blame the math. Mapping flat indices back to row and column is the same index discipline as any geometric transform, as the homography pixel-coordinate walkthrough shows.

Three routes to the same modes

All three decompositions below return the same leading modes for a clean stack; they differ in cost and in what they claim about noise.

Snapshot PCA and pixel-covariance PCA return identical nonzero eigenvalues, so the second row exists only as the thing you avoid. Factor analysis spends its iterations moving variance into the noise term; start with PCA and move on only when a run needs that separation.

Recording and computing hardware for a camera series

Match frame rate to behavior first. At 30 fps, 100 frames span 3.3 s — right for cooling, curing, drying or a slow process drift; a 50 ms event needs a high-speed model with a hardware trigger. The decomposition itself is indifferent to frame rate; it only needs frames that cover the full behavior. The camera finder filters area scan cameras by resolution, interface and shutter.

The LS-CMR-0166 covers slow processes with 4000 × 3000 frames at 30 fps over USB3, while the MV-CH120-11UC brings a global shutter and the high-speed range for fast events.

Now size the computer. A 12 MP series — 4000 × 3000 pixels over 100 frames — is 1.2 × 10⁹ values: 4.8 GB as float32, 9.6 GB as float64, and the scores add the same again. Budget roughly 10 GB of working set per series in float32: an 8 GB machine swaps, 16 GB handles one series at a time, and 32 GB leaves room for float64 or post-processing in the same box.

For the controller side, the vision controllers category lists machines by memory, LAN ports and PoE. The MV-VC5525-256G70-NN comes with 32 GB of memory, PoE camera ports and eight built-in light channels; the MV-VC2140-128G40-1T comes with 8 GB of memory and four built-in light channels. Pair those figures with the budget above to pick per case.

Decision per case: a slow process at modest resolution runs on any modern controller with 8 GB of memory; a 12 MP, 100-frame series wants 16 GB as a floor and sits comfortably at 32 GB; a fast event swaps in the high-speed camera and changes nothing in the math.

When PCA on time is the wrong tool

• Rare, short events: variance ranking buries them in trailing modes. Use frame differencing or background subtraction instead.
• Statistically independent, non-Gaussian sources mixed inside single pixels: independent component analysis separates them where PCA's orthogonality cannot.
• Regions rather than per-pixel scores: segment or mask first, then decompose inside each region.
• Non-linear pixel dynamics: when the response itself is curved, a per-pixel physics fit beats linear modes.

One sanity check survives every case: if mode 1 has a flat loading and its map is your lamp, the decomposition worked and the illumination dominated it — stabilize the light and re-acquire. More applied diagnostics live on the applications topic page.

Checks before you trust the maps

1. Round-trip the math: reconstruct the centered stack from scores times transposed eigenvectors. The maximum absolute error must be 0 in float64, about 10⁻⁶ relative in float32; anything larger means the flatten and unflatten orders disagree.
2. Cut the modes where cumulative explained variance passes 90%, and check the increment per added mode is small past the cut.
3. Read the loading shapes: constant means illumination, ramp means drift, oscillation means periodic, one isolated spike means a glitch in a single frame.
4. Apply one sign convention to every run so maps compare across days.
5. Convert to a float type before centering; integer types silently clip the negative values that centering creates. The refresher on why 10-bit sensors ship 8-bit images covers what the camera actually hands you.
6. Overlay one map on a known edge or fiducial; a correct overlay proves the reshape and the geometry together.

To size a recording rig for your own series, start the system configurator at the camera step, or bring your stack parameters to the contact page and we will spec the camera-and-controller pair.