From 1069948ef63bd2975a673fe0749196e4f7608481 Mon Sep 17 00:00:00 2001 From: Taylor Salo Date: Fri, 7 Jun 2019 12:03:42 +0200 Subject: [PATCH 1/4] Change default thresholding to bi-sided. We should treat positively and negatively weighted voxels separately when performing cluster-extent thresholding. --- tedana/utils.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/tedana/utils.py b/tedana/utils.py index dbfb9b4bc..ab663b950 100644 --- a/tedana/utils.py +++ b/tedana/utils.py @@ -233,7 +233,7 @@ def get_spectrum(data: np.array, tr: float = 1.0): def threshold_map(img, min_cluster_size, threshold=None, mask=None, - binarize=True, sided='two'): + binarize=True, sided='bi'): """ Cluster-extent threshold and binarize image. @@ -250,10 +250,10 @@ def threshold_map(img, min_cluster_size, threshold=None, mask=None, Boolean array for masking resultant data array. Default is None. binarize : bool, optional Default is True. - sided : {'two', 'one', 'bi'}, optional + sided : {'bi', 'two', 'one'}, optional How to apply thresholding. One-sided thresholds on the positive side. Two-sided thresholds positive and negative values together. Bi-sided - thresholds positive and negative values separately. Default is 'two'. + thresholds positive and negative values separately. Default is 'bi'. """ if not isinstance(img, np.ndarray): arr = img.get_data() From c1e7229bf2c24128f920d314ce4511e4389700d9 Mon Sep 17 00:00:00 2001 From: Taylor Salo Date: Fri, 7 Jun 2019 13:56:11 +0200 Subject: [PATCH 2/4] Draft apply_xforms function. --- tedana/utils.py | 20 +++++++++++ tedana/workflows/post_aroma.py | 65 ++++++++++++++++++++++++++++++++++ 2 files changed, 85 insertions(+) create mode 100644 tedana/workflows/post_aroma.py diff --git a/tedana/utils.py b/tedana/utils.py index ab663b950..3066c8829 100644 --- a/tedana/utils.py +++ b/tedana/utils.py @@ -14,6 +14,26 @@ LGR = logging.getLogger(__name__) +def apply_xforms(in_file, out_file, xforms): + """ + Apply some transformations to an input file. + """ + from nipype.interfaces.ants.resampling import ApplyTransforms + if not isinstance(xforms, list): + xforms = [xforms] + + for xform in xforms: + assert op.isfile(xform) + + at = ApplyTransforms() + at.inputs.input_image = in_file + at.inputs.reference_image = in_file + at.inputs.transforms = xforms + at.inputs.invert_transform_flags = [False for _ in range(len(xforms))] + at.inputs.output_image = out_file + at.run() + + def load_image(data): """ Takes input `data` and returns a sample x time array diff --git a/tedana/workflows/post_aroma.py b/tedana/workflows/post_aroma.py new file mode 100644 index 000000000..d0453b4ee --- /dev/null +++ b/tedana/workflows/post_aroma.py @@ -0,0 +1,65 @@ +""" +Apply AROMA to tedana output folder. +""" +import os.path as op +import logging + +import aroma +import numpy as np +import nibabel as nib + +from tedana.utils import apply_xforms + +LGR = logging.getLogger(__name__) + + +def aroma_workflow(tedana_dir, motpars_file, t_r, xforms=None): + """ + Apply AROMA to a tedana output directory to further identify components + that are motion-related from the original ICA. Overwrites original + component table with updated component classifications. + + Parameters + ---------- + tedana_dir : :obj:`str` + Directory containing tedana run outputs. + motpars_file : :obj:`str` + Motion parameters file. + t_r : :obj:`float` + Repetition time, in seconds. + xforms : :obj:`list` of :obj:`str`, optional + List of transforms to apply to warp component maps in tedana_dir into + standard space. + """ + betas_file = op.join(tedana_dir, 'betas_OC.nii') + mix_file = op.join(tedana_dir, 'meica_mix.1D') + betas_file2 = op.join(tedana_dir, 'betas_OC_std.nii') + apply_xforms(betas_file, betas_file2, xforms) + betas_data = nib.load(betas_file2) + # variance normalize betas to make pseudo-z-values + z_data = betas_data / np.std(betas_data, axis=-1) + z_thresh_idx = np.abs(z_data) > np.median(np.abs(z_data)) + z_thresh_data = z_data * z_thresh_idx + clf_df = aroma.run_aroma(mix_file, z_thresh_data, motpars_file, t_r) + clf_df.loc[clf_df['classification'] == 'rejected'] + comptable = pd.read_csv(op.join(tedana_dir, 'comp_table_ica.txt'), + sep='\t', index_col='component') + comptable.rename(columns={'classification': 'original_classification', + 'rationale': 'original_rationale'}, + inplace=True) + comptable = pd.merge(comptable, clf_df, on='component', how='outer') + # AROMA is only used for rejection, so original classifications for any + # components not rejected by AROMA from tedana. + LGR.info('Overriding original classifications for AROMA-rejected ' + 'components. All other components will retain their original ' + 'classifications.') + not_rej_idx = comptable.loc[comptable['classification'] != 'rejected'].index.values + comptable.loc[not_rej_idx, 'classification'] = comptable.loc[ + not_rej_idx, 'original_classification'] + comptable.loc[not_rej_idx, 'rationale'] = comptable.loc[ + not_rej_idx, 'original_rationale'] + + # Overwrite original component table + LGR.info('Overwriting original component table.') + comptable.to_csv(op.join(tedana_dir, 'comp_table_ica.txt'), sep='\t', + index_label='component') From 2276a8df8d00acc0538de4e3cb962a40bd1b172f Mon Sep 17 00:00:00 2001 From: Taylor Salo Date: Fri, 7 Jun 2019 13:58:24 +0200 Subject: [PATCH 3/4] Revert. Ugh. --- tedana/utils.py | 20 ----------- tedana/workflows/post_aroma.py | 65 ---------------------------------- 2 files changed, 85 deletions(-) delete mode 100644 tedana/workflows/post_aroma.py diff --git a/tedana/utils.py b/tedana/utils.py index 3066c8829..ab663b950 100644 --- a/tedana/utils.py +++ b/tedana/utils.py @@ -14,26 +14,6 @@ LGR = logging.getLogger(__name__) -def apply_xforms(in_file, out_file, xforms): - """ - Apply some transformations to an input file. - """ - from nipype.interfaces.ants.resampling import ApplyTransforms - if not isinstance(xforms, list): - xforms = [xforms] - - for xform in xforms: - assert op.isfile(xform) - - at = ApplyTransforms() - at.inputs.input_image = in_file - at.inputs.reference_image = in_file - at.inputs.transforms = xforms - at.inputs.invert_transform_flags = [False for _ in range(len(xforms))] - at.inputs.output_image = out_file - at.run() - - def load_image(data): """ Takes input `data` and returns a sample x time array diff --git a/tedana/workflows/post_aroma.py b/tedana/workflows/post_aroma.py deleted file mode 100644 index d0453b4ee..000000000 --- a/tedana/workflows/post_aroma.py +++ /dev/null @@ -1,65 +0,0 @@ -""" -Apply AROMA to tedana output folder. -""" -import os.path as op -import logging - -import aroma -import numpy as np -import nibabel as nib - -from tedana.utils import apply_xforms - -LGR = logging.getLogger(__name__) - - -def aroma_workflow(tedana_dir, motpars_file, t_r, xforms=None): - """ - Apply AROMA to a tedana output directory to further identify components - that are motion-related from the original ICA. Overwrites original - component table with updated component classifications. - - Parameters - ---------- - tedana_dir : :obj:`str` - Directory containing tedana run outputs. - motpars_file : :obj:`str` - Motion parameters file. - t_r : :obj:`float` - Repetition time, in seconds. - xforms : :obj:`list` of :obj:`str`, optional - List of transforms to apply to warp component maps in tedana_dir into - standard space. - """ - betas_file = op.join(tedana_dir, 'betas_OC.nii') - mix_file = op.join(tedana_dir, 'meica_mix.1D') - betas_file2 = op.join(tedana_dir, 'betas_OC_std.nii') - apply_xforms(betas_file, betas_file2, xforms) - betas_data = nib.load(betas_file2) - # variance normalize betas to make pseudo-z-values - z_data = betas_data / np.std(betas_data, axis=-1) - z_thresh_idx = np.abs(z_data) > np.median(np.abs(z_data)) - z_thresh_data = z_data * z_thresh_idx - clf_df = aroma.run_aroma(mix_file, z_thresh_data, motpars_file, t_r) - clf_df.loc[clf_df['classification'] == 'rejected'] - comptable = pd.read_csv(op.join(tedana_dir, 'comp_table_ica.txt'), - sep='\t', index_col='component') - comptable.rename(columns={'classification': 'original_classification', - 'rationale': 'original_rationale'}, - inplace=True) - comptable = pd.merge(comptable, clf_df, on='component', how='outer') - # AROMA is only used for rejection, so original classifications for any - # components not rejected by AROMA from tedana. - LGR.info('Overriding original classifications for AROMA-rejected ' - 'components. All other components will retain their original ' - 'classifications.') - not_rej_idx = comptable.loc[comptable['classification'] != 'rejected'].index.values - comptable.loc[not_rej_idx, 'classification'] = comptable.loc[ - not_rej_idx, 'original_classification'] - comptable.loc[not_rej_idx, 'rationale'] = comptable.loc[ - not_rej_idx, 'original_rationale'] - - # Overwrite original component table - LGR.info('Overwriting original component table.') - comptable.to_csv(op.join(tedana_dir, 'comp_table_ica.txt'), sep='\t', - index_label='component') From e459168be28094384d47f32bb79d439d87a08c01 Mon Sep 17 00:00:00 2001 From: Taylor Salo Date: Fri, 13 Dec 2019 15:34:30 -0500 Subject: [PATCH 4/4] Add returns section to threshold_map docstring. --- tedana/utils.py | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/tedana/utils.py b/tedana/utils.py index 7a385c422..d3954c4bb 100644 --- a/tedana/utils.py +++ b/tedana/utils.py @@ -259,6 +259,11 @@ def threshold_map(img, min_cluster_size, threshold=None, mask=None, How to apply thresholding. One-sided thresholds on the positive side. Two-sided thresholds positive and negative values together. Bi-sided thresholds positive and negative values separately. Default is 'bi'. + + Returns + ------- + clust_thresholded : (M) :obj:`numpy.ndarray` + Cluster-extent thresholded (and optionally binarized) map. """ if not isinstance(img, np.ndarray): arr = img.get_data()