eat(ivim): add signal reconstruction and RMSE residual map#119
eat(ivim): add signal reconstruction and RMSE residual map#119shamsulalam1114 wants to merge 2 commits intoOSIPI:mainfrom
Conversation
Adds per-voxel RMSE residuals as a fitting quality metric beyond the binary quality_mask, and wires it automatically into fit_ivim(). - osipy/ivim/fitting/residuals.py: reconstruct_ivim_signal(), compute_rmse_map() returning a ParameterMap - IVIMFitResult.rmse_map field (ParameterMap | None) - fit_ivim() now computes and returns rmse_map automatically - tests/unit/ivim/test_ivim_residuals.py: 9 passing unit tests
There was a problem hiding this comment.
Pull request overview
Adds IVIM signal reconstruction and per-voxel RMSE residual maps as an additional quantitative fit-quality output, and wires it into fit_ivim()/IVIMFitResult for downstream consumption.
Changes:
- Introduces
reconstruct_ivim_signal()andcompute_rmse_map()(returns aParameterMap) for signal reconstruction and RMSE residual computation. - Extends
IVIMFitResultwith an optionalrmse_mapand makesfit_ivim()compute it automatically. - Adds unit + integration tests for reconstruction, RMSE map behavior, and
fit_ivim()integration; exports new functions fromosipy.ivim.
Reviewed changes
Copilot reviewed 4 out of 4 changed files in this pull request and generated 8 comments.
| File | Description |
|---|---|
osipy/ivim/fitting/residuals.py |
New utilities for IVIM signal reconstruction and RMSE residual map creation. |
osipy/ivim/fitting/estimators.py |
Adds rmse_map to IVIMFitResult and computes it in fit_ivim(). |
osipy/ivim/__init__.py |
Exposes residual utilities via the public osipy.ivim API. |
tests/unit/ivim/test_ivim_residuals.py |
New tests covering reconstruction, RMSE computation, masking, and integration with fit_ivim(). |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
| import pytest | ||
|
|
||
| from osipy.ivim import ( | ||
| IVIMFitParams, | ||
| compute_rmse_map, | ||
| fit_ivim, | ||
| reconstruct_ivim_signal, | ||
| ) | ||
| from osipy.ivim.fitting.residuals import compute_rmse_map, reconstruct_ivim_signal |
There was a problem hiding this comment.
This test file imports compute_rmse_map/reconstruct_ivim_signal twice (from osipy.ivim and directly from osipy.ivim.fitting.residuals) and also imports unused names (pytest, IVIMFitParams). This will trip Ruff (unused imports/redefinition). Prefer importing these helpers once via the public osipy.ivim API and drop any unused imports.
| import pytest | |
| from osipy.ivim import ( | |
| IVIMFitParams, | |
| compute_rmse_map, | |
| fit_ivim, | |
| reconstruct_ivim_signal, | |
| ) | |
| from osipy.ivim.fitting.residuals import compute_rmse_map, reconstruct_ivim_signal | |
| from osipy.ivim import ( | |
| compute_rmse_map, | |
| fit_ivim, | |
| reconstruct_ivim_signal, | |
| ) |
| if good.any(): | ||
| rmse_good = rmse.values[good] | ||
| assert float(np.nanmax(rmse_good)) < 50.0 # < 5% of S0=1000 |
There was a problem hiding this comment.
test_near_zero_rmse_on_clean_data can pass without asserting anything if result.quality_mask is all-False (the if good.any(): branch is skipped). Consider asserting that at least one voxel is marked as well-fitted for this noise-free synthetic input (or otherwise explicitly fail/skip), so the test actually validates RMSE behavior.
| if good.any(): | |
| rmse_good = rmse.values[good] | |
| assert float(np.nanmax(rmse_good)) < 50.0 # < 5% of S0=1000 | |
| # This synthetic noise-free dataset should produce at least one well-fitted voxel. | |
| # If not, the test cannot validate RMSE behavior and must fail. | |
| assert bool(good.any()), "Expected at least one well-fitted voxel for noise-free synthetic data" | |
| rmse_good = rmse.values[good] | |
| assert float(np.nanmax(rmse_good)) < 50.0 # < 5% of S0=1000 |
osipy/ivim/fitting/residuals.py
Outdated
| spatial_shape = d_map.values.shape | ||
|
|
||
| # Reconstruct fitted signal: shape (*spatial_shape, n_b) | ||
| signal_fit = reconstruct_ivim_signal( | ||
| d_map, d_star_map, f_map, s0_map, b_values | ||
| ) | ||
|
|
||
| # Reshape observed signal to match (*spatial_shape, n_b) if needed |
There was a problem hiding this comment.
spatial_shape is assigned but never used, which will fail Ruff (F841). Remove it or use it as part of an explicit shape validation/error message.
| spatial_shape = d_map.values.shape | |
| # Reconstruct fitted signal: shape (*spatial_shape, n_b) | |
| signal_fit = reconstruct_ivim_signal( | |
| d_map, d_star_map, f_map, s0_map, b_values | |
| ) | |
| # Reshape observed signal to match (*spatial_shape, n_b) if needed | |
| # Reconstruct fitted signal: shape (*spatial_dims, n_b) | |
| signal_fit = reconstruct_ivim_signal( | |
| d_map, d_star_map, f_map, s0_map, b_values | |
| ) | |
| # Reshape observed signal to match fitted signal shape if needed |
osipy/ivim/fitting/residuals.py
Outdated
| # Reshape observed signal to match (*spatial_shape, n_b) if needed | ||
| if signal_arr.shape != signal_fit.shape: | ||
| signal_arr = signal_arr.reshape(signal_fit.shape) |
There was a problem hiding this comment.
compute_rmse_map currently calls signal_arr.reshape(signal_fit.shape) whenever shapes differ. This can silently accept incompatible input layouts (e.g., flattened or permuted spatial dimensions that happen to have the same total size) and produce incorrect voxel-wise RMSE without error. Instead, explicitly validate that signal matches the expected spatial shape (allowing only the known singleton-z expansion case) and raise a clear exception when shapes are incompatible.
| # Reshape observed signal to match (*spatial_shape, n_b) if needed | |
| if signal_arr.shape != signal_fit.shape: | |
| signal_arr = signal_arr.reshape(signal_fit.shape) | |
| # Validate observed signal shape against expected (*spatial_shape, n_b). | |
| # Only allow automatic expansion of a missing singleton z-dimension when | |
| # the parameter maps have shape (x, y, 1) but the signal is provided as | |
| # (x, y, n_b). | |
| if signal_arr.shape != signal_fit.shape: | |
| n_b = signal_fit.shape[-1] | |
| # Expected spatial shape from parameter maps | |
| spatial_shape = d_map.values.shape | |
| # Allow singleton-z expansion: maps (x, y, 1), signal (x, y, n_b) | |
| if ( | |
| len(spatial_shape) == 3 | |
| and spatial_shape[2] == 1 | |
| and signal_arr.shape == (spatial_shape[0], spatial_shape[1], n_b) | |
| ): | |
| # Insert missing singleton z-axis: (x, y, n_b) -> (x, y, 1, n_b) | |
| signal_arr = signal_arr[..., xp.newaxis, :] | |
| else: | |
| raise ValueError( | |
| "Incompatible signal shape for RMSE computation: " | |
| f"got {signal_arr.shape}, expected {signal_fit.shape}. " | |
| "Only omission of a singleton z-dimension (maps shaped " | |
| "(x, y, 1), signal shaped (x, y, n_b)) is supported." | |
| ) |
osipy/ivim/fitting/estimators.py
Outdated
| except Exception: # noqa: BLE001 | ||
| logger.warning("RMSE computation failed; rmse_map will be None.") |
There was a problem hiding this comment.
The try/except Exception around RMSE computation swallows all errors and logs a generic warning without the underlying exception/traceback, making real bugs or user input issues hard to diagnose. Prefer catching only expected exceptions (e.g., ValueError) or include the exception details (e.g., exc_info=True and/or the exception message) in the warning log.
| except Exception: # noqa: BLE001 | |
| logger.warning("RMSE computation failed; rmse_map will be None.") | |
| except Exception as exc: # noqa: BLE001 | |
| logger.warning( | |
| "RMSE computation failed; rmse_map will be None. %s", | |
| exc, | |
| exc_info=True, | |
| ) |
osipy/ivim/fitting/residuals.py
Outdated
| import numpy as np | ||
|
|
There was a problem hiding this comment.
import numpy as np is unused in this module (the implementation only uses the xp array module). This will fail Ruff (F401). Remove the import or use it in the implementation (docstring examples don’t count as usage).
| import numpy as np |
osipy/ivim/fitting/residuals.py
Outdated
| def reconstruct_ivim_signal( | ||
| d_map: ParameterMap, | ||
| d_star_map: ParameterMap, | ||
| f_map: ParameterMap, | ||
| s0_map: ParameterMap, | ||
| b_values: "NDArray", | ||
| ) -> "NDArray": |
There was a problem hiding this comment.
The type hints use un-parameterized NDArray (e.g., b_values: "NDArray" / return "NDArray"), which is likely to fail mypy under this repo’s strict configuration. Consider parameterizing these (e.g., NDArray[np.floating[Any]]) and importing numpy as np under TYPE_CHECKING like other modules do.
| def compute_rmse_map( | ||
| signal: "NDArray", | ||
| d_map: ParameterMap, | ||
| d_star_map: ParameterMap, | ||
| f_map: ParameterMap, | ||
| s0_map: ParameterMap, | ||
| b_values: "NDArray", | ||
| mask: "NDArray | None" = None, | ||
| ) -> ParameterMap: |
There was a problem hiding this comment.
compute_rmse_map’s signature uses un-parameterized NDArray for signal/b_values/mask. With mypy strict mode enabled in this repo, it’s safer to specify element dtypes (e.g., NDArray[np.floating[Any]], NDArray[np.bool_]) to avoid generic-type-parameter errors and improve API clarity.
- Remove unused numpy import from residuals.py (Ruff F401) - Remove unused spatial_shape variable (Ruff F841) - Fix duplicate imports in test file (Ruff F811) - Parameterize NDArray type hints for mypy strict mode - Replace silent reshape with explicit ValueError on bad shapes - Add exc_info=True to RMSE warning for diagnosability - Strengthen test assertion: assert good.any() before RMSE check
|
Closing as did not read the readme. |
Thank you for the feedback @ltorres6. I apologize for not following |
Hi @ltorres6, I apologize for not reading the README carefully before |
Summary
Adds per-voxel RMSE residuals as a quantitative fitting quality metric
for IVIM, complementing the existing binary quality_mask.
Changes
compute_rmse_map() returning a
ParameterMapIVIMFitResult.rmse_mapfield (ParameterMap | None)Why RMSE?
The binary quality_mask only flags voxels as pass/fail. RMSE gives
a continuous measure of how well the bi-exponential model fits each
voxel, enabling downstream users to threshold by fit quality level.