[BUG]: Correct file reading for AFNI outputs #409

Merged
synchon merged 6 commits from update/afni-impls into main 2024-12-03 10:54:18 +00:00
9 changed files with 79 additions and 80 deletions

View file

@ -0,0 +1 @@
Use correct file output suffices for AFNI-based markers by `Synchon Mandal`_

View file

@ -0,0 +1 @@
Ease asset storage for AFNI-based markers by `Synchon Mandal`_

View file

@ -50,7 +50,7 @@ class AFNIALFF(metaclass=Singleton):
@lru_cache(maxsize=None, typed=True) @lru_cache(maxsize=None, typed=True)
def compute( def compute(
self, self,
data: "Nifti1Image", input_path: Path,
highpass: float, highpass: float,
lowpass: float, lowpass: float,
tr: Optional[float], tr: Optional[float],
@ -59,8 +59,8 @@ class AFNIALFF(metaclass=Singleton):
Parameters Parameters
---------- ----------
data : 4D Niimg-like object input_path : pathlib.Path
Images to process. Path to the input data.
highpass : positive float highpass : positive float
Highpass cutoff frequency. Highpass cutoff frequency.
lowpass : positive float lowpass : positive float
@ -82,19 +82,17 @@ class AFNIALFF(metaclass=Singleton):
""" """
logger.debug("Creating cache for ALFF computation via AFNI") logger.debug("Creating cache for ALFF computation via AFNI")
# Create component-scoped tempdir # Create element-scoped tempdir
tempdir = WorkDirManager().get_tempdir(prefix="afni_alff+falff") element_tempdir = WorkDirManager().get_element_tempdir(
prefix="afni_lff"
# Save target data to a component-scoped tempfile )
nifti_in_file_path = tempdir / "input.nii" # needs to be .nii
nib.save(data, nifti_in_file_path)
# Set 3dRSFC command # Set 3dRSFC command
alff_falff_out_path_prefix = tempdir / "alff_falff" lff_out_path_prefix = element_tempdir / "output"
bp_cmd = [ bp_cmd = [
"3dRSFC", "3dRSFC",
f"-prefix {alff_falff_out_path_prefix.resolve()}", f"-prefix {lff_out_path_prefix.resolve()}",
f"-input {nifti_in_file_path.resolve()}", f"-input {input_path.resolve()}",
f"-band {highpass} {lowpass}", f"-band {highpass} {lowpass}",
"-no_rsfa -nosat -nodetrend", "-no_rsfa -nosat -nodetrend",
] ]
@ -104,49 +102,48 @@ class AFNIALFF(metaclass=Singleton):
# Call 3dRSFC # Call 3dRSFC
run_ext_cmd(name="3dRSFC", cmd=bp_cmd) run_ext_cmd(name="3dRSFC", cmd=bp_cmd)
# Create element-scoped tempdir so that the ALFF and fALFF maps are # Read header to get output suffix
# available later as nibabel stores file path reference for niimg = nib.load(input_path)
# loading on computation header = niimg.header
element_tempdir = WorkDirManager().get_element_tempdir( sform_code = header.get_sform(coded=True)[1]
prefix="afni_alff_falff" if sform_code == 4:
) output_suffix = "tlrc"
else:
output_suffix = "orig"
# Set params suffix
params_suffix = f"_{highpass}_{lowpass}_{tr}" params_suffix = f"_{highpass}_{lowpass}_{tr}"
# Convert alff afni to nifti # Convert alff afni to nifti
alff_afni_to_nifti_out_path = ( alff_nifti_out_path = (
element_tempdir / f"alff{params_suffix}_output.nii" element_tempdir / f"output_alff{params_suffix}.nii"
) # needs to be .nii ) # needs to be .nii
convert_alff_cmd = [ convert_alff_cmd = [
"3dAFNItoNIFTI", "3dAFNItoNIFTI",
f"-prefix {alff_afni_to_nifti_out_path.resolve()}", f"-prefix {alff_nifti_out_path.resolve()}",
f"{alff_falff_out_path_prefix}_ALFF+orig.BRIK", f"{lff_out_path_prefix}_ALFF+{output_suffix}.BRIK",
] ]
# Call 3dAFNItoNIFTI # Call 3dAFNItoNIFTI
run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_alff_cmd) run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_alff_cmd)
# Convert falff afni to nifti # Convert falff afni to nifti
falff_afni_to_nifti_out_path = ( falff_nifti_out_path = (
element_tempdir / f"falff{params_suffix}_output.nii" element_tempdir / f"output_falff{params_suffix}.nii"
) # needs to be .nii ) # needs to be .nii
convert_falff_cmd = [ convert_falff_cmd = [
"3dAFNItoNIFTI", "3dAFNItoNIFTI",
f"-prefix {falff_afni_to_nifti_out_path.resolve()}", f"-prefix {falff_nifti_out_path.resolve()}",
f"{alff_falff_out_path_prefix}_fALFF+orig.BRIK", f"{lff_out_path_prefix}_fALFF+{output_suffix}.BRIK",
] ]
# Call 3dAFNItoNIFTI # Call 3dAFNItoNIFTI
run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_falff_cmd) run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_falff_cmd)
# Load nifti # Load nifti
alff_data = nib.load(alff_afni_to_nifti_out_path) alff_data = nib.load(alff_nifti_out_path)
falff_data = nib.load(falff_afni_to_nifti_out_path) falff_data = nib.load(falff_nifti_out_path)
# Delete tempdir
WorkDirManager().delete_tempdir(tempdir)
return ( return (
alff_data, alff_data,
falff_data, falff_data,
alff_afni_to_nifti_out_path, alff_nifti_out_path,
falff_afni_to_nifti_out_path, falff_nifti_out_path,
) # type: ignore )

View file

@ -47,7 +47,7 @@ class JuniferALFF(metaclass=Singleton):
@lru_cache(maxsize=None, typed=True) @lru_cache(maxsize=None, typed=True)
def compute( def compute(
self, self,
data: "Nifti1Image", input_path: Path,
highpass: float, highpass: float,
lowpass: float, lowpass: float,
tr: Optional[float], tr: Optional[float],
@ -56,8 +56,8 @@ class JuniferALFF(metaclass=Singleton):
Parameters Parameters
---------- ----------
data : 4D Niimg-like object input_path : pathlib.Path
Images to process. Path to the input data.
highpass : positive float highpass : positive float
Highpass cutoff frequency. Highpass cutoff frequency.
lowpass : positive float lowpass : positive float
@ -80,9 +80,10 @@ class JuniferALFF(metaclass=Singleton):
logger.debug("Creating cache for ALFF computation via junifer") logger.debug("Creating cache for ALFF computation via junifer")
# Get scan data # Get scan data
niimg_data = data.get_fdata().copy() niimg = nib.load(input_path)
niimg_data = niimg.get_fdata().copy()
if tr is None: if tr is None:
tr = float(data.header["pixdim"][4]) # type: ignore tr = float(niimg.header["pixdim"][4]) # type: ignore
logger.info(f"`tr` not provided, using `tr` from header: {tr}") logger.info(f"`tr` not provided, using `tr` from header: {tr}")
# Bandpass the data within the lowpass and highpass cutoff freqs # Bandpass the data within the lowpass and highpass cutoff freqs
@ -120,19 +121,17 @@ class JuniferALFF(metaclass=Singleton):
# Calculate ALFF # Calculate ALFF
alff = numerator / np.sqrt(niimg_data.shape[-1]) alff = numerator / np.sqrt(niimg_data.shape[-1])
alff_data = nimg.new_img_like( alff_data = nimg.new_img_like(
ref_niimg=data, ref_niimg=niimg,
data=alff, data=alff,
) )
falff_data = nimg.new_img_like( falff_data = nimg.new_img_like(
ref_niimg=data, ref_niimg=niimg,
data=falff, data=falff,
) )
# Create element-scoped tempdir so that the ALFF and fALFF maps are # Create element-scoped tempdir
# available later as nibabel stores file path reference for
# loading on computation
element_tempdir = WorkDirManager().get_element_tempdir( element_tempdir = WorkDirManager().get_element_tempdir(
prefix="junifer_alff+falff" prefix="junifer_lff"
) )
output_alff_path = element_tempdir / "output_alff.nii.gz" output_alff_path = element_tempdir / "output_alff.nii.gz"
output_falff_path = element_tempdir / "output_falff.nii.gz" output_falff_path = element_tempdir / "output_falff.nii.gz"

View file

@ -146,7 +146,7 @@ class ALFFBase(BaseMarker):
estimator = JuniferALFF() estimator = JuniferALFF()
# Compute ALFF + fALFF # Compute ALFF + fALFF
alff, falff, alff_path, falff_path = estimator.compute( # type: ignore alff, falff, alff_path, falff_path = estimator.compute( # type: ignore
data=input_data["data"], input_path=input_data["path"],
highpass=self.highpass, highpass=self.highpass,
fraimondo commented 2024-12-02 13:28:33 +00:00 (Migrated from github.com)

Is there any reason why we move to paths instead of in-memory data?

Is there any reason why we move to paths instead of in-memory data?
synchon commented 2024-12-02 14:03:48 +00:00 (Migrated from github.com)

Passing the path is easier to reason about inside the method and cheaper to hash for the LRU cache.

Passing the path is easier to reason about inside the method and cheaper to hash for the LRU cache.
lowpass=self.lowpass, lowpass=self.lowpass,
tr=self.tr, tr=self.tr,

View file

@ -69,7 +69,9 @@ def test_ALFFSpheres(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None:
# Fit transform marker on data # Fit transform marker on data
output = marker.fit_transform(element_data) output = marker.fit_transform(element_data)
assert "Creating cache" in caplog.text # Tests for ALFFParcels run before this with the same data and that
# should create the cache
assert "Calculating ALFF and fALFF" in caplog.text
# Get BOLD output # Get BOLD output
assert "BOLD" in output assert "BOLD" in output

View file

@ -50,7 +50,7 @@ class AFNIReHo(metaclass=Singleton):
@lru_cache(maxsize=None, typed=True) @lru_cache(maxsize=None, typed=True)
def compute( def compute(
self, self,
data: "Nifti1Image", input_path: Path,
nneigh: int = 27, nneigh: int = 27,
neigh_rad: Optional[float] = None, neigh_rad: Optional[float] = None,
neigh_x: Optional[float] = None, neigh_x: Optional[float] = None,
@ -65,8 +65,8 @@ class AFNIReHo(metaclass=Singleton):
Parameters Parameters
---------- ----------
data : 4D Niimg-like object input_path : pathlib.Path
Images to process. Path to the input data.
nneigh : {7, 19, 27}, optional nneigh : {7, 19, 27}, optional
Number of voxels in the neighbourhood, inclusive. Can be: Number of voxels in the neighbourhood, inclusive. Can be:
@ -128,19 +128,17 @@ class AFNIReHo(metaclass=Singleton):
""" """
logger.debug("Creating cache for ReHo computation via AFNI") logger.debug("Creating cache for ReHo computation via AFNI")
# Create component-scoped tempdir # Create element-scoped tempdir
tempdir = WorkDirManager().get_tempdir(prefix="afni_reho") element_tempdir = WorkDirManager().get_element_tempdir(
prefix="afni_reho"
# Save target data to a component-scoped tempfile )
nifti_in_file_path = tempdir / "input.nii" # needs to be .nii
nib.save(data, nifti_in_file_path)
# Set 3dReHo command # Set 3dReHo command
reho_out_path_prefix = tempdir / "reho" reho_out_path_prefix = element_tempdir / "output"
reho_cmd = [ reho_cmd = [
"3dReHo", "3dReHo",
f"-prefix {reho_out_path_prefix.resolve()}", f"-prefix {reho_out_path_prefix.resolve()}",
f"-inset {nifti_in_file_path.resolve()}", f"-inset {input_path.resolve()}",
] ]
# Check ellipsoidal / cuboidal volume arguments # Check ellipsoidal / cuboidal volume arguments
if neigh_rad: if neigh_rad:
@ -164,28 +162,28 @@ class AFNIReHo(metaclass=Singleton):
# Call 3dReHo # Call 3dReHo
run_ext_cmd(name="3dReHo", cmd=reho_cmd) run_ext_cmd(name="3dReHo", cmd=reho_cmd)
# Create element-scoped tempdir so that the ReHo map is # Read header to get output suffix
# available later as nibabel stores file path reference for niimg = nib.load(input_path)
# loading on computation header = niimg.header
element_tempdir = WorkDirManager().get_element_tempdir( sform_code = header.get_sform(coded=True)[1]
prefix="afni_reho" if sform_code == 4:
) output_suffix = "tlrc"
else:
output_suffix = "orig"
# Convert afni to nifti # Convert afni to nifti
reho_afni_to_nifti_out_path = ( reho_nifti_out_path = (
element_tempdir / "output.nii" # needs to be .nii element_tempdir / "output.nii" # needs to be .nii
) )
convert_cmd = [ convert_cmd = [
"3dAFNItoNIFTI", "3dAFNItoNIFTI",
f"-prefix {reho_afni_to_nifti_out_path.resolve()}", f"-prefix {reho_nifti_out_path.resolve()}",
f"{reho_out_path_prefix}+orig.BRIK", f"{reho_out_path_prefix}+{output_suffix}.BRIK",
] ]
# Call 3dAFNItoNIFTI # Call 3dAFNItoNIFTI
run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_cmd) run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_cmd)
# Load nifti # Load nifti
output_data = nib.load(reho_afni_to_nifti_out_path) output_data = nib.load(reho_nifti_out_path)
# Delete tempdir return output_data, reho_nifti_out_path
WorkDirManager().delete_tempdir(tempdir)
return output_data, reho_afni_to_nifti_out_path # type: ignore

View file

@ -48,15 +48,15 @@ class JuniferReHo(metaclass=Singleton):
@lru_cache(maxsize=None, typed=True) @lru_cache(maxsize=None, typed=True)
def compute( def compute(
self, self,
data: "Nifti1Image", input_path: Path,
nneigh: int = 27, nneigh: int = 27,
) -> tuple["Nifti1Image", Path]: ) -> tuple["Nifti1Image", Path]:
"""Compute ReHo map. """Compute ReHo map.
Parameters Parameters
---------- ----------
data : 4D Niimg-like object input_path : pathlib.Path
Images to process. Path to the input data.
nneigh : {7, 19, 27, 125}, optional nneigh : {7, 19, 27, 125}, optional
Number of voxels in the neighbourhood, inclusive. Can be: Number of voxels in the neighbourhood, inclusive. Can be:
@ -89,7 +89,8 @@ class JuniferReHo(metaclass=Singleton):
logger.debug("Creating cache for ReHo computation via junifer") logger.debug("Creating cache for ReHo computation via junifer")
# Get scan data # Get scan data
niimg_data = data.get_fdata() niimg = nib.load(input_path)
niimg_data = niimg.get_fdata().copy()
# Get scan dimensions # Get scan dimensions
n_x, n_y, n_z, _ = niimg_data.shape n_x, n_y, n_z, _ = niimg_data.shape
@ -119,7 +120,7 @@ class JuniferReHo(metaclass=Singleton):
# after #299 is merged # after #299 is merged
# Calculate whole brain mask # Calculate whole brain mask
mni152_whole_brain_mask = nmask.compute_brain_mask( mni152_whole_brain_mask = nmask.compute_brain_mask(
target_img=data, target_img=niimg,
threshold=0.5, threshold=0.5,
mask_type="whole-brain", mask_type="whole-brain",
) )
@ -227,7 +228,7 @@ class JuniferReHo(metaclass=Singleton):
# Create new image like target image # Create new image like target image
output_data = nimg.new_img_like( output_data = nimg.new_img_like(
ref_niimg=data, ref_niimg=niimg,
data=reho_map, data=reho_map,
copy_header=False, copy_header=False,
) )

View file

@ -125,7 +125,7 @@ class ReHoBase(BaseMarker):
estimator = JuniferReHo() estimator = JuniferReHo()
# Compute reho # Compute reho
reho_map, reho_map_path = estimator.compute( # type: ignore reho_map, reho_map_path = estimator.compute( # type: ignore
data=input_data["data"], input_path=input_data["path"],
**reho_params, **reho_params,
) )