Results
Every method per stage, every full-pipeline combination, and the patterns behind the numbers.
Combination matrix
Dipole inversion ↓ × background removal →, chained end-to-end · best methods top-left · click a cell to inspect.
Combined single-method algorithms (one step straight to a χ map)
No pipeline combinations for this field mapping yet.
Patterns the leaderboard tables don't surface, computed live from the same results. The recurring lesson: a method's standing depends heavily on what you score, and on the pipeline it sits in, not just on a single number.
The isolated ranking misleads about real pipelines
Dipole inversion ranked on the clean ground-truth field (left) versus chained after real background removal (right), by NRMSE, best at top. Green rises when composed (robust to upstream error); red falls. TV-family methods top the clean ranking, then collapse; several networks climb.
Accurate but fragile, or reliable but modest
Each dipole method by its mean NRMSE across background-removal partners (x) versus how much that error varies between partners (y). Bottom-left is the sweet spot: accurate and consistent. The networks cluster there, accurate and dependable across upstream partners; the iterative methods scatter, and while the best few rival the networks, many degrade sharply with a poor upstream stage.
Greedy pipeline-building is close, but not optimal
Picking the strongest method at each stage independently lands within NRMSE of the best full pipeline; composition is fairly forgiving. But the two share no methods: stages interact, so the actual best combo isn't the one you'd assemble greedily.
Top pipelines (of )
Which stage decides your result?
How much the choice at each stage reshuffles the full-pipeline ranking (rank-based, so one broken method can't skew it), across every pipeline. Background removal and dipole inversion matter about equally: each reshuffles the result far more than field mapping. The notable part: the background-removal choice is as influential as the dipole-inversion step the field usually fixates on.
methods compared at this stage
Rank-based on purpose: under a raw-NRMSE variance measure a single catastrophic method (AMP-PE fails on in-silico, ~4000% NRMSE) would make its stage look all-important. Ranking each pipeline first keeps the comparison robust.
Error compounds through the pipeline
The best achievable xSIM (blue) and the median (green box = interquartile range) at three pipeline depths: dipole inversion alone from the ground-truth local field, then with real background removal added (from the ground-truth total field), then a full pipeline from raw phase. Each real upstream stage roughly halves the median and pulls the ceiling down with it. Upstream error is never recovered downstream, so a method's isolated score badly overstates its end-to-end result.
Pick your dipole independently of background removal
For each dipole method, its rank among dipoles (dot = mean, bar = min–max) across the strong background-removal methods. Most bars are tight: the dipole ranking barely shifts with the upstream stage, so the stages are near-separable — the methods at the top lead wherever they sit, and those at the bottom trail wherever they sit. The few with wide bars swing; those are the genuinely upstream-sensitive methods. Coloured by family.
The best background-removal method depends on what you measure
Each background-removal method's rank (1 = best, green good → red poor) under five metrics, rows sorted by xSIM. The similarity metrics (xSIM, NRMSE, HFEN) are nearly interchangeable (their columns match). But calcification and iron linearity reshuffle the ranking: iron linearity is essentially uncorrelated with similarity (near-zero rank correlation), and the deep-learning bfrnet is worst-tier on similarity yet best on iron quantification. A single headline number hides which sources a method actually recovers.
The best method depends on what you measure
Isolated dipole methods across five metrics, each normalised so higher is better. Grey lines are all methods; a coloured line is the winner of one axis, showing it gives up ground elsewhere: the global-error (NRMSE / xSIM) leader slips down the calcification axis, and the fast direct methods under-quantify iron in deep grey matter. No single method wins every metric.
Where the error lives
Isolated dipole error by region and artefact (rows sorted by overall NRMSE). Green is good, red is poor, coloured within each column. Venous blood is universally hardest; the calcification columns reshuffle the ranking entirely.
Which metrics are redundant
Rank agreement (Spearman ρ) between metrics over the isolated dipole methods, sign-aligned so +1 means they reward the same methods. NRMSE, detrended NRMSE and correlation are near-interchangeable; DGM linearity, calcification and runtime are their own axes, so a single headline number can't stand in for them.
One region, every method
Reconstructed susceptibility inside a single segmented region, per isolated dipole method: dot = mean, bar = ±1 std, tick = median (ppm). The green band is ground truth (mean ± std; dashed line = median). Methods are sorted by how far their regional mean lands from the truth, a per-region view of quantification bias that the whole-brain metrics average away.
Every method, across every region
Every isolated dipole method's regional mean (one coloured dot per method) against the grey ground-truth band (mean ± std; tick = median), for all regions at once. Hover a dot or a legend entry to trace that method across every region; drag the plot sideways to pan the χ axis (calcification sits far to the left). Click a dot to open its page.
Statistics are computed over the full brain mask (a voxel an eroding method dropped counts as zero), so from the next full re-score the paired ground-truth values are identical for every run on a phantom; runs scored earlier used each run's own valid support and can differ slightly. Currently populated for isolated dipole runs; composed pipelines fill in at the next full re-score.
Accuracy vs compute
Isolated dipole inversion: structural similarity (xSIM, higher is better) against wall-clock runtime on a log scale. Green points are Pareto-optimal (nothing is both faster and more accurate); grey points are dominated. Some methods cost minutes for no accuracy gain.
What tuning buys
Isolated dipole NRMSE at default parameters (grey) versus tuned on the scoring phantom (green). Only a handful of methods are tunable (pretrained networks aren't), which is why defaults is the fair cross-method ranking. Tuning on the scored data is a mild form of overfitting.
Build your own comparison
Plot any two metrics against each other for any isolated stage. Points are methods, coloured by family; hover to identify, click to open a method. Runtime uses a log scale.
In silico · QSM Challenge 2.0 (2019)
The main dataset is the QSM Reconstruction Challenge 2.0 in-silico head phantom, forward-simulated to a BIDS multi-echo acquisition with qsm-forward. Its per-voxel susceptibility is known exactly, so QSM-CI scores each pipeline stage on its own (field mapping, background removal, dipole inversion) rather than only the final map. It backs the Full pipelines, Stages and Findings views.
Why a synthetic ground truth
QSM is an ill-posed dipole deconvolution, so different algorithms produce visibly different maps from the same phase. Judging them needs a true susceptibility map, which no in vivo acquisition provides, so every in vivo reference is a surrogate that biases the comparison: a COSMOS reconstruction (Liu et al., 2009) averages anisotropic susceptibility across orientations, and simple piecewise-constant phantoms (de Rochefort, 2008; Wharton, 2010) flatter the TV/TGV methods that reconstruct piecewise-constant maps by construction (Knoll, 2011). A ground truth is never neutral; its texture decides which assumptions get rewarded.
Challenge 2.0 (Marques et al., 2021) answered this with an in-silico head phantom: a realistic digital brain from a 0.64 mm 7 T acquisition, with an interhemispheric calcification, whose susceptibility is known exactly and from which a data-consistent gradient-echo acquisition is forward-simulated. It deliberately carries physiological texture so it does not reward piecewise-smooth solutions the way earlier phantoms did.
What QSM-CI scores, and what it leaves out
Because the total field, local field and tissue susceptibility are all known, QSM-CI scores each stage independently and reruns methods as they appear. The simulation carries realistic background fields from air/tissue interfaces and wrapped multi-echo phase, so field mapping and background removal genuinely have work to do, not just dipole inversion.
It stays idealised in known ways: susceptibility is an isotropic scalar (no white-matter anisotropy), there is no flow or motion, microstructure enters only through the χ-separation compartments, noise is complex Gaussian, and the phantom comes from a single subject. A QSM-CI ranking is a strong controlled test, not a guarantee of in vivo performance.
Reproduce the phantom
qsm-forward simulates the acquisition from the Challenge 2.0 tissue-property maps. Request those from the Donders repository (data.ru.nl/…/DSC_3015069.02_542, DOI 10.34973/m20r-jt17), then generate:
pip install qsm-forward
# 7 T, 4 echoes, 1 mm, fixed seed: matches the scored phantom
qsm-forward head ./qsm-challenge-2.0 bids_qsm --save-field--save-field emits the total field (fieldmap) and local field (fieldmap-local) that the isolated background-removal and dipole stages are scored against. Acquisition parameters (--B0, --TEs, --voxel-size, --peak-snr, --random-seed) default to the challenge phantom.
Citations
- Head phantom: Marques J.P., Meineke J., Milovic C., et al. “QSM reconstruction challenge 2.0: A realistic in silico head phantom for MRI data simulation and evaluation of susceptibility mapping procedures.” Magn Reson Med 2021;86(1):526–542. doi:10.1002/mrm.28716. Data: doi:10.34973/m20r-jt17.
- Challenge 2.0 design & results: QSM Challenge 2.0 Organization Committee. “QSM reconstruction challenge 2.0: Design and report of results.” Magn Reson Med 2021;86(3):1241–1255. doi:10.1002/mrm.28754.
- COSMOS reference: Liu T., Spincemaille P., de Rochefort L., et al. “Calculation of susceptibility through multiple orientation sampling (COSMOS).” Magn Reson Med 2009;61(1):196–204. doi:10.1002/mrm.21828.
- Piecewise-constant inversion: de Rochefort L., et al. Magn Reson Med 2008;60(4):1003–1009. doi:10.1002/mrm.21710. Wharton S., et al. Magn Reson Med 2010;63(5):1292–1304. doi:10.1002/mrm.22334.
- TGV regularisation: Knoll F., Bredies K., Pock T., Stollberger R. “Second order total generalized variation (TGV) for MRI.” Magn Reson Med 2011;65(2):480–491. doi:10.1002/mrm.22595.
- Forward simulation: qsm-forward (Stewart A., et al.).
| # | Algorithm | |
|---|---|---|
| GRE-native paper computing… DNF |
R2′ generators
Candidate R2′ substitutes, each scored directly against the phantom's true (spin-echo-derived) R2′. The best one (lowest detrended NRMSE) supplies R2′ to every pipeline above (named after each row's ⟵ arrow); the Findings tab compares how far each substitute carries every method.
| Generator | NRMSE | Detrended | Corr. | xSIM | Runtime |
|---|---|---|---|---|---|
| supplies R2′ above DNF |
Myelin (χ−) is harder to recover than iron (χ+)
Structural similarity (xSIM) of the recovered paramagnetic χ+ (iron) source versus the diamagnetic χ− (myelin / calcium) source. Every method drops from left to right, so the diamagnetic map is systematically the harder of the two.
Cross-contamination between the two sources
How much each source map bleeds into the other (regression slope, 0 = clean): iron leaking into the myelin map (x) versus myelin leaking into iron (y). The origin is perfect separation; distance from it is contamination that xSIM can miss.
What does dropping the spin-echo cost?
χ-separation wants R2′ = R2* − R2, and measuring R2 takes a separate (multi-echo) spin-echo acquisition. This figure shows each method's separation error (mean of the χ+ and χ− detrended NRMSE) when fed the true R2′ (● full MEGRE+MESE protocol) versus a GRE-derived substitute: ● R2PRIMEnet (learned R2*→R2′) or ● the fixed 0.52·R2* scaling. The gap is the price of GRE-only; methods that never read R2′ (◆) pay nothing by construction.
χ-separation (2026)
Susceptibility source separation splits the net susceptibility χ into a paramagnetic part χ⁺ (iron) and a diamagnetic part χ⁻ (myelin, calcium). Plain QSM measures only their signed sum, χtotal = χ⁺ + χ⁻, in which the two partly cancel, so it cannot tell them apart. The reversible relaxation rate R2′ responds instead to their combined magnitude, |χ⁺| + |χ⁻|. Having both the signed sum and the magnitude gives two independent constraints per voxel, which is what makes the sources recoverable. Each method is given χtotal, R2′, the local field and the multi-echo signals, and must produce χ⁺ and χ⁻, which are scored separately.
The GRE-only condition. Measuring R2′ = R2* − R2 needs a separate multi-echo spin-echo acquisition for R2, often unavailable in practice. The leaderboard’s MEGRE+MESE / GRE-only toggle compares each method fed the phantom’s true R2′ against the same method fed a GRE-derived substitute: an R2′ generator (a separate submission stage) estimates R2′ from the multi-echo GRE magnitude alone and is composed in front of every R2′-consuming method. The GRE-only leaderboard feeds every method the best-performing generator on that phantom (lowest detrended NRMSE against the true R2′, a scale-independent structural criterion, the “use the best available tool” policy a GRE-only protocol would follow); the full generator × method grid is scored regardless and compared on the Findings tab. Current generators are R2PRIMEnet (the learned R2*→R2′ conversion from the χ-sepnet pipeline, Kim et al., 2025, where composing it with χ-sepnet reproduces the paper’s χ-sepnet-R2* variant) and the fixed 0.52·R2* scaling evaluated by Ji et al. (2024) and Oliveira Assunção et al. (2026). Methods that never read R2′ (R2*-QSM, DECOMPOSE) are natively GRE-only and appear unchanged under both conditions. The Findings tab charts what the missing spin-echo costs each method.
The phantoms
Four phantoms are scored independently, all built on the Challenge 2.0 head model (Marques et al., 2021) via qsm-forward, and switchable with the leaderboard’s phantom selector. The first three reproduce the Susceptibility-Separation-Phantom of Ridani et al. (2026), verified voxelwise against the authors’ own ground truth (isotropic χ⁻ matches to < 10⁻⁵ ppm; per-ROI means to ≤ 5×10⁻⁴ ppm); the fourth, the Stewart multicompartment phantom (2026), is a QSM-CI extension that makes white-matter anisotropy recoverable from the signal itself.
| Phantom | Field / echoes | χ⁻ anisotropy | Dr | Signal |
|---|---|---|---|---|
| Ridani 3T isotropic | 3 T / 6 | none (constant Dr⁻) | scaled (theoretical) | mono-exponential |
| Ridani 3T anisotropic | 3 T / 6 | Dr⁻(θ) ∝ sin²θ | scaled (theoretical) | mono-exponential |
| Ridani 7T anisotropic | 7 T / 4 | Dr⁻(θ) ∝ sin²θ | scaled (theoretical) | mono-exponential |
| Stewart 7T multicompartment | 7 T / 8 | θ-encoded in the signal | fixed (Dr⁺ = 320) | hollow-cylinder 3-pool |
The three Ridani sets are noiseless and gradient-echo-only, with a mono-exponential white-matter decay in which the fibre angle θ modulates only the scalar R2′. They deliberately omit a fibre-orientation input so orientation cannot be handed to methods; a closed-form solver given the fibre angle separates the anisotropic phantom almost exactly.
The Stewart multicompartment phantom adds a matched spin-echo and replaces that mono-exponential decay with a three-pool hollow-cylinder white-matter signal model (Wharton & Bowtell, 2012). Its orientation-dependent pool interference produces a non-mono-exponential GRE decay whose beat encodes the fibre angle, so θ can be recovered from the signal itself rather than needing DTI. It is the minimal model in which fibre angle leaves any mark on the signal: mono-exponential decay carries none, and fuller models (fibre dispersion, exchange, diffusion) add realism at the cost of recoverability. The hollow-cylinder sits at the point where anisotropy is realistic enough to matter yet still cleanly recoverable.
Across all four, χ⁺ and χ⁻ are built per tissue and summed to give χtotal, so the net susceptibility shares the same anatomy as the In silico dataset. All are generated at 1 mm, chosen so that every method (including the memory-hungry deep-learning ones) runs within standard CI limits.
The forward model
Building χ⁺ and χ⁻
Each tissue is assigned literature χ⁺ and χ⁻ values (χ⁺lit, χ⁻lit), with within-tissue texture from the acquired R1 and R2* maps:
where the bar denotes a tissue-mean (\(\overline{R_2^*}\) and \(\overline{R_1}\) are the average R2* and R1 over that tissue, so each texture term is this voxel’s deviation from the tissue average), and the modulation coefficients a and b are set per tissue and per source. Summing the two gives χtotal.
χtotal is one of the two constraints from above; the other is R2′, which relates to the source magnitudes through a relaxivity Dr:
If Dr were a single known constant, the two would solve directly, a closed-form null solution:
This null solution is a useful baseline: it adds nothing beyond the standard assumption, so any method should at least match it. Real separation is harder because the true relaxivity is not a single known constant. Its main complication, described next, is that the myelin relaxivity depends on fibre orientation, which the field and R2′ alone cannot reveal; it also differs between iron and myelin, but that is a smaller, calibratable offset. That orientation dependence is why χ⁻ is the harder map to recover, and it is what the standard separation model (which assumes isotropic sources) leaves out.
Methods improve on this baseline by relaxing one of its assumptions: estimating the relaxivity rather than fixing it to a single constant, fitting the multi-echo signal decay directly instead of relying on the scalar R2′, or imposing spatial and learned priors that the voxel-wise closed form has none of.
White-matter anisotropy
Myelin is magnetically anisotropic, so in white matter χ⁻ depends on the angle θ between the nerve fibre and B₀, taken from the diffusion V1 eigenvector:
where Δχ is the myelin susceptibility anisotropy (susceptibility along the fibre minus across it) and χ₀ the orientation-independent part, both per tract. A single-orientation acquisition cannot observe θ, so the myelin relaxivity Dr⁻(θ) is unknown per voxel, which is what breaks the null solution for χ⁻. Iron is modelled as isotropic spheres and carries no such term.
R2′ and the relaxivity Dr
The true relaxivity is source-dependent, so R2′ is built as
Dr comes from the static-dephasing model (Yablonskiy & Haacke, 1994): spheres for χ⁺, cylinders for χ⁻, giving Dr⁺ = 107.84·B₀ for iron and Dr⁻ = 133.77·sin²θ·B₀ for myelin (755 and 936·sin²θ Hz/ppm at 7 T). These grow with field and run about 5× higher than measured in vivo. This is not an error, but the idealisation the model makes: static dephasing assumes perfectly static spins around ideal spheres/cylinders, whereas real tissue has diffusion averaging and heterogeneous microstructure that lower the effective relaxivity, so empirical calibrations (e.g. Shin et al.) land well below the theoretical value.
The three Ridani phantoms use the scaled setting: the static-dephasing values above, field-consistent and matching the reference exactly. The multicompartment phantom uses a fixed setting: both relaxivities are rescaled by a single factor so the iron value is Dr⁺ = 320 Hz/ppm (the field-scaled empirical calibration of Shin et al., 2021 at 7 T), field-independent, with the same χ⁺/χ⁻ ground truth.
Scaled vs fixed. Scaled is field-physical and reproduces Ridani; on such a phantom, methods that assume a lower empirical Dr carry a scale offset (an informative, real model error). Fixed anchors Dr⁺ to the empirical value so the score reflects separation rather than Dr calibration.
Design notes and caveats
Why not simply the physical split relaxivities everywhere. An earlier single χ-separation phantom shipped a noiseless R2′ built with two relaxivities (Dr⁺ = 114, Dr⁻ = 30 Hz/ppm). The 30 had no published basis, and every available separation method assumes one kernel, so a perfectly split noise-free R2′ was an internally inconsistent constraint that broke all of them (iLSQR scored xSIM ≈ 0.03) for a reason unrelated to reconstruction quality; on real data the R2′ noise and the R2 baseline buffer that mismatch. That phantom moved to the single empirical kernel Dr = 137 Hz/ppm (Shin et al., 2021). The current four phantoms keep the physically realistic split and address the same problem with the scaled / fixed settings above, and with the multicompartment phantom’s matched spin-echo, which lets a method estimate R2′ itself.
- χ⁺ and χ⁻ are scored separately. Every method recovers χ⁻ far worse than χ⁺ (in the thesis that informed this design, χ⁻ NRMSE was 78–87% across all seven methods), so a single combined score would hide where a method actually fails. χ⁻ is stored and shown as a positive magnitude.
- Learning-based versus iterative rankings deserve caution here. The deep-learning methods were trained on in vivo labels derived from multi-orientation acquisitions; a single-orientation synthetic phantom cannot reproduce that advantage, so on this data they do not necessarily beat the iterative methods, and that says little about in vivo performance.
- Known simplifications. Iron and myelin are taken as linear in χ (a post-mortem relation), R2′ follows the static-dephasing regime (valid for low source concentration and low diffusivity), and the Ridani phantoms carry no noise.
Several of these choices follow Paula Stoll, Development of a Deep Learning Framework for Iron and Myelin Mapping from QSM, MSc thesis, ETH Zürich / University of Queensland, 2025 (advisors S. Bollmann, A. Stewart), which extended the QSM forward model to χ-separation, benchmarked the same methods, and confirmed the xSIM constants used here (K1 = 0.01, K2 = 0.001).
The signal
The simulated multi-echo GRE magnitude carries the separated sources through R2′ rather than a single lumped R2*, so its decay is written directly in terms of χ⁺ and χ⁻:
so the effective decay rate is R2* = R2 + R2′. On the multicompartment phantom the matched spin-echo removes the reversible term and decays with R2 alone (\(S(\mathrm{TE}) = M_0\,e^{-\mathrm{TE}\,R_2}\)), letting a signal-domain method estimate R2′ = R2* − R2 itself instead of relying on the provided R2′ map.
The multicompartment signal (χ-sep MC)
On the Ridani phantoms the fibre angle θ enters only as a scalar that modulates the white-matter R2′, so θ is not recoverable from a single-orientation acquisition, and white-matter χ⁻ is unbeatable without an external DTI input. The multicompartment phantom removes that ceiling: white matter is modelled as three water pools (myelin, axonal, extra-axonal) with orientation-dependent frequency offsets (a hollow-cylinder model, Wharton & Bowtell, 2012). The pools beat against each other, giving a non-mono-exponential magnitude whose shape encodes θ and the myelin-water fraction, so a method can recover both from the GRE signal alone. This makes it the only phantom where anisotropy-aware separation can win without being handed the fibre orientation.
Scoring
The headline metric is the per-ROI MSPE (mean squared percentage error of ROI means) of Ridani et al. (2026), so the numbers are directly comparable to the paper: χ⁺ is averaged over the deep-gray-matter iron nuclei and cortical gray matter, χ⁻ over the white-matter fibre bundles of a co-registered atlas. Detrended NRMSE (scale-independent, voxelwise) and xSIM are also reported per source. For the anisotropic phantoms a MEV diagnostic (maximum error variation) reports how much the χ⁻ error falls from fibres parallel to B₀ to perpendicular. Each phantom also carries a null baseline (the closed-form solve with the phantom’s own assumptions) that any real method must beat.
Differences from the reference phantom
The three Ridani phantoms reproduce the χ⁺/χ⁻ construction, white-matter anisotropy and acquisition of Ridani et al. (2026), and their ground truth is verified voxelwise against the authors’ own OSF data. A few implementation choices differ:
- Orientation withheld. No fibre-angle (DTI) input is shipped, so a method cannot be handed the orientation that would make the anisotropic phantom nearly closed-form solvable. It must recover white-matter χ⁻ from the field and R2′ alone (or, on the multicompartment phantom, from the signal beat).
- Seeded noise. The reference’s white-matter texture noise is seeded here for reproducibility, so the ground truth is bit-reproducible; the authors’ unseeded realisation is matched only in distribution.
- Sign enforced. χ⁺ is clamped ≥ 0 and χ⁻ ≤ 0 (a source map cannot flip sign).
- Multicompartment (χ-sep MC only). The 7 T fixed-Dr phantom additionally simulates the hollow-cylinder white-matter signal and a matched spin-echo, a QSM-CI extension beyond the reference’s mono-exponential model.
The diffusion V1 eigenvector is stored mirrored (L-A-S) relative to the segmentation; qsm-forward flips it to register the fibres, which reproduces the reference’s published fibre-angle map exactly. The public reference code omits this flip.
Reproduce the phantoms
qsm-forward builds the phantoms from a data directory holding the Challenge 2.0 tissue-property maps and diffusion V1 eigenvector. Request those from the Donders repository (data.ru.nl/…/DSC_3015069.02_542, DOI 10.34973/m20r-jt17), and add the white-matter tract mask from the Susceptibility-Separation-Phantom repository (Ridani et al.). QSM-CI’s scripts/gen_chisep.py then generates each phantom from a preset:
# The fibre-bundle atlas: Masks/WM_fibers_seg.nii.gz from the Ridani et al. OSF project
# (osf.io/9xwhz, 168 KB, md5 06b1aec8b2a8a001e7e856c6c1fddcf0). Every shipped phantom carries it
# resampled as groundtruth/wm_rois.nii.gz — the label map the χ- per-ROI MSPE averages over — so
# --wm-rois-src is required to reproduce them fully.
ROIS=/path/to/WM_fibers_seg.nii.gz
gen_chisep.py --preset ridani-3t-iso --wm-rois-src $ROIS --out data/ridani-3t-iso
gen_chisep.py --preset ridani-3t-aniso --wm-rois-src $ROIS --out data/ridani-3t-aniso
gen_chisep.py --preset ridani-7t-aniso --wm-rois-src $ROIS --out data/ridani-7t-aniso
gen_chisep.py --multicompartment --dr-model fixed --dr-fixed 320 \
--tes 3,9,15,21,27,33,39,45 --wm-rois-src $ROIS --out data/chisep-mc # χ-sep MCEach preset sets the field, echo train, isotropy and relaxivity, simulates the signals and flattens the result into the dataset layout (public inputs/, scoring groundtruth/); point qsm-ci run straight at the NIfTIs.
Citations
- Separation phantom: Ridani D., De Leener B., Alonso-Ortiz E. “A realistic in-silico brain phantom for quantifying susceptibility anisotropy-induced error in susceptibility separation.” Magn Reson Med 2026. doi:10.1002/mrm.70468.
- Head model: Marques J.P., Meineke J., Milovic C., et al. “QSM reconstruction challenge 2.0: a realistic in silico head phantom for MRI data simulation and evaluation of susceptibility mapping procedures.” Magn Reson Med 2021;86(1):526–542. doi:10.1002/mrm.28716.
- χ-separation model & Dr: Shin H.G., Lee J., Yun Y.H., et al. “χ-separation: Magnetic susceptibility source separation toward iron and myelin mapping in the brain.” NeuroImage 2021;240:118371. doi:10.1016/j.neuroimage.2021.118371.
- Static-dephasing regime: Yablonskiy D.A., Haacke E.M. “Theory of NMR signal behavior in magnetically inhomogeneous tissues: the static dephasing regime.” Magn Reson Med 1994;32(6):749–763. doi:10.1002/mrm.1910320610.
- DECOMPOSE: Chen J., et al. “Decompose QSM to sub-voxel diamagnetic and paramagnetic components based on gradient-echo MRI data.” NeuroImage 2021. doi:10.1016/j.neuroimage.2021.118735.
- SUSEP-Net: Li, Gao, Sun, et al. “SUSEP-Net: simulation-supervised contrastive learning for susceptibility source separation via a sub-voxel QSM network.” arXiv:2506.13293, 2025. doi:10.48550/arXiv.2506.13293.
- WaveSep: Fang Z., Shin H.G., van Zijl P., Li X., Sulam J. “WaveSep: A Flexible Wavelet-Based Approach for Source Separation in Susceptibility Imaging.” MLCN, MICCAI 2023. doi:10.1007/978-3-031-44858-4_6.
- R2*-QSM: Dimov A.V., et al. “Magnetic Susceptibility Source Separation Solely from Gradient Echo Data: Histological Validation.” Tomography / J Neuroimaging 2022. doi:10.1111/jon.13014.
- Multicompartment white-matter signal: Wharton S., Bowtell R. “Fiber orientation-dependent white matter contrast in gradient echo MRI.” PNAS 2012;109(45):18559–18564. doi:10.1073/pnas.1211075109.
| # | Algorithm | vs COSMOS | vs STI χ33 | Runtime |
|---|---|---|---|---|
| DNF |
No in-vivo results yet.
COSMOS vs STI χ33: the harder reference
Each method's structural similarity against the two references. Every point sits below the diagonal: STI χ33 is a systematically harder target than COSMOS (roughly +20 NRMSE across the board), yet the two references rank methods almost identically. Coloured by family; hover for values.
Do NRMSE and HFEN agree?
The headline error metric (NRMSE) against the fine-detail metric (HFEN, the 2016 challenge's signature), both vs COSMOS. A high Spearman ρ means they reward the same methods; where a point strays off the trend, a method wins on bulk error but loses fine structure, or vice-versa.
Runtime vs accuracy
In-vivo structural similarity (xSIM vs COSMOS) against wall-clock runtime on a log scale. Green points are Pareto-optimal (nothing is both faster and more accurate); grey points are dominated.
Tuning rescues classical methods in vivo
Isolated dipole xSIM vs COSMOS at default parameters (grey) versus tuned on the in-vivo data (green), ordered by gain. At their defaults the classical iterative methods look far weaker than the pretrained networks; tuned, several leap most of the way to the best learned method (dashed). The payoff dwarfs the phantom's (medi and l1-qsm each gain roughly +0.3 xSIM) because their defaults are mistuned for real-data noise. Tuning on the scored data is a mild form of overfitting, and only some methods are tunable, so defaults stays the fair cross-method ranking.
In vivo · 2016 Reconstruction Challenge
This dataset is real acquired brain data from the 2016 QSM Reconstruction Challenge: a single-orientation 3 T scan of one volunteer. Each method inverts the ground-truth local tissue field to a susceptibility map, scored against two reference reconstructions rather than a true susceptibility (none exists in vivo). Only the dipole-inversion stage is scored, so the ranking reads as agreement with a reference, not accuracy against ground truth.
The challenge and its references
The 2016 Reconstruction Challenge (Langkammer et al., 2018) took the χ33 component of a susceptibility-tensor (STI) reconstruction as its reference, hoping to sidestep the orientation bias of COSMOS. Twenty-seven submissions were scored on RMSE, HFEN, SSIM and regional error. The lesson was as much about metrics as algorithms: tuned to minimise those errors, the winning maps were over-smoothed and lost fine structure, and near-identical scores hid visibly different maps. A follow-up (Milovic et al., 2020) then showed the χ33 reference itself did not match the measured phase once anisotropy, microstructure and noise were accounted for.
QSM-CI scores each method against both available references: the multi-orientation COSMOS reconstruction and the STI χ33 map. Source: the 2016 challenge (qsm.neuroimaging.at); the local tissue field and both reference χ maps are in ppm.
Citations
- 2016 Challenge report: Langkammer C., Schweser F., Shmueli K., et al. “Quantitative susceptibility mapping: Report from the 2016 reconstruction challenge.” Magn Reson Med 2018;79(3):1661–1673. doi:10.1002/mrm.26830.
- 2016 Challenge lessons: Milovic C., Tejos C., Acosta-Cabronero J., et al. “The 2016 QSM Challenge: Lessons learned and considerations for a future challenge design.” Magn Reson Med 2020;84(3):1624–1637. doi:10.1002/mrm.28185.
- COSMOS reference: Liu T., Spincemaille P., de Rochefort L., Kressler B., Wang Y. “Calculation of susceptibility through multiple orientation sampling (COSMOS).” Magn Reson Med 2009;61(1):196–204. doi:10.1002/mrm.21828.
How reproducible is each acquisition protocol, marginalised over reconstruction pipelines? Slope deviation |a−1| across all pipelines, as median (mean) [5th–95th percentile].
How wide is the spread, protocol by protocol?
Where the mean sits well right of the median, a minority of pipelines is dragging the average out: the protocol works for most pipelines and fails for a few.
Hover a dot for its pipeline, or click to open it.
| Protocol | Test–retest |a−1| | Inter-scanner |a−1| | BA bias | Pipelines |
|---|---|---|---|---|
Reproducibility results are not computed yet; the reconstruction matrix over all 23 acquisitions is pending.
pipeline(s) with in-silico xSIM ≥ , ranked by reproducibility.
How reproducible is each reconstruction pipeline across acquisitions of the same head? |a−1| over the pairwise per-ROI orthogonal fits in each comparison class. A pipeline that harmonizes well has a small inter-scanner slope deviation without giving up test–retest precision. Both statistics are shown in every cell — where they disagree sharply the pipeline is bimodal, reproducible on some protocols and not others.
Combination matrix
Dipole inversion ↓ × background removal → · best combinations top-left · bar length = in-silico accuracy · click a cell to inspect the pipeline.
Combined & end-to-end methods
These span background removal and dipole inversion in one step (or the whole reconstruction), so they have no cell on the axes above. Same colour scale.
No pipelines match the current field mapping and accuracy floor.
| # | Pipeline | Test–retest |a−1| | Inter-scanner |a−1| | vs bridge |a−1| | In-silico xSIM | Pairs |
|---|---|---|---|---|---|---|
Reproducibility results are not computed yet; the reconstruction matrix over all 23 acquisitions is pending.
Which stage decides reproducibility?
How much each stage's method choice reshuffles the full-pipeline |a−1| ranking, across every pipeline in the factorial (rank-based η², so one degenerate fit can't dominate). A larger share means that stage's choice matters more for reproducibility.
methods compared at this stage
Rank-based η² on each pipeline's median |a−1|: the share of the outcome ordering each stage's choice explains. Because the design is the complete factorial, these shares aren't confounded by the other stages.
Which methods are most reproducible, stage by stage?
For each method, its median |a−1| marginalised over every pipeline it appears in (dot), with the interquartile spread (bar) showing how much its partner stages swing it. Lower is more reproducible. Coloured by family; hover for the numbers.
Reproducible and accurate? Every method, both axes
Each point is a method, placed by its median |a−1| across every pipeline it appears in (y, lower is more reproducible) against its median in-silico (2019) xSIM over those same pipelines (x, higher is more accurate). Coloured by stage: field mapping, background removal, and the final χ-producing step (dipole inversion plus combined bf-removal+dipole and end-to-end methods). The bottom-right corner is the sweet spot: reproducible and accurate. Hover any point for its numbers. method(s) with a reproducibility fit but no in-silico run are omitted (no accuracy value to place them).
The reproducibility frontier
Each point is a pipeline: test–retest |a−1| (run-to-run precision) on x, inter-scanner |a−1| (Prisma↔Cima harmonizability) on y. Bottom-left is both precise and harmonizable. Pick a method to light up every pipeline that uses it; click any point to open it.
: median inter-scanner across pipelines that use it, vs for those that don't.
Where does harmonization fail? Region ×
Median |Δχ| for each nucleus, per , marginalised over the other pipeline stages and protocols. Green = the region's mean χ agrees across the selected acquisition pair; red = it doesn't. Columns sorted best→worst.
Does one scanner read a nucleus systematically higher?
Forest plot of Cima − Prisma mean-χ bias per region (matched protocol + run), median ± IQR across all pipelines. Deep-grey nuclei are bold. A dot off the dashed zero line is a systematic cross-scanner offset.
χ–iron calibration: do the scanners agree?
Median χ per deep-grey nucleus vs literature non-heme iron, one regression per scanner. The slope is χ-sensitivity to iron; matching slopes and overlapping points mean the scanners are harmonized on the iron axis. A tight linear χ–iron relation is expected (it’s one of QSM’s foundational validations), so this is a physical-plausibility & agreement check, not a precise calibration (n = 5 nuclei, and the iron values are fixed population constants, not this subject’s). Iron (µg Fe/g wet wt) from Hallgren & Sourander 1958: pallidum 213, putamen 133, caudate 93, ventral DC (SN/RN) ~190, thalamus ~55.
Which protocol shifts which region?
Each protocol's per-region mean-χ deviation from the cross-protocol consensus (same scanner/run/pipeline), averaged over everything else. This is a signed map: pale/near-white = the protocol agrees with the consensus (good); saturated blue = it reads that region lower, saturated red = higher. So both strong colours mean the protocol shifts that region, just in opposite directions (e.g. whether Pulseq pulls a nucleus up or down versus bridge/local).
Which stage drives deep-grey variability?
For 3-stage pipelines, swap the method at one stage while holding the other two plus the acquisition fixed, and measure the spread of deep-grey mean χ. A taller bar means that stage moves the clinically-relevant χ more, attributing reproducibility to unwrapping vs background removal vs dipole inversion.
Is inter-nucleus contrast preserved across scanners?
The χ difference between neighbouring nuclei on Cima (y) vs Prisma (x), one dot per pipeline × protocol × run. On the dashed identity line the contrast (e.g. pallidum−putamen) reproduces across scanners; scatter off it means the contrast a clinician reads changes with the scanner even when global χ looks fine.
Does the headline metric reflect the clinical regions?
Each pipeline's whole-brain inter-scanner slope-deviation (today's harmonization headline) vs its deep-grey regional χ agreement. A tight correlation means the global number is a good proxy; scatter means pipelines can look harmonized globally yet disagree in the nuclei that matter (or vice versa).
Harmonization · two scanners, four protocols, one head
One volunteer, scanned on two Siemens 3 T systems at MGH (bays 4–5: a MAGNETOM Prisma Fit and a MAGNETOM Cima.X), with four multi-echo GRE protocols, three runs each, for 23 acquisitions (one Pulseq offline run could not be reconstructed). There is no ground truth and none is possible; the question is reproducibility: how consistent are a pipeline's susceptibility values across runs, protocols and scanners of the same brain?
The dataset comes from the QSM harmonization initiative of the ISMRM Electro-Magnetic Tissue Properties Study Group (Costagli, Lancione, Bilgic, Robinson, Schweser, Li, and colleagues). Data were acquired by Berkin Bilgic; the analysis method (per-ROI orthogonal fits χB = a·χA + b between acquisitions, with |a−1| as the headline) follows Marta Lancione's two-pipeline analysis, which QSM-CI extends across its full pipeline matrix.
The four protocols
| Protocol | Sequence | TEs (ms) | TR / FA |
|---|---|---|---|
| Bridge | Siemens product GRE implementing the consensus protocol (GRAPPA R=2, 1 mm iso) | 5 / 11 / 17 / 23 / 29 | 35 ms / 15° |
| Local | Each site’s optimised protocol: Prisma is GRE 2×2, sum-of-squares; Cima is 3D-EPI, CAIPIRINHA 3×2 | Prisma: 7.04 / 13.4 / 19.76 / 26.12 Cima: 6.5 / 12 / 17.5 / 23 | 31 ms / 15° |
| Pulseq (online) | Vendor-neutral Pulseq consensus sequence (HarmonizedMRI/megre_label), Siemens ICE reconstruction | 5 / 11 / 17 / 23 / 29 | 35 ms / 15° |
| Pulseq (offline) | Same acquisition, offline GRAPPA reconstruction (removes the vendor recon from the chain) | 5 / 11 / 17 / 23 / 29 | 35 ms / 15° |
Echo times were recovered from the scanners’ exam-card PDFs and the Pulseq sequence source (the exported DICOM sidecars carry none, and the Pulseq exam cards show placeholder timings). The full recovered acquisition manifest ships with the dataset as protocols.json.
How it is evaluated
Every QSM-CI pipeline combination (field mapping × background removal × dipole inversion, plus single-step methods) reconstructs all acquisitions from raw phase. Each χ map is rigid-registered to a common target acquisition, referenced to its whole-brain mean, and reduced to per-ROI means over one SynthSeg segmentation (FreeSurfer mri_synthseg, which is license-free and emits standard FreeSurfer aseg labels (the deep-grey nuclei used in the regional Findings figures below), so ROI definitions cannot vary between acquisitions). Every acquisition pair then gets an orthogonal (total-least-squares) fit χB = a·χA + b: |a−1| measures scale disagreement, b the bias, and Bland–Altman limits the ROI-level scatter. Test–retest pairs bound what each protocol/pipeline could possibly achieve across scanners.
A fit is rejected outright when its slope is non-physical (|a−1| > , i.e. a outside [, ]), which would otherwise let one blown-up number dominate an average. This is rare: of pairs (%), and not a borderline call — rejected slopes run to a median |a−1| of and a maximum of . They are not spread evenly: and % are comparisons, , where a handful of solvers lose χ contrast entirely. A further pairs are missing because pipelines failed to reconstruct an acquisition. In total of possible pairs (%) survive.
Data availability
The dataset is on OSF at osf.io/gkemr: raw/ holds the original exports per scanner/protocol (with corrected sidecars, exam-card PDFs and the recovered-parameter manifest), and qsmci-inputs/ one prepacked zip per acquisition in the QSM-CI canonical layout (4D magnitude/phase in radians, SynthStrip mask, params.json).
Dipole inversion across datasets
The same isolated dipole-inversion methods, scored on the in-silico phantom (2019, true χ) and the in-vivo challenge (2016, COSMOS reference). Does a ranking on the idealised phantom carry over to real acquired data?
Phantom rank does not predict in-vivo rank
Each method's xSIM on the in-silico phantom (x) against its xSIM in vivo (y), coloured by family. The two are essentially uncorrelated: classical iterative and direct methods that win on the clean phantom collapse in vivo, while deep-learning methods hold up. Hover for values; click to open a method.
How methods move between datasets
Each method's rank on the in-silico phantom (left) linked to its rank in vivo (right), by xSIM. Lines coloured by family: deep-learning methods tend to climb left-to-right; classical methods fall.
Domain robustness by family
Mean xSIM per method family on each dataset. Deep-learning methods barely drop from phantom to in-vivo; classical iterative and direct methods lose far more: tuned to the idealised phantom, they generalise worse to real 3 T data.
The tuned optimum moves with the data
For every dipole method tuned on both datasets, its optimal regularisation strength on the in-silico phantom (left) versus in vivo (right), on a log scale. The in-vivo optimum sits roughly 10–30× higher almost everywhere: real 3 T data (with acquisition noise and forward-model mismatch) needs far heavier regularisation than the idealised phantom, so a single “recommended” value can't serve both. Coloured by family; hover for the shift.