longitudinal/REGISTRATION_BENCHMARK_DEFORMABLE.md

# Registration Pipeline Benchmark — CT → T1c (M6) — Deformable

Companion to [REGISTRATION_BENCHMARK_RIGID.md](REGISTRATION_BENCHMARK_RIGID.md):
the **deformable** variant of the same six-pipeline benchmark. The rigid
benchmark asked "which toolkit aligns the head best with 6 DOF?"; this one
asks "which toolkit's native deformable workflow aligns the head (and the
tissue) best?" on the same 10 M6 CT→T1c test cases, the same unified
contract, and the same metrics.

Run date: 2026-10-02 (main run 05:43–07:52 local; serial `synthmorph`
re-run 09:58–10:07 local; accelerated `simpleitk` re-run 15:31–16:26
local). Follow-up ANTs-variant runs 2026-10-03 (08:44–09:12 local: five
stock antspyx presets — direct `SyNRA` / `SyN`, `Repro[sr|br]`,
`QuickRepro[sr|br]`; 11:00–11:26 local: the legacy two-call
`Rigid → SyNOnly` workflow; 12:44–12:49 local: re-run of the stock `SyN`
preset as `ants_syn`, replacing the initial `SyNOnly` row) on the same
10 cases (see the *ANTs variants* section). 10 M6 test cases, 3× NVIDIA A100 80GB PCIe, 96-core host, conda
env `longitudinal-313`. Main run: 6 workers (24 threads each, 2 workers
per GPU). `synthmorph` was re-run serially (1 worker, one container per
GPU at a time) because the concurrent run lost 6/9 cases to cuDNN
scratch-space allocation failures (see the SynthMorph note under
Interpretation); the re-run registered 9/9, matching the concurrent
run's scores to 4 decimals on the 3 cases it completed. `simpleitk`
was re-run after its spline grid was accelerated from 16 mm to 24 mm
(see the simpleitk row of the pipeline table; the original 16 mm results
are archived as `simpleitk_16mm_legacy.json`).

## Unified contract

Identical to the rigid benchmark: the returned transform `T` must satisfy
`sitk.Resample(ct, t1c, T, ...)` = CT on the T1c grid, i.e. **T maps T1c →
CT coordinates**. fixed = T1c, moving = CT, full volumes (no head box, no
mask, no crop, no pad), each toolkit's native automatic transform
initialization, CT clipped to [0, 80] HU only for the pipelines where the
rigid A/B showed clipping helps (niftyreg, synthmorph).

One deformable-specific change: every pipeline's output (elastix BSpline
combination, ANTs warp + affine composite, NiftyReg affine + F3D spline,
SimpleITK Euler + BSpline, SynthMorph joint field, IMPACT RecursiveBSpline)
is **rasterised to a SimpleITK `DisplacementFieldTransform` on the T1c
grid** before scoring, so all six pipelines return the same transform type
and the shared scoring / GTV warping / QA code is used verbatim.

## The six pipelines (deformable configurations)

Each pipeline runs its toolkit's **native deformable workflow** — a linear
init stage + a deformable stage, or a single deformable stage where that is
the toolkit's canonical setup:

| pipeline | deformable configuration | metric |
|---|---|---|
| **antspyx** | `type_of_transform='SyN'` (ANTs' native deformable): Affine stage (2100×1200×1200×0, shrink 4×2×2×1) + SyN[0.2,3,0] deformable stage (40×20×0, shrink 4×2×1), center-of-mass auto init | Mattes MI (32 bins), both stages |
| **elastix** | two-stage (the toolkit's canonical deformable workflow): EulerTransform rigid stage (2000 iters × 3 resolutions, CenterOfGravity auto init) → BSplineTransform stage (16 mm grid, spline order 3, 2000 iters × 3 resolutions) initialized from stage 1 via `-t0` | Mattes MI (32 bins), 5000 samples/level |
| **impact** | the ImpactReg `MR_CT_TS` **reference map in its native deformable form**: RecursiveBSplineTransform (14 mm final grid, order 3) driven by the MIND + TotalSegmentator static feature combination (4 per-level scales 400/300/200/200 iters, 2000 samples) with the reference's three-metric combination: IMPACT (w 1.0) + Mattes MI (w 0.4) + TransformBendingEnergyPenalty (w 10). Two M6 adaptations, same as the rigid benchmark: CenterOfGravity auto init (the reference's GeometricalCenter start leaves the heads ~100 mm apart and aborts) and RequiredRatioOfValidSamples 0.1 | IMPACT semantic (L1 + Dice) + Mattes MI + bending energy |
| **niftyreg** | reg_aladin rigid init (3-level pyramid, 15 TLS iters/level) → reg_f3d non-rigid (F3D B-spline Free-Form deformation, Modat et al. 2010; NMI metric, symmetric bending/linear regularisation, 3-level pyramid, 150 iters at the final level, 8 mm final grid) | NCC block matching (init) + NMI (F3D) |
| **simpleitk** | two-stage (the toolkit's canonical deformable demo): Euler3D stage (2000 iters/level, pyramid [4,2,1]) → BSplineTransform3D stage (24 mm grid, 50 iters/level on pyramid [4,2]), the BSpline's control points initialised to reproduce the Euler stage exactly. The 24 mm grid (5577 parameters) is the accelerated setting (the 2026-10-02 original used 16 mm / 13005 parameters, ~61 min/case under the harness): this build's spline per-iteration cost scales ~linearly with the control-point count (~2.8 s/iter at the shrink-4 level and ~20 s/iter at shrink-2 for 24 mm vs ~6.4 / ~45.8 s for 16 mm — independent of the metric sample count, which is fixed in this SimpleITK build; LBFGS2 was tested and fails to move the transform on this noisy MI), and the full re-run keeps the 16 mm quality within 0.003–0.029 MI per case (median 0.698 vs 0.705, ImpactLoss median unchanged) while retaining simpleitk's `m6_26TYQI3C` MI/ImpactLoss win — see the acceleration note under Interpretation. A full-resolution polish pass (10 iters at shrink-1) was tested and dropped the achieved MI (0.6891 at 32 mm), so the pyramid stays [4,2] | Mattes MI (32 bins), both stages |
| **synthmorph** | `joint` model (affine + deformable in one symmetric network pass — the toolkit's canonical deformable; the `deformable` model requires a pre-existing field and would bias the comparison) | anatomy-aware contrast-invariant network (pretrained) |

## Test cases

First 10 M6 labelled cases of `m6.discover()` (all with a GTV label):
`m6_22GD5VL2_20230724`, `m6_22IRRKU6_20240115`, `m6_23XIMO4M_20230222`,
`m6_24V62QQM_20230223`, `m6_24V62QQM_20231017`, `m6_24V62QQM_20240820`,
`m6_26TYQI3C_20221027`, `m6_2BBN5XA4_20250411`, `m6_2C4TKI6H_20240506`,
`m6_2CDNN7BH_20240102`.

(`m6_22GD5VL2_20230724` has an incomplete CT head — a known data issue,
not a pipeline failure; it is scored as ct_incomplete in every pipeline,
as in the rigid benchmark.)

(Ground-truth reference transforms are **not** part of this benchmark —
scoring is by achieved image similarity only, as in the rigid benchmark.)

## Metrics (identical for every pipeline)

- **Mattes MI** (nats) between the warped CT and the T1c over the MR head
  (largest connected component of the T1c mask) — achieved alignment, higher
  = better.
- **IMPACT loss** (vboussot/ImpactLoss; pretrained TotalSegmentator
  `TS/M730` deep features — weighted L1 over the two deepest layers, full
  volumes; the same scorer as the rigid benchmark) between the warped CT
  and the T1c — semantic-anatomy similarity, lower = better.
- **Registration wall time** (seconds, excludes I/O + scoring).

Per-case JSONs additionally retain the production skull-ring QA (gap +
concentricity, and the acceptance verdict) for reference.

### Result

Alignment is scored by the achieved **Mattes MI** (nats, warped-CT ~ T1c
over the MR head; higher = better) **and** the **ImpactLoss** semantic score
(warped-CT ~ T1c over the full volumes; lower = better), on the same
unclipped scoring volumes as the rigid benchmark, so all rows are comparable
on the same scales and against the rigid numbers. Full volumes, no mask /
crop / pad, fixed = T1c / moving = CT (shared rule). `niftyreg` and
`synthmorph` clip the CT (the moving image) to [0, 80] HU *before
registration* (per the rigid benchmark's clipping A/B); `antspyx`,
`elastix`, `impact` and `simpleitk` register the CT as-is. All six
transforms are deformable (each toolkit's native workflow) and are
rasterised to T1c-grid displacement fields before scoring.

All six pipelines registered 9/9 cases (the 10th is the known
incomplete-CT data issue, unscored in every pipeline). **Bold** = best
for the column (MI higher is better, ImpactLoss lower is better).

| pipeline | toolkit | registered | s/reg (med) | ImpactLoss (med/mean) | MI (med/mean, nats) |
|---|---|---|---|---|---|
| **antspyx** | ANTsPy 0.6.3 | 9 | 31.7 | 0.158 / 0.156 | **0.723 / 0.729** |
| **elastix** | ITK-Elastix 0.25.4 | 9 | 45.8 | 0.164 / 0.164 | 0.622 / 0.627 |
| **impact** | Docker `vboussot/elastix_impact` (Elastix 5.2.0 + IMPACT, deformable MIND+TS reference map) | 9 | 273.2 | 0.159 / 0.159 | 0.553 / 0.547 |
| **niftyreg** | NiftyReg 2.1.1 (reg_aladin + reg_f3d, Docker + CUDA) | 9 | 38.1 | 0.161 / 0.161 | 0.666 / 0.686 |
| **simpleitk** | SimpleITK 2.5.6 | 9 | 1369.9 | 0.158 / 0.159 | 0.698 / 0.704 |
| **synthmorph** | Docker `freesurfer/synthmorph` (TensorFlow, pretrained) | 9 | 49.1 | **0.156 / 0.156** | 0.627 / 0.640 |

Per-case MI (nats, unclipped score; **bold** = best for the case). The
`synthmorph` column is from the serial re-run (see the SynthMorph note
under Interpretation):

| case | antspyx | elastix | impact | niftyreg | simpleitk | synthmorph |
|---|---|---|---|---|---|---|
| m6_22GD5VL2_20230724 | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT |
| m6_22IRRKU6_20240115 | **0.679** | 0.580 | 0.548 | 0.637 | 0.654 | 0.620 |
| m6_23XIMO4M_20230222 | **0.723** | 0.622 | 0.573 | 0.711 | 0.693 | 0.618 |
| m6_24V62QQM_20230223 | **0.789** | 0.570 | 0.635 | 0.729 | 0.738 | 0.701 |
| m6_24V62QQM_20231017 | **0.720** | 0.622 | 0.580 | 0.665 | 0.692 | 0.628 |
| m6_24V62QQM_20240820 | **0.738** | 0.674 | 0.462 | 0.707 | 0.700 | 0.666 |
| m6_26TYQI3C_20221027 | 0.799 | 0.694 | 0.537 | 0.793 | **0.801** | 0.709 |
| m6_2BBN5XA4_20250411 | **0.666** | 0.585 | 0.436 | 0.612 | 0.641 | 0.569 |
| m6_2C4TKI6H_20240506 | **0.719** | 0.635 | 0.553 | 0.656 | 0.698 | 0.627 |
| m6_2CDNN7BH_20240102 | **0.726** | 0.660 | 0.602 | 0.666 | 0.718 | 0.624 |

Per-case ImpactLoss (semantic, unclipped score; **bold** = best (lowest)
for the case):

| case | antspyx | elastix | impact | niftyreg | simpleitk | synthmorph |
|---|---|---|---|---|---|---|
| m6_22GD5VL2_20230724 | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT |
| m6_22IRRKU6_20240115 | 0.1541 | 0.1595 | 0.1532 | 0.1530 | 0.1572 | **0.1512** |
| m6_23XIMO4M_20230222 | 0.1581 | 0.1657 | **0.1567** | 0.1605 | 0.1603 | 0.1573 |
| m6_24V62QQM_20230223 | 0.1473 | 0.1634 | 0.1497 | 0.1541 | 0.1502 | **0.1470** |
| m6_24V62QQM_20231017 | **0.1518** | 0.1591 | 0.1547 | 0.1620 | 0.1548 | **0.1518** |
| m6_24V62QQM_20240820 | **0.1545** | 0.1646 | 0.1694 | 0.1706 | 0.1584 | 0.1562 |
| m6_26TYQI3C_20221027 | 0.1580 | 0.1638 | 0.1656 | 0.1672 | **0.1562** | 0.1602 |
| m6_2BBN5XA4_20250411 | 0.1624 | 0.1663 | 0.1669 | 0.1648 | 0.1631 | **0.1616** |
| m6_2C4TKI6H_20240506 | 0.1594 | 0.1640 | 0.1585 | **0.1552** | 0.1626 | 0.1564 |
| m6_2CDNN7BH_20240102 | 0.1590 | 0.1662 | 0.1597 | 0.1605 | 0.1633 | **0.1583** |

### Interpretation

- **`antspyx` (Affine + SyN) leads on achieved MI** (median 0.723; best on
  8 of 9 cases and within 0.002 of best on the other —
  `m6_26TYQI3C`, 0.799 vs `simpleitk`'s 0.801) at the best median
  ImpactLoss of the gradient-based tools (0.158, within 0.002 of
  `synthmorph`'s 0.156), and it is the fastest full-workflow pipeline
  (31.7 s/case) — the best time/quality trade of the six.
- **`simpleitk` is second on MI** (0.698; it keeps its per-case win on
  `m6_26TYQI3C` — 0.801, still the highest MI of the whole benchmark,
  though now only 0.001 ahead of `antspyx` — and cedes its
  `m6_2CDNN7BH` win to `antspyx`, 0.718 vs 0.726) at a cost of ~23
  min/case under the harness — ~43× the fastest pipeline after the
  spline-grid acceleration below (this build's BSpline per-iteration
  cost scales ~linearly with the control-point count, independent of the
  metric sample count — see the pipeline table).
- **`simpleitk` acceleration (2026-10-02):** the original setting (16 mm
  spline grid, 13005 parameters) ran 3658.9 s/case (median) under the
  harness. The spline stage's per-iteration cost in this build scales
  ~linearly with the control-point count (16 mm: ~6.4 s/iter at the
  shrink-4 level, ~45.8 s at shrink-2; 24 mm: ~2.8 / ~20 s; 32 mm:
  ~1.2 / ~9.6 s), and the metric sample count is fixed in this SimpleITK
  build (`SetMetricAsMattesMutualInformation` exposes only the histogram
  bins); LBFGS2 was tested and its line search rejects every step on
  this noisy MI (the transform does not move). Switching to the **24 mm
  grid (5577 parameters) cut the wall to 1369.9 s median (2.7×)** at a
  cost of −0.003…−0.029 MI on 8 of 9 cases (+0.017 on `m6_2C4TKI6H`;
  median 0.698 vs 0.705), with the median ImpactLoss unchanged (0.158)
  and `simpleitk`'s `m6_26TYQI3C` MI/ImpactLoss wins intact. A 32 mm
  grid (a further ~1.9×) was rejected: it drops below `antspyx`'s 0.7994
  on `m6_26TYQI3C` (0.7902, and 0.7963 even with 100 fine-level iters).
  A full-resolution polish pass (10 iters at shrink-1) was also tested
  and dropped the achieved MI (0.6891 at 32 mm). The original 16 mm
  results are archived as `simpleitk_16mm_legacy.json`.
- **`niftyreg`'s `reg_f3d` stage delivers the largest deformable gain of
  any pipeline** (0.292 → 0.666 median MI, +0.374): it recovers the two
  gross-offset cases where the rigid block-matching locked onto local
  minima (`m6_24V62QQM_20240820` 0.118 → 0.707, `m6_2BBN5XA4` 0.069 →
  0.612), lifts the toolkit to third on MI, and it is the most
  concentric of the six (median skull-ring gap 5.5 mm mid / 3.8 mm std,
  6/9 accepted).
- **The deformable stage does not uniformly pay off.** Against the rigid
  benchmark's per-case MI it improves `antspyx` (0.652 → 0.723 median),
  `simpleitk` (0.647 → 0.698, accelerated 24 mm setting) and `niftyreg`
  (0.292 → 0.666), but
  **degrades** `elastix` (0.646 → 0.622; mixed per case, worst regression
  `m6_24V62QQM_20230223` 0.698 → 0.570) and `impact` (0.652 → 0.553 —
  lower than its own rigid init on **all 9 cases**). The two degradations
  are the coarse-grid BSpline configurations (14–16 mm) with heavy
  bending regularisation — for `impact` the reference map's
  TransformBendingEnergyPenalty at weight 10 — so the deformable stage
  smooths rather than refines, and (for `impact`) its IMPACT-dominant
  objective is simply not the achieved-MI objective.
- **`impact`'s deformable reference map is the weakest on the intensity
  scale** (MI 0.553) even though its *semantic* score is competitive
  (ImpactLoss 0.159, within 0.003 of the best median) — a
  metric-objective/achieved-MI divergence the rigid benchmark did not
  show (there the two scores agreed, with `impact` in the top tier on
  both). At 273 s/case it is the worst time/quality trade of the six.
- **`synthmorph` (serial re-run, 9/9) splits the leaderboards: best on
  the semantic scale, mid-pack on the intensity scale.** Median
  ImpactLoss 0.156 — the best of the six — with the per-case best on 4 of
  9 (plus a 4-decimal tie with `antspyx` on `m6_24V62QQM_20231017`),
  while its MI (0.627) is fourth, statistically tied with `elastix`
  (0.622). The pretrained contrast-invariant network optimises learned
  anatomy rather than raw CT/MR intensity synchrony, so it wins where
  anatomy matters and loses where the intensity histogram matters. Fast
  (49.1 s/case) and the read-back cross-check against the tool's own
  resampling is NCC 1.0 on all 9 cases.
- **The two score families discriminate differently here than in the
  rigid benchmark.** The ImpactLoss scale is nearly saturated (medians
  0.156–0.164, ~5 % spread; per-case bests spread across five tools,
  `elastix` with 0/9), whereas the MI scale spans 0.553–0.723
  (0.17 nats absolute) and the per-case bests go 8/9 `antspyx`, 1/9
`simpleitk`, 0/9 the rest. The deformable stage's main effect is fine
   intensity alignment on top of the gross anatomy the rigid stages
   already place — exactly what MI measures and ImpactLoss is less
   sensitive to. For ranking the deformable workflows, achieved MI is
   the more discriminative score.
- **`synthmorph` infrastructure note.** In the main concurrent run (6
  workers, 2 containers per GPU) 6 of 9 cases failed with `CUDNN failed
  to allocate the scratch space ... [Op:Conv3D]` inside the joint model's
  256-channel 128³ hyperconvolution. Each case is a fresh
  `docker run --rm` container, so the failures are cuDNN scratch-space
  contention when two SynthMorph containers share one A100, not a
  pipeline fault: the serial re-run registered 9/9, reproduced the
  concurrent run's MI/ImpactLoss to 4 decimals on its 3 completed cases
  (deterministic model), and ran ~17 s/case faster without the
  contention (49.1 s vs 65.7 s median). The concurrent-run JSON is archived as
  `synthmorph_concurrent_run.json`.
- **Production skull-ring QA (reference only, kept in the per-case
  JSONs):** acceptance 6/9 `niftyreg`, 5/9 `elastix` / `synthmorph`,
  4/9 `antspyx` / `simpleitk`, 2/9 `impact`. Concentricity: `niftyreg`
  best (median gap 5.5 mm mid / 3.8 mm std), `simpleitk` (19.1 mm std,
  accelerated setting) and `impact` (18.9 mm) worst on the rotational
  term. The same weak correlation with perceived alignment noted in the
  rigid benchmark, so it is reported for reference only.

## ANTs variants: direct SyNRA / SyN and the Repro / QuickRepro workflows

The `ants_syn` row below is antspyx's stock **`SyN` preset** (Affine
2100×1200×1200×0 + SyN[0.2,3,0] 40×20×0, Mattes MI 32, center-of-mass
init) — the same configuration as the `antspyx` row of the main table
above, re-run in this follow-up batch and reproduced case-for-case within
0.002 nats (results in `ants_syn.json`), listed under its own name to
avoid confusing it with the `ants_synra` / `ants_rigid_synonly` variants.
As a follow-up (2026-10-03, 6 workers as the main run) six further ANTs
configurations were benchmarked on the same 10 cases and the same
contract — the question: is there a better stock ANTs deformable workflow
for CT→T1c than the `SyN` preset?

| pipeline | antspyx preset | stages (all COM init unless noted) |
|---|---|---|
| **ants_rigid_synonly** | legacy two-call (`Rigid` + `SyNOnly`) | Rigid[0.25] (Mattes 32, 500×500×200, 4×2×1, 3×2×0) → SyNOnly SyN[0.2,3,0] (Mattes 32, 100×100×0, 4×2×1, 2×1×0) initialized from the Rigid stage; `Composite.h5` read-back |
| **ants_synra** | `SyNRA` | Rigid[0.25] + Affine[0.25] (each Mattes 32, 2100×1200×1200×0, 4×2×2×1, 3×2×1×0) + SyN[0.2,3,0] (Mattes 32, 40×20×0, 4×2×1) |
| **ants_synreprosr** | `antsRegistrationSyNRepro[sr]` | Rigid[0.1] (GC 1/1/Regular/0.25, 1000×500×250×100, 8×4×2×1, 3×2×1×0vox) + SyN[0.2,3,0] (CC 1/2, 50×20, 2×1, 1×0vox) |
| **ants_synreprobr** | `antsRegistrationSyNRepro[br]` | the `Repro` Rigid[0.1] stage + BSplineSyN[0.2,26,0,3] (CC 1/2, 50×20, 2×1, 1×0vox) |
| **ants_synquickreprosr** | `antsRegistrationSyNQuickRepro[sr]` | Rigid[0.1] (GC, 1000×500×250×0, 8×4×2×1) + SyN[0.2,3,0] (CC 1/2, 50×0, 2×1) |
| **ants_synquickreprobr** | `antsRegistrationSyNQuickRepro[br]` | Rigid[0.1] (GC, 1000×500×250×0, 8×4×2×1) + BSplineSyN[0.2,26,0,3] (CC 1/2, 50×0, 2×1) |

`sr`/`br` select the deformable stage (SyN vs BSplineSyN, antspyx preset
naming — not a brain mask); these are the same stage configurations the
`antsRegistrationSyNRepro*.sh` / `antsRegistrationSyNQuickRepro*.sh`
workflows use. `ants_rigid_synonly` is the two-call workflow that an early
version of this benchmark used: a Rigid stage, then a
*deeper* SyNOnly stage (100 iterations at every level incl. full
resolution, vs the `SyN` preset's 40×20×0) initialized from the Rigid
stage's transform. Full volumes, no CT clip, each toolkit preset's own
metrics (GC/CC/Mattes as above). As for the other rows, the stage files
are rasterised to the T1c-grid displacement field and scored with the
shared code; the composite read-back was cross-checked per case against
ANTs' own resampling (**NCC 1.0 on all 63 scored cases**,
`own_warp_ncc_head` in the JSONs). One read-back subtlety: when ANTs runs Rigid + Affine as separate
linear stages (`SyNRA`), it writes the merged linear stage as a *binary*
`0GenericAffine.mat` that SimpleITK cannot parse; the pipeline recovers it
by point-probing the transform (exact for affine maps, verified to
2×10⁻⁴ mm).

All seven rows registered 9/9 (the 10th is the known incomplete-CT data
issue). **Bold** = best for the column.

| pipeline | registered | s/reg (med) | ImpactLoss (med/mean) | MI (med/mean, nats) |
|---|---|---|---|---|
| **ants_rigid_synonly** | 9 | 208.4 | 0.156 / 0.155 | **0.753 / 0.763** |
| **ants_syn** | 9 | 33.8 | 0.158 / 0.156 | 0.722 / 0.729 |
| **ants_synra** | 9 | 43.6 | 0.158 / 0.156 | 0.722 / 0.729 |
| **ants_synreprosr** | 9 | 294.7 | **0.151 / 0.153** | 0.530 / 0.529 |
| **ants_synreprobr** | 9 | 370.3 | 0.156 / 0.157 | 0.541 / 0.537 |
| **ants_synquickreprosr** | 9 | 83.7 | 0.154 / 0.156 | 0.542 / 0.534 |
| **ants_synquickreprobr** | 9 | 96.0 | 0.158 / 0.158 | 0.549 / 0.538 |

Per-case MI (nats, unclipped score; **bold** = best for the case):

| case | ants_syn | ants_rigid_synonly | ants_synra | ants_synreprosr | ants_synreprobr | ants_synquickreprosr | ants_synquickreprobr |
|---|---|---|---|---|---|---|---|
| m6_22GD5VL2_20230724 | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT |
| m6_22IRRKU6_20240115 | 0.681 | **0.741** | 0.681 | 0.623 | 0.646 | 0.630 | 0.627 |
| m6_23XIMO4M_20230222 | 0.722 | **0.760** | 0.722 | 0.530 | 0.544 | 0.529 | 0.530 |
| m6_24V62QQM_20230223 | 0.787 | **0.824** | 0.789 | 0.520 | 0.528 | 0.542 | 0.549 |
| m6_24V62QQM_20231017 | **0.720** | 0.716 | **0.720** | 0.473 | 0.475 | 0.469 | 0.462 |
| m6_24V62QQM_20240820 | 0.738 | **0.753** | 0.738 | 0.524 | 0.525 | 0.517 | 0.524 |
| m6_26TYQI3C_20221027 | 0.799 | **0.857** | 0.798 | 0.557 | 0.541 | 0.599 | 0.632 |
| m6_2BBN5XA4_20250411 | 0.666 | **0.696** | 0.665 | 0.230 | 0.303 | 0.235 | 0.297 |
| m6_2C4TKI6H_20240506 | 0.719 | **0.751** | 0.719 | 0.672 | 0.649 | 0.654 | 0.621 |
| m6_2CDNN7BH_20240102 | 0.725 | **0.765** | 0.726 | 0.629 | 0.622 | 0.631 | 0.599 |

Per-case ImpactLoss (semantic, unclipped score; **bold** = best (lowest) for
the case):

| case | ants_syn | ants_rigid_synonly | ants_synra | ants_synreprosr | ants_synreprobr | ants_synquickreprosr | ants_synquickreprobr |
|---|---|---|---|---|---|---|---|
| m6_22GD5VL2_20230724 | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT | incomplete-CT |
| m6_22IRRKU6_20240115 | 0.1540 | 0.1530 | 0.1542 | **0.1442** | 0.1482 | 0.1474 | 0.1493 |
| m6_23XIMO4M_20230222 | 0.1581 | 0.1565 | 0.1581 | **0.1505** | 0.1549 | 0.1535 | 0.1585 |
| m6_24V62QQM_20230223 | 0.1472 | **0.1464** | 0.1472 | 0.1515 | 0.1553 | 0.1535 | 0.1560 |
| m6_24V62QQM_20231017 | 0.1519 | **0.1512** | 0.1516 | 0.1559 | 0.1575 | 0.1577 | 0.1581 |
| m6_24V62QQM_20240820 | 0.1546 | **0.1545** | 0.1547 | 0.1553 | 0.1578 | 0.1586 | 0.1586 |
| m6_26TYQI3C_20221027 | 0.1579 | 0.1562 | 0.1578 | **0.1561** | 0.1574 | 0.1602 | 0.1610 |
| m6_2BBN5XA4_20250411 | 0.1624 | **0.1611** | 0.1624 | 0.1680 | 0.1695 | 0.1690 | 0.1707 |
| m6_2C4TKI6H_20240506 | 0.1593 | 0.1582 | 0.1594 | **0.1490** | 0.1526 | 0.1517 | 0.1545 |
| m6_2CDNN7BH_20240102 | 0.1589 | 0.1582 | 0.1589 | **0.1493** | 0.1560 | 0.1538 | 0.1580 |

### Interpretation (ANTs variants)

- **`ants_rigid_synonly` is the best ANTs configuration of the seven** —
  the only one that beats the current `ants_syn` row: median MI **0.753**
  (per-case best on 8 of 9; its 0.857 on `m6_26TYQI3C` is the highest MI
  of the whole benchmark, ahead of the previous 0.801) at ImpactLoss
  0.156 (second-best median after `synthmorph`'s 0.156). The gain over
  the `SyN` preset comes from the *deeper* SyN stage: 100 iterations at
  every level **including full resolution** (100×100×0 vs 40×20×0)
  refining the Rigid stage's pose. The cost is 208.4 s/case (~6.2× the
  `ants_syn` time) — the full-resolution SyN iterations dominate.
- **`ants_synra` is indistinguishable from the `ants_syn` row** on
  quality (median MI 0.722 vs 0.722, ImpactLoss 0.158 vs 0.158; per-case
  deltas within 0.002 nats, identical on 6 of 9) at **+29 % wall time**
  (43.6 s vs 33.8 s): the explicit Rigid stage ahead of the Affine stage
  buys nothing — on these pairs the 12-DOF stage from the center-of-mass
  init already finds the gross pose, so the extra 6-DOF stage just
  re-solves it.
- **The Repro / QuickRepro family is semantically competitive but weaker on
  intensity:** median MI 0.530–0.549 (vs 0.722 for the `SyN` preset, i.e.
  `ants_syn`),
  driven by the GC-based Rigid stage and the short 2-level SyN/BSplineSyN
  refinement (50×20 or 50×0 iters, 2×1 pyramid, CC 2-bin metric) — a
  workflow tuned for intra-modality T1w-T1w. `ants_synreprosr` still takes
  the **best ImpactLoss median of all thirteen pipelines here** (0.151;
  per-case wins on 4 of 9), so the semantic/intensity split of the main
  table carries over to the stock ANTs workflows.
- **Time tracks budget, not quality:** the full Repro Rigid stage
  (1000×500×250×100, 8×4×2×1, GC) dominates — `ants_synreprobr` is the
  slowest ANTs configuration (370 s/case) at the same (low) MI as the
  84–96 s `QuickRepro` variants (which drop the 4th rigid level and the
  2nd SyN level, changing the MI by < 0.02). Only `ants_rigid_synonly`
  challenges `ants_syn`'s quality (better MI at ≈6× the time); every
  other variant is worse for more or for the same.
- **Production skull-ring QA (reference only):** acceptance 4/9
  `ants_syn` / `ants_synra` / `ants_rigid_synonly`, 2/9
  `ants_synreprobr`, 1/9 the rest; rotational dispersion worst in the
  `Repro` / `QuickRepro` family (24–29 mm std) vs 14–15 mm for
  `ants_syn` / `ants_synra` / `ants_rigid_synonly` — consistent with
  their alignment.

## Transform conventions (deformable-specific)

Unlike the rigid benchmark (where the elastix-style toolkits came out in the
opposite direction 8/9 cases and needed the best-of-T vs T⁻¹ orientation
gate), every deformable output here is consumed in the documented
fixed→moving = T1c→CT direction; the gate is replaced by an analytic
read-back plus a skull-ring sanity check. Two non-obvious conventions were
established empirically (validated against each tool's own resampling, NCC
≥ 0.998 on the head):

- **NiftyReg and surfa (FreeSurfer's NIfTI stack) decode RAS+ NIfTI qforms
  naively, i.e. they work in an x/y-mirrored world.** A displacement delta
  (dx, dy, dz) stored by such a tool is (−dx, −dy, dz) in the true RAS
  world; the field grid (file indices) is unaffected. The read-back flips
  the x/y components of the dense field. (ANTs writes correct RAS fields —
  no flip needed.)
- **NiftyReg reg_f3d writes only the control-point grid** (outputCPP.nii);
  `reg_transform -disp` evaluates it to a dense field on the reference grid.
  The reg_aladin affine is folded into the control points before
  optimization (reg_f3d initialises the spline with the affine's
  deformation field), so the dense field **is** the full transform — no
  separate composition step. (The `-comp` subcommand was tested and does
  not reproduce reg_f3d's internal transform; it is not used.)
- **SynthMorph joint**: driven as `register -m joint moving fixed`
  (moving = clipped CT first, fixed = T1c second); the `-t` output is the
  displacement field on the fixed (T1c) grid — the contract transform
  directly. The tool's own `-o` resampling is written as a cross-check
  (per-case `own_resample_ncc_head` in the JSONs).
- **ANTs SyN**: the two stage outputs (001Warp.nii.gz + 000GenericAffine.mat)
  are composed as warp-then-affine (verified to reproduce ANTs' own
  `warpedmovout` to < 0.001 MI) and rasterised.
- **Elastix / IMPACT**: the combination transform (auto-init folded in) is
  rasterised — elastix in-process via `itk.TransformToDisplacementFieldFilter`
  (exact, machine precision on affine fields), IMPACT by evaluating the
  BSpline control-point file directly (SimpleITK cannot read elastix's
  BSpline .txt format; scipy's order-3 B-spline = ITK's kernel).

## Reproducing

```bash
source /opt/conda/etc/profile.d/conda.sh && conda activate longitudinal-313
cd /mnt/b4/xfr/git26/longitudinal

# all six deformable pipelines (results go to results/bench_reg/bench_registration_deformable/)
# the 2026-10-02 main run used --workers 6 --threads 24 (~2 h 9 min wall,
# dominated by simpleitk's ~111 min phase)
python scripts/bench_registration.py --deformable --workers 6 --threads 24 --qa-images

# single pipeline / case subset
python scripts/bench_registration.py --deformable --pipelines simpleitk impact --n-cases 3

# the ANTs-variant follow-up (direct SyNRA / SyN + Repro / QuickRepro
# sr/br; results go to the same results/bench_reg/bench_registration_deformable/ dir,
# 2026-10-03 run, 6 workers, ~28 min wall; `ants_syn` is the re-run of the
# stock SyN preset that replaced the initial SyNOnly row)
python scripts/bench_registration.py --deformable --pipelines "ants_synra ants_syn ants_synreprosr ants_synreprobr ants_synquickreprosr ants_synquickreprobr" --workers 6 --threads 24 --qa-images

# the legacy two-call Rigid -> SyNOnly(100,100,100) workflow (2026-10-03,
# 6 workers, ~25 min wall; dominated by the full-resolution SyN iters)
python scripts/bench_registration.py --deformable --pipelines "ants_rigid_synonly" --workers 6 --threads 24 --qa-images

# serial SynthMorph re-run (the 2026-10-02 final numbers): one container per
# GPU at a time; two concurrent containers on one A100 hit cuDNN
# scratch-space OOM in the joint model's 256-ch 128^3 hyperconvolution
python scripts/bench_registration.py --deformable --pipelines synthmorph --workers 1 --threads 24 --qa-images

# combine the per-pipeline JSONs into summary.json / summary.md / plot
python scripts/bench_registration.py --deformable --merge
```

## Visual QA images

`--qa-images` writes a `checkreg`-style PNG per (pipeline, case) to
`results/bench_reg/bench_registration_deformable/qa/<pipeline>__<case>.png`
(9 per pipeline; the incomplete-CT case is skipped) — same format as the
rigid benchmark: a 3×2 grid of the mid-axial / mid-coronal / mid-sagittal
slices through the head center, **left column** checkerboard blend of the
warped CT vs the T1c (clean registration = seamless mosaic), **right
column** additive color overlay (CT red, T1c green; agreement = yellow,
misregistration = red/green fringing). CT is brain-windowed 0–80 HU, T1c
2–98 % head window; the title carries the pipeline, case key and the
achieved Mattes MI.