[BUG]: Correct file reading for AFNI outputs #409
9 changed files with 79 additions and 80 deletions
1
docs/changes/newsfragments/409.bugfix
Normal file
1
docs/changes/newsfragments/409.bugfix
Normal file
|
|
@ -0,0 +1 @@
|
||||||
|
Use correct file output suffices for AFNI-based markers by `Synchon Mandal`_
|
||||||
1
docs/changes/newsfragments/409.enh
Normal file
1
docs/changes/newsfragments/409.enh
Normal file
|
|
@ -0,0 +1 @@
|
||||||
|
Ease asset storage for AFNI-based markers by `Synchon Mandal`_
|
||||||
|
|
@ -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
|
)
|
||||||
|
|
|
||||||
|
|
@ -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"
|
||||||
|
|
|
||||||
|
|
@ -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,
|
||||||
|
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,
|
||||||
|
|
|
||||||
|
|
@ -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
|
||||||
|
|
|
||||||
|
|
@ -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
|
|
||||||
|
|
|
||||||
|
|
@ -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,
|
||||||
)
|
)
|
||||||
|
|
|
||||||
|
|
@ -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,
|
||||||
)
|
)
|
||||||
|
|
||||||
|
|
|
||||||
Loading…
Reference in a new issue
Is there any reason why we move to paths instead of in-memory data?