From c8e4cce8e1f75a9991053312a29a737fc5e69f26 Mon Sep 17 00:00:00 2001 From: Kenya Otsuka Date: Wed, 24 Jun 2026 04:36:06 +0000 Subject: [PATCH 1/2] refactor: replace nipy image I/O with nibabel and deprecate nipy-based GLM paradigm Work toward dropping the unmaintained nipy dependency (KamitaniLab/bdpy#68). - Replace nipy.load_image with nibabel.load and img.coordmap.affine with img.affine in load_mri, load_epi, and fmriprep volume loading. These paths no longer use nipy; behavior is unchanged (nibabel's affine matches nipy's coordmap.affine for the supported images). - Emit a DeprecationWarning from make_paradigm, which still relies on nipy's paradigm objects, pointing users toward nilearn. nipy is retained as a dependency so the GLM path keeps working until a later release removes it. Co-Authored-By: Claude Opus 4.8 --- bdpy/mri/fmriprep.py | 13 ++++++------- bdpy/mri/glm.py | 9 +++++++++ bdpy/mri/load_epi.py | 6 +++--- bdpy/mri/load_mri.py | 8 ++++---- 4 files changed, 22 insertions(+), 14 deletions(-) diff --git a/bdpy/mri/fmriprep.py b/bdpy/mri/fmriprep.py index 6567a69..70ca5ff 100644 --- a/bdpy/mri/fmriprep.py +++ b/bdpy/mri/fmriprep.py @@ -11,7 +11,6 @@ from collections import OrderedDict import nibabel -import nipy import numpy as np import pandas as pd @@ -406,17 +405,17 @@ def __load_volume(self): - Returns voxel ijk indexes (3 x voxel) - Data, xyz, and ijk are flattened by Fortran-like index order """ - img = nipy.load_image(self.__dpath) + img = nibabel.load(self.__dpath) data = img.get_fdata() if data.ndim == 4: data = data.reshape(-1, data.shape[-1], order='F').T i_len, j_len, k_len, t = img.shape - affine = np.delete(np.delete(img.coordmap.affine, 3, axis=0), 3, axis=1) + affine = np.delete(np.delete(img.affine, 3, axis=0), 3, axis=1) elif data.ndim == 3: data = data.flatten(order='F') i_len, j_len, k_len = img.shape - affine = img.coordmap.affine + affine = img.affine else: raise ValueError('Invalid shape.') @@ -838,17 +837,17 @@ def __load_mri(fpath): - Returns voxel ijk indexes (3 x voxel) - Data, xyz, and ijk are flattened by Fortran-like index order """ - img = nipy.load_image(fpath) + img = nibabel.load(fpath) data = img.get_fdata() if data.ndim == 4: data = data.reshape(-1, data.shape[-1], order='F').T i_len, j_len, k_len, t = img.shape - affine = np.delete(np.delete(img.coordmap.affine, 3, axis=0), 3, axis=1) + affine = np.delete(np.delete(img.affine, 3, axis=0), 3, axis=1) elif data.ndim == 3: data = data.flatten(order='F') i_len, j_len, k_len = img.shape - affine = img.coordmap.affine + affine = img.affine else: raise ValueError('Invalid shape.') diff --git a/bdpy/mri/glm.py b/bdpy/mri/glm.py index 2916cac..c775c96 100644 --- a/bdpy/mri/glm.py +++ b/bdpy/mri/glm.py @@ -2,6 +2,7 @@ import csv +import warnings import numpy as np from nipy.modalities.fmri.experimental_paradigm import ( @@ -45,6 +46,14 @@ def make_paradigm(event_files, num_vols, tr=2., cond_col=2, label_col=None, regr run_regressors : nuisance regressors for runs run_regressors_label : labels for the run regressors """ + warnings.warn( + "make_paradigm depends on nipy, which is unmaintained. nipy support " + "will be removed in a future release of bdpy. Consider migrating to " + "nilearn (e.g. nilearn.glm.first_level.make_first_level_design_matrix).", + DeprecationWarning, + stacklevel=2, + ) + onset = [] duration = [] conds = [] diff --git a/bdpy/mri/load_epi.py b/bdpy/mri/load_epi.py index 7e1d87b..95ec1d4 100644 --- a/bdpy/mri/load_epi.py +++ b/bdpy/mri/load_epi.py @@ -6,7 +6,7 @@ import itertools as itr -import nipy +import nibabel import numpy as np @@ -35,7 +35,7 @@ def load_epi(datafiles): print("Loading %s" % df) # Load an EPI image - img = nipy.load_image(df) + img = nibabel.load(df) xyz = _check_xyz(xyz, img) data_list.append(np.array(img.get_fdata().flatten(), dtype=np.float64)) @@ -49,7 +49,7 @@ def _check_xyz(xyz, img): """Check voxel xyz consistency.""" # Use the image shape metadata; avoid get_fdata() here, which would read the # whole volume just to check consistency (load_epi reads the data later). - xyz_current = _get_xyz(img.coordmap.affine, img.shape[:3]) + xyz_current = _get_xyz(img.affine, img.shape[:3]) if xyz.size == 0: xyz = xyz_current diff --git a/bdpy/mri/load_mri.py b/bdpy/mri/load_mri.py index 772d6b0..e0b5b26 100644 --- a/bdpy/mri/load_mri.py +++ b/bdpy/mri/load_mri.py @@ -1,7 +1,7 @@ """load_mri""" -import nipy +import nibabel import numpy as np @@ -13,17 +13,17 @@ def load_mri(fpath): - Returns voxel ijk indexes (3 x voxel) - Data, xyz, and ijk are flattened by Fortran-like index order """ - img = nipy.load_image(fpath) + img = nibabel.load(fpath) data = img.get_fdata() if data.ndim == 4: data = data.reshape(-1, data.shape[-1], order='F').T i_len, j_len, k_len, t = img.shape - affine = np.delete(np.delete(img.coordmap.affine, 3, axis=0), 3, axis=1) + affine = np.delete(np.delete(img.affine, 3, axis=0), 3, axis=1) elif data.ndim == 3: data = data.flatten(order='F') i_len, j_len, k_len = img.shape - affine = img.coordmap.affine + affine = img.affine else: raise ValueError('Invalid shape.') From 91bcb2b8d573b719bb79e7949b5775504f2bff70 Mon Sep 17 00:00:00 2001 From: Kenya Otsuka Date: Wed, 24 Jun 2026 05:13:10 +0000 Subject: [PATCH 2/2] fix: handle 4D affine in load_mri and BrainData volume loading The 4D code path reduced the affine to a 3x3 matrix before applying it to homogeneous (4 x N) voxel indices, raising a shape-mismatch ValueError whenever a 4D image was loaded. This bug predates the nipy->nibabel swap (nipy's coordmap.affine and nibabel's affine are both 4x4, so both backends hit it). Use the full 4x4 affine for the 4D path, matching the 3D path. Fixed in the active loaders: load_mri.load_mri and fmriprep.BrainData.__load_volume. The latter is the real volume path: create_bdata_fmriprep treats BrainData.data.shape[0] as the number of volumes (time points), which only holds for the (T, N) shape the 4D branch produces. The unused, name-private fmriprep.__load_mri and __get_xyz are left untouched to keep the change minimal. Adds test_load_mri_4d and test_braindata_load_volume_4d covering the data/xyz/ijk shapes and values for a 4D volume with a non-trivial affine. Co-Authored-By: Claude Opus 4.8 --- bdpy/mri/fmriprep.py | 2 +- bdpy/mri/load_mri.py | 2 +- tests/test_mri.py | 88 ++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 90 insertions(+), 2 deletions(-) diff --git a/bdpy/mri/fmriprep.py b/bdpy/mri/fmriprep.py index 70ca5ff..7e786d0 100644 --- a/bdpy/mri/fmriprep.py +++ b/bdpy/mri/fmriprep.py @@ -411,7 +411,7 @@ def __load_volume(self): if data.ndim == 4: data = data.reshape(-1, data.shape[-1], order='F').T i_len, j_len, k_len, t = img.shape - affine = np.delete(np.delete(img.affine, 3, axis=0), 3, axis=1) + affine = img.affine elif data.ndim == 3: data = data.flatten(order='F') i_len, j_len, k_len = img.shape diff --git a/bdpy/mri/load_mri.py b/bdpy/mri/load_mri.py index e0b5b26..94c31b1 100644 --- a/bdpy/mri/load_mri.py +++ b/bdpy/mri/load_mri.py @@ -19,7 +19,7 @@ def load_mri(fpath): if data.ndim == 4: data = data.reshape(-1, data.shape[-1], order='F').T i_len, j_len, k_len, t = img.shape - affine = np.delete(np.delete(img.affine, 3, axis=0), 3, axis=1) + affine = img.affine elif data.ndim == 3: data = data.flatten(order='F') i_len, j_len, k_len = img.shape diff --git a/tests/test_mri.py b/tests/test_mri.py index bee53ae..e1e27fb 100644 --- a/tests/test_mri.py +++ b/tests/test_mri.py @@ -47,6 +47,94 @@ def test_load_mri_3d(self) -> None: self.assertEqual(xyz.shape, (3, arr.size)) self.assertEqual(ijk.shape, (3, arr.size)) + def test_load_mri_4d(self) -> None: + """Test load_mri on a 4D volume. + + Guards a bug in the 4D path where the affine was reduced to 3x3 before + being applied to homogeneous (4 x N) voxel indices, raising a shape + mismatch. The 4D path must use the full 4x4 affine, just like the 3D + path, and return one row per time point. + """ + import os + import tempfile + + import nibabel + + i_len, j_len, k_len, t_len = 2, 3, 4, 5 + n_vox = i_len * j_len * k_len + arr = np.arange(n_vox * t_len, dtype=np.float32).reshape( + i_len, j_len, k_len, t_len) + affine = np.array([[2., 0., 0., -10.], + [0., 2., 0., -12.], + [0., 0., 2., -8.], + [0., 0., 0., 1.]]) + + with tempfile.TemporaryDirectory() as tmpdir: + fpath = os.path.join(tmpdir, 'test4d.nii.gz') + nibabel.save(nibabel.Nifti1Image(arr, affine=affine), fpath) + data, xyz, ijk = bmr.load_mri(fpath) # type: ignore + + # Data is (sample x voxel): one sample per time point, voxels flattened + # in Fortran order. + self.assertEqual(data.shape, (t_len, n_vox)) + np.testing.assert_array_equal( + data, arr.reshape(-1, t_len, order='F').T) + + # ijk are the voxel indices in Fortran order. + expected_ijk = np.array(np.unravel_index( + np.arange(n_vox), (i_len, j_len, k_len), order='F')) + np.testing.assert_array_equal(ijk, expected_ijk) + + # xyz are the world coordinates from the full 4x4 affine. + ijk_b = np.vstack([expected_ijk, np.ones((1, n_vox))]) + expected_xyz = np.dot(affine, ijk_b)[:3] + np.testing.assert_allclose(xyz, expected_xyz) + self.assertEqual(xyz.shape, (3, n_vox)) + + def test_braindata_load_volume_4d(self) -> None: + """Test BrainData volume loading on a 4D volume. + + Guards the same 4D affine bug as test_load_mri_4d in the duplicated + loader inside fmriprep.BrainData.__load_volume. The surrounding + create_bdata_fmriprep code treats ``data.shape[0]`` as the number of + volumes (time points), so the 4D path must return ``(T, N)`` data and + spatial xyz from the full 4x4 affine. + """ + import os + import tempfile + + import nibabel + + from bdpy.mri.fmriprep import BrainData + + i_len, j_len, k_len, t_len = 2, 3, 4, 5 + n_vox = i_len * j_len * k_len + arr = np.arange(n_vox * t_len, dtype=np.float32).reshape( + i_len, j_len, k_len, t_len) + affine = np.array([[2., 0., 0., -10.], + [0., 2., 0., -12.], + [0., 0., 2., -8.], + [0., 0., 0., 1.]]) + + with tempfile.TemporaryDirectory() as tmpdir: + fpath = os.path.join(tmpdir, 'bold4d.nii.gz') + nibabel.save(nibabel.Nifti1Image(arr, affine=affine), fpath) + brain = BrainData(fpath, dtype='volume') + + # data.shape[0] must be the number of time points (volumes). + self.assertEqual(brain.data.shape, (t_len, n_vox)) + np.testing.assert_array_equal( + brain.data, arr.reshape(-1, t_len, order='F').T) + + expected_ijk = np.array(np.unravel_index( + np.arange(n_vox), (i_len, j_len, k_len), order='F')) + np.testing.assert_array_equal(brain.index, expected_ijk) + + ijk_b = np.vstack([expected_ijk, np.ones((1, n_vox))]) + expected_xyz = np.dot(affine, ijk_b)[:3] + np.testing.assert_allclose(brain.xyz, expected_xyz) + self.assertEqual(brain.xyz.shape, (3, n_vox)) + def test_get_roiflag_pass0002(self) -> None: """Test for get_roiflag (pass case 0002).""" roi_xyz = [np.array([[1, 2, 3],