Skip to content

Commit

Permalink
[ENH] Use different masking thresholds for denoising and classificati…
Browse files Browse the repository at this point in the history 
…on (#736)

* Use two different adaptive masks.

masksum_denoise: Threshold of 1, used for denoising and optimal combination.
masksum_clf: Threshold of 3, used for decomposition and component selection.

* Fix long line.

* Derive classifier masks directly from denoising ones.

* Use denoising mask for component figures.

* Swap full and limited T2*/S0 maps.

Now that the full maps are used throughout the pipeline, they are the primary T2*/S0 outputs. The limited maps are not used for anything, and are only written out if the workflow is run in verbose mode.

* Apply the same logic to t2smap.

* Fix tests.

* Fix documentation about full vs. limited maps.

* Fix bad merge.

* Log updated masking procedure.

* Fix!
  • Loading branch information
@tsalo
tsalo committed Aug 13, 2021
1 parent 7e1cf2c commit 3147e01
Show file tree
Hide file tree
Showing 9 changed files with 113 additions and 88 deletions.
5 changes: 1 addition & 4 deletions docs/approach.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -108,12 +108,9 @@ this voxel), so the line is fit to all available data.
``tedana`` actually performs and uses two sets of :math:`T_{2}^*`/:math:`S_0` model fits.
In one case, ``tedana`` estimates :math:`T_{2}^*` and :math:`S_0` for voxels with good signal in at
least two echoes.
The resulting "limited" :math:`T_{2}^*` and :math:`S_0` maps are used throughout
most of the pipeline.
In the other case, ``tedana`` estimates :math:`T_{2}^*` and :math:`S_0` for voxels
with good data in only one echo as well, but uses the first two echoes for those voxels.
The resulting "full" :math:`T_{2}^*` and :math:`S_0` maps are used to generate the
optimally combined data.
The resulting "full" :math:`T_{2}^*` and :math:`S_0` maps are used throughout the rest of the pipeline.

.. image:: /_static/a05_loglinear_regression.png

Expand Down
36 changes: 20 additions & 16 deletions docs/outputs.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,19 +13,19 @@ Outputs of the tedana workflow
Filename Content
================================================ =====================================================
dataset_description.json Top-level metadata for the workflow.
T2starmap.nii.gz Limited estimated T2* 3D map.
T2starmap.nii.gz Full estimated T2* 3D map.
Values are in seconds.
The difference between the limited and full maps
is that, for voxels affected by dropout where
only one echo contains good data, the full map
uses the single echo's value while the limited
map has a NaN.
S0map.nii.gz Limited S0 3D map.
only one echo contains good data, the full map uses
the T2* estimate from the first two echoes, while the
limited map has a NaN.
S0map.nii.gz Full S0 3D map.
The difference between the limited and full maps
is that, for voxels affected by dropout where
only one echo contains good data, the full map
uses the single echo's value while the limited
map has a NaN.
only one echo contains good data, the full map uses
the S0 estimate from the first two echoes, while the
limited map has a NaN.
desc-optcom_bold.nii.gz Optimally combined time series.
desc-optcomDenoised_bold.nii.gz Denoised optimally combined time series. Recommended
dataset for analysis.
Expand Down Expand Up @@ -79,15 +79,19 @@ If ``verbose`` is set to True:
============================================================== =====================================================
Filename Content
============================================================== =====================================================
desc-full_T2starmap.nii.gz Full T2* map/time series.
desc-limited_T2starmap.nii.gz Limited T2* map/time series.
Values are in seconds.
The difference between the limited and full maps is
that, for voxels affected by dropout where only one
echo contains good data, the full map uses the
single echo's value while the limited map has a NaN.
Only used for optimal combination.
desc-full_S0map.nii.gz Full S0 map/time series. Only used for optimal
combination.
The difference between the limited and full maps
is that, for voxels affected by dropout where
only one echo contains good data, the full map uses
the S0 estimate from the first two echoes, while the
limited map has a NaN.
desc-limited_S0map.nii.gz Limited S0 map/time series.
The difference between the limited and full maps
is that, for voxels affected by dropout where
only one echo contains good data, the full map uses
the S0 estimate from the first two echoes, while the
limited map has a NaN.
echo-[echo]_desc-[PCA|ICA]_components.nii.gz Echo-wise PCA/ICA component weight maps.
echo-[echo]_desc-[PCA|ICA]R2ModelPredictions_components.nii.gz Component- and voxel-wise R2-model predictions,
separated by echo.
Expand Down
16 changes: 8 additions & 8 deletions tedana/resources/config/outputs.json
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,11 +4,11 @@
"bidsv1.5.0": "desc-adaptiveGoodSignal_mask"
},
"t2star img": {
"orig": "t2sv",
"orig": "t2svG",
"bidsv1.5.0": "T2starmap"
},
"s0 img": {
"orig": "s0v",
"orig": "s0vG",
"bidsv1.5.0": "S0map"
},
"combined img": {
Expand Down Expand Up @@ -47,13 +47,13 @@
"orig": "lowk_ts_OC",
"bidsv1.5.0": "desc-optcomRejected_bold"
},
"full t2star img": {
"orig": "t2svG",
"bidsv1.5.0": "desc-full_T2starmap"
"limited t2star img": {
"orig": "t2sv",
"bidsv1.5.0": "desc-limited_T2starmap"
},
"full s0 img": {
"orig": "s0vG",
"bidsv1.5.0": "desc-full_S0map"
"limited s0 img": {
"orig": "s0v",
"bidsv1.5.0": "desc-limited_S0map"
},
"whitened img": {
"orig": "ts_OC_whitened",
Expand Down
4 changes: 2 additions & 2 deletions tedana/tests/data/fiu_four_echo_outputs.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,10 +20,10 @@ desc-PCA_mixing.tsv
desc-PCA_stat-z_components.nii.gz
desc-T1likeEffect_min.nii.gz
desc-adaptiveGoodSignal_mask.nii.gz
desc-full_S0map.nii.gz
desc-full_T2starmap.nii.gz
desc-globalSignal_map.nii.gz
desc-globalSignal_timeseries.tsv
desc-limited_S0map.nii.gz
desc-limited_T2starmap.nii.gz
desc-optcomAcceptedMIRDenoised_bold.nii.gz
desc-optcomAccepted_bold.nii.gz
desc-optcomDenoised_bold.nii.gz
Expand Down
4 changes: 2 additions & 2 deletions tedana/tests/data/nih_five_echo_outputs_t2smap.txt
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
dataset_description.json
desc-full_S0map.nii.gz
desc-full_T2starmap.nii.gz
desc-limited_S0map.nii.gz
desc-limited_T2starmap.nii.gz
desc-optcom_bold.nii.gz
S0map.nii.gz
T2starmap.nii.gz
4 changes: 2 additions & 2 deletions tedana/tests/data/nih_five_echo_outputs_verbose.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,8 +18,8 @@ desc-PCA_metrics.tsv
desc-PCA_mixing.tsv
desc-PCA_stat-z_components.nii.gz
desc-adaptiveGoodSignal_mask.nii.gz
desc-full_S0map.nii.gz
desc-full_T2starmap.nii.gz
desc-limited_S0map.nii.gz
desc-limited_T2starmap.nii.gz
desc-optcomAccepted_bold.nii.gz
desc-optcomDenoised_bold.nii.gz
desc-optcomPCAReduced_bold.nii.gz
Expand Down
20 changes: 10 additions & 10 deletions tedana/tests/test_t2smap.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,9 +31,9 @@ def test_basic_t2smap1(self):
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_T2starmap.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_T2starmap.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_S0map.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-optcom_bold.nii.gz'))
assert len(img.shape) == 4
Expand All @@ -57,9 +57,9 @@ def test_basic_t2smap2(self):
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'S0map.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-full_T2starmap.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_T2starmap.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-full_S0map.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_S0map.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-optcom_bold.nii.gz'))
assert len(img.shape) == 4
Expand All @@ -83,9 +83,9 @@ def test_basic_t2smap3(self):
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_T2starmap.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_T2starmap.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_S0map.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-optcom_bold.nii.gz'))
assert len(img.shape) == 4
Expand All @@ -109,9 +109,9 @@ def test_basic_t2smap4(self):
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'S0map.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-full_T2starmap.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_T2starmap.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-full_S0map.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_S0map.nii.gz'))
assert len(img.shape) == 4
img = nib.load(op.join(out_dir, 'desc-optcom_bold.nii.gz'))
assert len(img.shape) == 4
Expand All @@ -135,9 +135,9 @@ def test_t2smap_cli(self):
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_T2starmap.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_T2starmap.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-full_S0map.nii.gz'))
img = nib.load(op.join(out_dir, 'desc-limited_S0map.nii.gz'))
assert len(img.shape) == 3
img = nib.load(op.join(out_dir, 'desc-optcom_bold.nii.gz'))
assert len(img.shape) == 4
Expand Down
46 changes: 26 additions & 20 deletions tedana/workflows/t2smap.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -171,21 +171,27 @@ def t2smap_workflow(data, tes, out_dir='.', mask=None,
-----
This workflow writes out several files, which are described below:
========================== =================================================
============================= =================================================
Filename Content
========================== =================================================
T2starmap.nii.gz Limited estimated T2* 3D map or 4D timeseries.
============================= =================================================
T2starmap.nii.gz Estimated T2* 3D map or 4D timeseries.
Will be a 3D map if ``fitmode`` is 'all' and a
4D timeseries if it is 'ts'.
S0map.nii.gz Limited S0 3D map or 4D timeseries.
desc-full_T2starmap.nii.gz Full T2* map/timeseries. The difference between
S0map.nii.gz S0 3D map or 4D timeseries.
desc-limited_T2starmap.nii.gz Limited T2* map/timeseries. The difference between
the limited and full maps is that, for voxels
affected by dropout where only one echo contains
good data, the full map uses the single echo's
value while the limited map has a NaN.
desc-full_S0map.nii.gz Full S0 map/timeseries.
good data, the full map uses the T2* estimate
from the first two echos, while the limited map
will have a NaN.
desc-limited_S0map.nii.gz Limited S0 map/timeseries. The difference between
the limited and full maps is that, for voxels
affected by dropout where only one echo contains
good data, the full map uses the S0 estimate
from the first two echos, while the limited map
will have a NaN.
desc-optcom_bold.nii.gz Optimally combined timeseries.
========================== =================================================
============================= =================================================
"""
out_dir = op.abspath(out_dir)
if not op.isdir(out_dir):
Expand Down Expand Up @@ -234,36 +240,36 @@ def t2smap_workflow(data, tes, out_dir='.', mask=None,

# set a hard cap for the T2* map/timeseries
# anything that is 10x higher than the 99.5 %ile will be reset to 99.5 %ile
cap_t2s = stats.scoreatpercentile(t2s_limited.flatten(), 99.5,
cap_t2s = stats.scoreatpercentile(t2s_full.flatten(), 99.5,
interpolation_method='lower')
cap_t2s_sec = utils.millisec2sec(cap_t2s * 10.)
LGR.debug('Setting cap on T2* map at {:.5f}s'.format(cap_t2s_sec))
t2s_limited[t2s_limited > cap_t2s * 10] = cap_t2s
t2s_full[t2s_full > cap_t2s * 10] = cap_t2s

LGR.info('Computing optimal combination')
# optimally combine data
OCcatd = combine.make_optcom(catd, tes, masksum, t2s=t2s_full,
combmode=combmode)

# clean up numerical errors
for arr in (OCcatd, s0_limited, t2s_limited):
for arr in (OCcatd, s0_full, t2s_full):
np.nan_to_num(arr, copy=False)

s0_limited[s0_limited < 0] = 0
t2s_limited[t2s_limited < 0] = 0
s0_full[s0_full < 0] = 0
t2s_full[t2s_full < 0] = 0

io_generator.save_file(
utils.millisec2sec(t2s_limited),
utils.millisec2sec(t2s_full),
't2star img',
)
io_generator.save_file(s0_limited, 's0 img')
io_generator.save_file(s0_full, 's0 img')
io_generator.save_file(
utils.millisec2sec(t2s_full),
'full t2star img',
utils.millisec2sec(t2s_limited),
'limited t2star img',
)
io_generator.save_file(
s0_full,
'full s0 img',
s0_limited,
'limited s0 img',
)
io_generator.save_file(OCcatd, 'combined img')

Expand Down
Loading

0 comments on commit 3147e01

Please sign in to comment.