# NOTES — MRI 2-slice observation geometry vs low-rank DVF reconstruction

## What this runs

`run.sh` -> `experiment.py` (data load, axis-calibration gate, PCA/low-rank
fit + reconstruction across 9 leave-one-cycle-out folds x 3 seeds x 5 arms x
4 candidate slice spacings) -> `analyze.py` (pools per-seed metrics, runs
the falsifiable-prediction test on the 27 paired fold/seed samples, writes
`results.json` and `figs/*.png`).

Total wall-clock: ~6s (way under the 30-min budget — the whole computation is
small linear algebra: point clouds up to 10125 points, PCA via one SVD per
fold on a ~24657x40-90 matrix, K<=3, everything vectorized in numpy).

## Data sanity-check (axis calibration gate)

Before touching the ds5 data, `axis_calibration_gate()` reloads the ds15
subject-A cyc1-4 files and checks against the plan's anchor numbers:
- position span per axis (AP=120.1mm match; SI=150-162mm range brackets the
  ~156mm anchor; LR=75-84mm range brackets the ~77mm anchor)
- SI (axis1) is the axis with by far the largest span and largest max
  displacement (20.17mm vs ~10.5mm on AP/LR) — matches respiratory-motion
  physiology (liver moves most superior-inferior) and confirms the x,y,z ->
  AP,SI,LR axis mapping is right, not just "approximately numerically close."
- cyc1-own-t0 vs cyc3-ref-t0 offset: mean 3.81mm, max 9.46mm — matches the
  plan's anchor numbers almost exactly, confirming the reference-frame
  definition (`p_ref` = cyc3 t=0, not each cycle's own t=0) is implemented
  correctly. This was the most important check since a wrong reference would
  silently inflate all downstream displacement numbers.
- The gate hard-fails (assert) if any of this doesn't hold; it ran and
  passed on the first try after the reference-offset numbers matched to two
  decimal places.

## Design decisions not fully pinned down by the plan

- **Reference (mean) subtraction is per-subject fixed geometry.** Slice
  planes (sagittal LR-position, coronal AP-position) are computed once from
  each subject's cyc3-t0 reference point cloud and reused across all folds
  and seeds of that subject — this matches how a real treatment-planning
  scan would fix slice positions once, not re-locate them every breathing
  cycle.
- **27 paired samples** = 9 leave-one-cycle-out folds (A: 4 cycles, B: 5
  cycles) x 3 seeds — this reproduces the plan's explicit "27" number
  exactly and is the natural way to get 27 without inventing an
  unspecified extra sweep axis.
- **Best spacing selection (PASS 1)** picks the single d in {10,20,30,50}mm
  that minimizes the *pooled* RMSE (sum of squared error / total count)
  across all 9 folds x 3 seeds, separately for parallel-2sag and
  parallel-2cor. This avoids per-fold cherry-picking (the spacing choice is
  a fixed geometry decision, not adapted per test sample) — see
  `figs/spacing_sweep.png`. Both arms picked d=50mm, the largest spacing
  tested; RMSE is still monotonically decreasing at d=50, so 50mm is
  probably a boundary optimum of the swept range, not a true interior
  optimum — flagged in `results.json.deviations`, not chased further
  (no >50mm points meaningfully exist within this ROI-local point cloud).
- **`per_seed` in results.json** (required to have exactly 3 entries) is the
  RMSE *pooled over all 9 folds* for that seed, not a single fold. The
  fold-level detail (all 27 samples) is what the prediction test consumes
  internally and is visible in `figs/paired_diff_27.png`; the schema only
  asks for one number per seed so pooling across folds was the natural
  collapse.
- **ROI sites**: `liver_dome`/`liver_inferior`/`liver_lateral` are defined
  geometrically (nearest actual grid point to the 90th/10th SI percentile,
  or the lateral LR extreme, each near central position on the other two
  axes) since the ETH dataset has no anatomical labels — this is explicit
  in `out_of_scope`. All three sites are computed and stored in
  `arms[].metrics.roi_centroid_err_mm__<site>`; `liver_dome` isn't singled
  out anywhere except informally in comments — no site is part of the
  headline claim.

## What was NOT run (see `results.json.deviations` for the authoritative list)

- No exhaustive grid over K / sigma_sag / alpha / interleave / roi_site.
  Only the plan's stated default (K=3, sigma_sag=0.5, alpha=1,
  interleave=false) was used for the headline/prediction test. This is a
  scope cut, not a compute-budget cut — the compute itself was so cheap
  (whole pipeline in ~6s) that we could have easily afforded the sweep; the
  cut was made to keep the falsifiable-prediction test focused on exactly
  the settings the plan specified as decisive, and to leave enough of the
  30-minute human/agent review budget for double-checking the result below
  rather than generating more numbers.
- `interleave=true` (slices acquired in alternating frames, requiring a
  temporal-interpolation model) was not implemented at all.

## Result: the prediction was refuted, and not for the reason we guessed

Headline (dvf_rmse3d_mm, pooled, mean of 3 per-seed values):
parallel-2sag-best(d=50mm) = 0.719mm, orthogonal-sag-cor = 0.736mm — i.e.
orthogonal is *not* better overall (11/27 paired samples even favor it,
below the 19 threshold, and the mean difference is tiny: -0.009mm / -1.1%,
bootstrap 95% CI [-0.017, -0.003]mm). That much matches what the prediction
expected.

But the prediction's mechanistic story — "parallel-2sag can't see LR motion
directly, so its LR-axis error should be >=1.5x worse than orthogonal's" —
is wrong: measured ratio is 0.975 (parallel-2sag's LR RMSE, 0.565mm, is
actually *slightly lower* than orthogonal's, 0.579mm), which trips
refutation condition (c) on its own regardless of the other two conditions.

**Why, once checked against the oracle arm:** `oracle-full-dvf` (fits the
same K=3 PCA coefficients directly against the *full* true 3D field, no
slice restriction, no noise at all) still has `rmse_lr_mm = 0.529mm` —
almost as large as either slice-based arm. See `figs/axis_breakdown.png`:
LR error (~0.53-0.58mm) dwarfs AP/SI error (~0.31-0.33mm) *for every arm,
including the oracle*. That means LR reconstruction error here is dominated
by K=3 PCA-truncation capacity, not by which slice geometry observes it —
respiratory motion is coherent enough across axes that a couple of AP/SI
in-plane observations already pin down the same low-rank coefficient that
also determines LR motion, so adding an explicit (noisy, single-slice)
direct LR reading via the coronal plane doesn't help, and here even
slightly hurts (fewer total AP/SI point observations than 2 stacked
sagittal slices: 2093 vs 2071 observations — close, not the explanation by
itself, but consistent with 2 spatially-separated sagittal slices giving
more diverse spatial information about the shared mode shape than 1
sagittal + 1 coronal).

This was checked for a bug before accepting it (n_observations counts are
close between arms — 2071 vs 2093, not wildly imbalanced; the oracle-arm
comparison independently explains the effect via a mechanism unrelated to
slice geometry) and is reported as-is per the no-cherry-picking rule: no
parameters, seeds, or subsets were changed after seeing this result.

## Practical implication (bounded by scope)

Within this specific low-rank (K=3 PCA) motion-model setup, on 2 liver
subjects / 9 breathing cycles, geometry choice between parallel-2-sagittal
and orthogonal sagittal+coronal barely matters for overall 3D DVF accuracy,
and the "orthogonal sees the missing axis so it must reconstruct it better"
intuition does not hold once the model is low-rank — the model's own
inter-axis correlation structure already recovers most of the missing-axis
information from the observed axes. No claim is made beyond this scope
(2 subjects, K=3, linear PCA model, no image-domain simulation — see
`results.json.deviations` and the original plan's `out_of_scope`).