dmipy-fit: Diffusion Microstructure Imaging in Python

dmipy-fit is a Python toolbox for biophysical modelling of diffusion MRI data. Given an acquisition scheme and a multi-compartment tissue model, it fits the model parameters voxel-by-voxel and returns interpretable microstructure maps — axon density and dispersion, diffusivities, cell size, volume fractions, T2, and more.

It models the diffusion signal plus the transverse-relaxation contrasts that ride on it — T2 and surface relaxivity — as composable, occupancy-gated factors that attach to any compartment. Every analytical model is validated effect-by-effect against the dmipy-sim Monte-Carlo simulator (the two share one free-waveform sequence/substrate interface; the dependency is one-directional, fit → sim).

from dmipy_fit.signal_models.gaussian_models import G1Ball
from dmipy_fit.signal_models.cylinder_models import C1Stick
from dmipy_fit.core.modeling_framework import MultiCompartmentModel

ball_stick = MultiCompartmentModel([G1Ball(), C1Stick()])
fit = ball_stick.fit(scheme, data, solver="jax")     # whole slice on GPU
intra_fraction = fit.fitted_parameters["partial_volume_1"]

The signal model

dmipy-fit models the measured signal as a sum of compartments, each carrying its diffusion attenuation times its (transverse) relaxation contrasts:

$$ S ;=; S_0 \sum_i f_i ; \underbrace{E^{\mathrm{diff}}i(b)}{\text{diffusion}} ; \underbrace{e^{-\mathrm{TE}/T_{2,i}}}_{T_2} ; \underbrace{\hat B^{\mathrm{surf}}i(\mathrm{TE})}{\text{surface relaxivity}} $$

Magnetisation is treated as fully transverse (ideal instantaneous pulses). Each non-diffusion effect is an occupancy-gated factor — a multiplicative attenuation that attaches to any compartment via OccupancyGatedModel, so any subset composes by listing more factors.

What's here

• Signal models (signal_models/) — sticks/cylinders, sphere, Gaussian (ball/zeppelin), plane, capped cylinder, tissue-response; b-tensor / free-waveform aware.
• Occupancy-gated factors (signal_models/attenuation.py) — OccupancyGatedModel + TransverseRelaxation, IntraPoreSurfaceRelaxivity, ExteriorSurfaceRelaxivity. Relaxation and surface relaxivity as composable add-ons to any compartment.
• Fitting framework (core/) — MultiCompartmentModel, spherical-mean and spherical-harmonics frameworks, fitted-model properties.
• Distributions (distributions/) — Watson / Bingham dispersion, Gamma diameter.
• Optimizers — brute2fine, MIX, multi-tissue NNLS; CSD (Tournier / cvxpy / OSQP-JAX); IVIM / impulsed / 3-tissue custom optimizers.
• GPU fitting (jax/) — diffusion signal models, multicompartment, CSD, DTI, fractions, vmap_fit; whole-slice fitting in seconds with solver="jax".
• White matter (white_matter/):
  • build_white_matter_model() — a decoupled, diffusion-only canonical WM model (stick + zeppelin + stuck-myelin dot + optional CSF ball), each compartment gated with T2 + surface-relaxivity. Because the intra-pore and exterior surface factors differ per compartment, surface relaxivity introduces a b-independent signal weighting between intra and extra (true spin fractions ≠ apparent fractions) — a physical effect that plain stick+zeppelin models omit.
  • t2_spectrum_mwf() — standard regularised NNLS T2-spectrum myelin-water fraction (Whittall–MacKay / DECAES). Estimate MWF directly from a multi-echo T2 decay, including one generated by a dmipy-sim CPMG-style multi-TE series.
• Rician noise model, b-tensor / free-waveform schemes, and a Monte-Carlo bridge to dmipy-sim.
• Flagship example (examples/flagship_canonical_wm/) — the canonical WM signal computed both ways (this analytical model and the dmipy-sim Monte-Carlo), agreeing for diffusion and surface relaxivity; fully reproducible (regenerates the cached MC bit-exactly). The clearest entry point to what the two engines are for.

Installation

pip install dmipy-fit                    # core only
pip install "dmipy-fit[jax]"             # + JAX GPU fitting
pip install "dmipy-fit[jax,data,viz]"    # + bundled-data loaders + plotting
pip install "dmipy-fit[all]"             # everything

Myelin water fraction from a Monte-Carlo CPMG decay

dmipy_sim.simulate_cpmg returns the full echo train from a single walk; feed it to the NNLS T2-spectrum estimator:

import numpy as np
import dmipy_sim as ds
from dmipy_fit.white_matter.mwf import t2_spectrum_mwf

geom = ds.MyelinatedCylinder(inner_radius=2.5e-6, outer_radius=3.57e-6, orientation=(0, 0, 1),
    D_intra=1.7e-9, D_myelin_radial=0.1e-9, D_myelin_tangential=0.5e-9, D_extra=1.7e-9,
    T2_intra=0.080, T2_myelin=0.015, T2_extra=0.080)

TE, n_echoes = 8e-3, 32
wf = ds.cpmg(n_echoes=n_echoes, TE=TE, G_magnitude=0.0, bvecs=[[1, 0, 0]], n_t_per_echo=60)
S = np.asarray(ds.simulate_cpmg(40000, None, wf, geom)).ravel()      # (n_echoes,), one walk
echo_times = np.arange(1, n_echoes + 1) * TE
mwf, T2_grid, spectrum = t2_spectrum_mwf(S / S[0], echo_times)

Testing

pytest -q                       # analytical + numerical (dipy/MISST references)
pytest -q -m "not slow"         # skip the heavy GPU battery

Tests assert against analytical results, MISST/dipy references, or analytic ↔ Monte-Carlo parity against dmipy-sim.

License

Dual-licensed: GNU AGPL-3.0 for open-source use, or a commercial license for proprietary/closed use. See LICENSE and LICENSING.md (commercial: rutger.fick@dmrai-lab.org).