Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

[REF, DOC] Document PAID combination method #264

Merged
merged 6 commits into from
Apr 23, 2019
Merged

[REF, DOC] Document PAID combination method #264

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
ceb3359
Rename ste to paid combination method.
tsalo Apr 20, 2019
b086261
Add note box about PAID to documentation.
tsalo Apr 20, 2019
5492777
Update docs/approach.rst
emdupre Apr 23, 2019
e4b9dc5
Update docs/approach.rst
emdupre Apr 23, 2019
f9e7bb4
Newline for each sentence.
tsalo Apr 23, 2019
9eb250d
[skip ci] Address review.
tsalo Apr 23, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
18 changes: 15 additions & 3 deletions docs/approach.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,9 +46,9 @@ calculated below, each voxel's values are only calculated from the first :math:`
echoes, where :math:`n` is the value for that voxel in the adaptive mask.

.. note::
``tedana`` allows users to provide their own mask. The adaptive mask will
be computed on this explicit mask, and may reduce it further based on the
data.
``tedana`` allows users to provide their own mask.
The adaptive mask will be computed on this explicit mask, and may reduce
it further based on the data.
If a mask is not provided, ``tedana`` runs `nilearn.masking.compute_epi_mask`_
on the first echo's data to derive a mask prior to adaptive masking.
The workflow does this because the adaptive mask generation function
Expand Down Expand Up @@ -140,6 +140,17 @@ of the other echoes (which it is).
.. image:: /_static/10_optimal_combination_timeseries.png
:align: center

.. note::
tsalo marked this conversation as resolved.
Show resolved Hide resolved
An alternative method for optimal combination that
does not use :math:`T_{2}^*`, is the parallel-acquired inhomogeneity
desensitized (PAID) ME-fMRI combination method (`Poser et al., 2006`_).
This method specifically assumes that noise in the acquired echoes is "isotopic and
homogeneous throughout the image," meaning it should be used on smoothed data.
As we do not recommend performing tedana denoising on smoothed data,
we discourage using PAID within the tedana workflow.
We do, however, make it accessible as an alternative combination method
in the t2smap workflow.

TEDPCA
``````
The next step is to identify and temporarily remove Gaussian (thermal) noise
Expand Down Expand Up @@ -223,3 +234,4 @@ robust PCA.

.. _nilearn.masking.compute_epi_mask: https://nilearn.github.io/modules/generated/nilearn.masking.compute_epi_mask.html
.. _Power et al. (2018): http://www.pnas.org/content/early/2018/02/07/1720985115.short
.. _Poser et al., 2006: https://onlinelibrary.wiley.com/doi/full/10.1002/mrm.20900
35 changes: 18 additions & 17 deletions tedana/combine.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -50,11 +50,13 @@ def _combine_t2s(data, tes, ft2s):


@due.dcite(Doi('10.1002/mrm.20900'),
description='STE method of combining data across echoes using just '
description='PAID method of combining data across echoes using just '
'SNR/signal and TE.')
def _combine_ste(data, tes):
def _combine_paid(data, tes):
"""
Combine data across echoes using SNR/signal and TE.
Combine data across echoes using SNR/signal and TE via the
parallel-acquired inhomogeneity desensitized (PAID) ME-fMRI combination
method.

Parameters
----------
Expand Down Expand Up @@ -90,9 +92,9 @@ def make_optcom(data, tes, mask, t2s=None, combmode='t2s', verbose=True):
t2s : (S [x T]) :obj:`numpy.ndarray` or None, optional
Estimated T2* values. Only required if combmode = 't2s'.
Default is None.
combmode : {'t2s', 'ste'}, optional
How to combine data. Either 'ste' or 't2s'. If 'ste', argument 't2s' is
not required. Default is 't2s'.
combmode : {'t2s', 'paid'}, optional
How to combine data. Either 'paid' or 't2s'. If 'paid', argument 't2s'
is not required. Default is 't2s'.
verbose : :obj:`bool`, optional
Whether to print status updates. Default is True.

Expand Down Expand Up @@ -127,32 +129,31 @@ def make_optcom(data, tes, mask, t2s=None, combmode='t2s', verbose=True):
'voxels/samples: {0} != {1}'.format(mask.shape[0],
data.shape[0]))

if combmode not in ['t2s', 'ste']:
raise ValueError("Argument 'combmode' must be either 't2s' or 'ste'")
if combmode not in ['t2s', 'paid']:
raise ValueError("Argument 'combmode' must be either 't2s' or 'paid'")
elif combmode == 't2s' and t2s is None:
raise ValueError("Argument 't2s' must be supplied if 'combmode' is "
"set to 't2s'.")
elif combmode == 'ste' and t2s is not None:
LGR.warning("Argument 't2s' is not required if 'combmode' is 'ste'. "
elif combmode == 'paid' and t2s is not None:
LGR.warning("Argument 't2s' is not required if 'combmode' is 'paid'. "
"'t2s' array will not be used.")

data = data[mask, :, :] # mask out empty voxels/samples
tes = np.array(tes)[np.newaxis, ...] # (1 x E) array_like

if t2s is not None:
if combmode == 'paid':
LGR.info('Optimally combining data with parallel-acquired inhomogeneity '
'desensitized (PAID) method')
combined = _combine_paid(data, tes)
else:
if t2s.ndim == 1:
msg = 'Optimally combining data with voxel-wise T2 estimates'
else:
msg = ('Optimally combining data with voxel- and volume-wise T2 '
'estimates')
t2s = t2s[mask, ..., np.newaxis] # mask out empty voxels/samples

if verbose:
LGR.info(msg)

if combmode == 'ste':
combined = _combine_ste(data, tes)
else:
LGR.info(msg)
combined = _combine_t2s(data, tes, t2s)

combined = unmask(combined, mask)
Expand Down
18 changes: 9 additions & 9 deletions tedana/decomposition/eigendecomp.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -165,7 +165,7 @@ def kundu_tedpca(comptable, n_echos, kdaw, rdaw, stabilize=False):


def tedpca(catd, OCcatd, combmode, mask, t2s, t2sG,
ref_img, tes, method='mle', ste=-1, kdaw=10., rdaw=1.,
ref_img, tes, method='mle', source_tes=-1, kdaw=10., rdaw=1.,
verbose=False):
"""
Use principal components analysis (PCA) to identify and remove thermal
Expand All @@ -177,9 +177,9 @@ def tedpca(catd, OCcatd, combmode, mask, t2s, t2sG,
Input functional data
OCcatd : (S x T) array_like
Optimally combined time series data
combmode : {'t2s', 'ste'} str
combmode : {'t2s', 'paid'} str
How optimal combination of echos should be made, where 't2s' indicates
using the method of Posse 1999 and 'ste' indicates using the method of
using the method of Posse 1999 and 'paid' indicates using the method of
Poser 2006
mask : (S,) array_like
Boolean mask array
Expand All @@ -193,7 +193,7 @@ def tedpca(catd, OCcatd, combmode, mask, t2s, t2sG,
Dimensionality augmentation weight for Rho calculations
method : {'mle', 'kundu', 'kundu-stabilize'}, optional
Method with which to select components in TEDPCA. Default is 'mle'.
ste : :obj:`int` or :obj:`list` of :obj:`int`, optional
source_tes : :obj:`int` or :obj:`list` of :obj:`int`, optional
Which echos to use in PCA. Values -1 and 0 are special, where a value
of -1 will indicate using the optimal combination of the echos
and 0 will indicate using all the echos. A list can be provided
Expand Down Expand Up @@ -260,17 +260,17 @@ def tedpca(catd, OCcatd, combmode, mask, t2s, t2sG,
"""

n_samp, n_echos, n_vols = catd.shape
ste = np.array([int(ee) for ee in str(ste).split(',')])
source_tes = np.array([int(ee) for ee in str(source_tes).split(',')])

if len(ste) == 1 and ste[0] == -1:
if len(source_tes) == 1 and source_tes[0] == -1:
LGR.info('Computing PCA of optimally combined multi-echo data')
d = OCcatd[mask, :][:, np.newaxis, :]
elif len(ste) == 1 and ste[0] == 0:
elif len(source_tes) == 1 and source_tes[0] == 0:
LGR.info('Computing PCA of spatially concatenated multi-echo data')
d = catd[mask, ...]
else:
LGR.info('Computing PCA of echo #%s' % ','.join([str(ee) for ee in ste]))
d = np.stack([catd[mask, ee, :] for ee in ste - 1], axis=1)
LGR.info('Computing PCA of echo #%s' % ','.join([str(ee) for ee in source_tes]))
d = np.stack([catd[mask, ee, :] for ee in source_tes - 1], axis=1)

eim = np.squeeze(eimask(d))
d = np.squeeze(d[eim])
Expand Down
4 changes: 2 additions & 2 deletions tedana/model/fit.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -39,9 +39,9 @@ def fitmodels_direct(catd, mmix, mask, t2s, t2s_full, tes, combmode, ref_img,
estimate using the first two echoes.
tes : list
List of echo times associated with `catd`, in milliseconds
combmode : {'t2s', 'ste'} str
combmode : {'t2s', 'paid'} str
How optimal combination of echos should be made, where 't2s' indicates
using the method of Posse 1999 and 'ste' indicates using the method of
using the method of Posse 1999 and 'paid' indicates using the method of
Poser 2006
ref_img : str or img_like
Reference image to dictate how outputs are saved to disk
Expand Down
10 changes: 5 additions & 5 deletions tedana/tests/test_combine.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,15 +27,15 @@ def test__combine_t2s():
assert comb.shape == (n_voxels, n_trs)


def test__combine_ste():
def test__combine_paid():
"""
Test tedana.model.combine._combine_ste
Test tedana.model.combine._combine_paid
"""
np.random.seed(0)
n_voxels, n_echos, n_trs = 20, 3, 10
data = np.random.random((n_voxels, n_echos, n_trs))
tes = np.array([[10, 20, 30]]) # 1 x E
comb = combine._combine_ste(data, tes)
comb = combine._combine_paid(data, tes)
assert comb.shape == (n_voxels, n_trs)


Expand All @@ -62,9 +62,9 @@ def test_make_optcom():
assert comb.shape == (n_voxels, n_trs)

# STE with erroneously included T2* argument
comb = combine.make_optcom(data, tes, mask, t2s=t2s, combmode='ste')
comb = combine.make_optcom(data, tes, mask, t2s=t2s, combmode='paid')
assert comb.shape == (n_voxels, n_trs)

# Normal STE call
comb = combine.make_optcom(data, tes, mask, t2s=None, combmode='ste')
comb = combine.make_optcom(data, tes, mask, t2s=None, combmode='paid')
assert comb.shape == (n_voxels, n_trs)
8 changes: 4 additions & 4 deletions tedana/tests/test_t2smap.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -67,13 +67,13 @@ def test_basic_t2smap2(self):
def test_basic_t2smap3(self):
"""
A very simple test, to confirm that t2smap creates output
files when combmode is set to 'ste'.
files when combmode is set to 'paid'.
"""
data_dir = get_test_data_path()
data = [op.join(data_dir, 'echo1.nii.gz'),
op.join(data_dir, 'echo2.nii.gz'),
op.join(data_dir, 'echo3.nii.gz')]
workflows.t2smap_workflow(data, [14.5, 38.5, 62.5], combmode='ste',
workflows.t2smap_workflow(data, [14.5, 38.5, 62.5], combmode='paid',
fitmode='all', label='t2smap')
out_dir = 'TED.echo1.t2smap'

Expand All @@ -93,15 +93,15 @@ def test_basic_t2smap3(self):
def test_basic_t2smap4(self):
"""
A very simple test, to confirm that t2smap creates output
files when combmode is set to 'ste' and fitmode is set to 'ts'.
files when combmode is set to 'paid' and fitmode is set to 'ts'.

Not sure why this fails.
"""
data_dir = get_test_data_path()
data = [op.join(data_dir, 'echo1.nii.gz'),
op.join(data_dir, 'echo2.nii.gz'),
op.join(data_dir, 'echo3.nii.gz')]
workflows.t2smap_workflow(data, [14.5, 38.5, 62.5], combmode='ste',
workflows.t2smap_workflow(data, [14.5, 38.5, 62.5], combmode='paid',
fitmode='ts', label='t2smap')
out_dir = 'TED.echo1.t2smap'

Expand Down
12 changes: 6 additions & 6 deletions tedana/workflows/t2smap.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -67,9 +67,9 @@ def _get_parser():
optional.add_argument('--combmode',
dest='combmode',
action='store',
choices=['t2s', 'ste'],
choices=['t2s', 'paid'],
help=('Combination scheme for TEs: '
't2s (Posse 1999, default), ste (Poser)'),
't2s (Posse 1999, default), paid (Poser)'),
default='t2s')
optional.add_argument('--label',
dest='label',
Expand Down Expand Up @@ -110,8 +110,8 @@ def t2smap_workflow(data, tes, mask=None, fitmode='all', combmode='t2s',
'all' means that the model is fit, per voxel, across all timepoints.
'ts' means that the model is fit, per voxel and per timepoint.
Default is 'all'.
combmode : {'t2s', 'ste'}, optional
Combination scheme for TEs: 't2s' (Posse 1999, default), 'ste' (Poser).
combmode : {'t2s', 'paid'}, optional
Combination scheme for TEs: 't2s' (Posse 1999, default), 'paid' (Poser).
label : :obj:`str` or :obj:`None`, optional
Label for output directory. Default is None.

Expand Down Expand Up @@ -162,9 +162,9 @@ def t2smap_workflow(data, tes, mask=None, fitmode='all', combmode='t2s',
LGR.debug('Resulting data shape: {}'.format(catd.shape))

try:
ref_label = os.path.basename(ref_img).split('.')[0]
ref_label = op.basename(ref_img).split('.')[0]
except (TypeError, AttributeError):
ref_label = os.path.basename(str(data[0])).split('.')[0]
ref_label = op.basename(str(data[0])).split('.')[0]

if label is not None:
out_dir = 'TED.{0}.{1}'.format(ref_label, label)
Expand Down
19 changes: 10 additions & 9 deletions tedana/workflows/tedana.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -90,7 +90,7 @@ def _get_parser():
'accepted components'),
default=None)
optional.add_argument('--sourceTEs',
dest='ste',
dest='source_tes',
type=str,
help=('Source TEs for models. E.g., 0 for all, '
'-1 for opt. com., and 1,2 for just TEs 1 and '
Expand All @@ -99,9 +99,9 @@ def _get_parser():
optional.add_argument('--combmode',
dest='combmode',
action='store',
choices=['t2s', 'ste'],
choices=['t2s'],
help=('Combination scheme for TEs: '
't2s (Posse 1999, default), ste (Poser)'),
't2s (Posse 1999, default)'),
default='t2s')
optional.add_argument('--verbose',
dest='verbose',
Expand Down Expand Up @@ -184,7 +184,7 @@ def _get_parser():

def tedana_workflow(data, tes, mask=None, mixm=None, ctab=None, manacc=None,
tedort=False, gscontrol=None, tedpca='mle',
ste=-1, combmode='t2s', verbose=False, stabilize=False,
source_tes=-1, combmode='t2s', verbose=False, stabilize=False,
out_dir='.', fixed_seed=42, maxit=500, maxrestart=10,
debug=False, quiet=False, png=False, png_cmap='coolwarm'):
"""
Expand Down Expand Up @@ -219,11 +219,11 @@ def tedana_workflow(data, tes, mask=None, mixm=None, ctab=None, manacc=None,
is None.
tedpca : {'mle', 'kundu', 'kundu-stabilize'}, optional
Method with which to select components in TEDPCA. Default is 'mle'.
ste : :obj:`int`, optional
source_tes : :obj:`int`, optional
Source TEs for models. 0 for all, -1 for optimal combination.
Default is -1.
combmode : {'t2s', 'ste'}, optional
Combination scheme for TEs: 't2s' (Posse 1999, default), 'ste' (Poser).
combmode : {'t2s'}, optional
Combination scheme for TEs: 't2s' (Posse 1999, default).
verbose : :obj:`bool`, optional
Generate intermediate and additional files. Default is False.
png : obj:'bool', optional
Expand Down Expand Up @@ -352,15 +352,16 @@ def tedana_workflow(data, tes, mask=None, mixm=None, ctab=None, manacc=None,
# Identify and remove thermal noise from data
n_components, dd = decomposition.tedpca(catd, data_oc, combmode, mask,
t2s, t2sG, ref_img,
tes=tes, method=tedpca, ste=ste,
tes=tes, method=tedpca,
source_tes=source_tes,
kdaw=10., rdaw=1.,
verbose=verbose)
mmix_orig = decomposition.tedica(n_components, dd, fixed_seed,
maxit, maxrestart)

if verbose:
np.savetxt(op.join(out_dir, '__meica_mix.1D'), mmix_orig)
if ste == -1:
if source_tes == -1:
io.filewrite(utils.unmask(dd, mask),
op.join(out_dir, 'ts_OC_whitened.nii'), ref_img)

Expand Down