ENH: Add `smooth` input to `RegridToZooms` #549

oesteban · 2020-08-12T09:03:43Z

Resampling (up- or down-) typically requires a gussian smoothing to be
run before, so that outliers are not introduced (esp. when
downsampling).

This PR adds this option, relying on scipy.

This is the first element of a series of PRs deconstructing #506 into digestible, independent bits.

pull-assistant · 2020-08-12T09:03:46Z

Best reviewed: commit by commit

Optimal code review plan (2 commits squashed)

ENH: Add smooth input to RegridToZooms

Update niworkflows/interfaces/images.py

fix: address code re... ... enh: eliminate casti...

Squashed 2 commits:

fix: address code review comments

enh: eliminate casting of smoothed data

Update niworkflows/utils/images.py

Powered by Pull Assistant. Last update 15d1aa0 ... f52c7fb. Read the comment docs.

codecov · 2020-08-12T09:15:59Z

Codecov Report

Merging #549 into master will decrease coverage by 4.85%.
The diff coverage is 37.50%.

@@            Coverage Diff             @@
##           master     #549      +/-   ##
==========================================
- Coverage   64.65%   59.80%   -4.86%     
==========================================
  Files          44       43       -1     
  Lines        5396     5386      -10     
  Branches      791      789       -2     
==========================================
- Hits         3489     3221     -268     
- Misses       1747     2028     +281     
+ Partials      160      137      -23

Flag	Coverage Δ
#documentation	`?`
#masks	`?`
#reportlettests	`?`
#travis	`59.80% <37.50%> (-0.06%)`	⬇️
#unittests	`?`

Flags with carried forward coverage won't be shown. Click here to find out more.

Impacted Files	Coverage Δ
niworkflows/utils/images.py	`74.28% <28.57%> (-3.94%)`	⬇️
niworkflows/interfaces/images.py	`66.14% <100.00%> (+0.25%)`	⬆️
niworkflows/anat/ants.py	`12.35% <0.00%> (-58.43%)`	⬇️
niworkflows/func/util.py	`22.35% <0.00%> (-57.65%)`	⬇️
niworkflows/anat/freesurfer.py	`39.13% <0.00%> (-52.18%)`	⬇️
niworkflows/anat/skullstrip.py	`30.00% <0.00%> (-50.00%)`	⬇️
niworkflows/interfaces/itk.py	`26.92% <0.00%> (-12.31%)`	⬇️
niworkflows/interfaces/fixes.py	`41.17% <0.00%> (-11.77%)`	⬇️
niworkflows/interfaces/bids.py	`87.58% <0.00%> (-8.97%)`	⬇️
niworkflows/interfaces/ants.py	`57.81% <0.00%> (-7.82%)`	⬇️
... and 6 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 43adaf0...f52c7fb. Read the comment docs.

effigies · 2020-08-12T13:05:39Z

niworkflows/utils/images.py

+        data = gaussian_filter(in_file.get_fdata(),
+                               2 if smooth is True else smooth).astype(dtype)


Should we default to something besides 2? Scikit-image uses (max downsampling factor - 1) / 2.

https://github.com/scikit-image/scikit-image/blob/1188c33c3defefdc0f672870b5bcc36beaec967f/skimage/transform/_warps.py#L120

mgxd

minor comments

niworkflows/interfaces/images.py

mgxd · 2020-08-12T13:16:15Z

niworkflows/utils/images.py

@@ -195,6 +194,13 @@ def resample_by_spacing(in_file, zooms, order=3, clip=True):
        new_card.dot(np.vstack((new_grid, np.ones((1, new_grid.shape[1])))))
    )

+    if smooth:
+        from scipy.ndimage import gaussian_filter
+        data = gaussian_filter(in_file.get_fdata(),


why the difference of in_file.get_fdata() and in_file.dataobj?

gaussian_filter wants float, and nibabel does it conveniently. Please note the astype in the next line

You'll probably want to round if re-casting to int, otherwise you truncate 2.9 to 2.

Resampling (up- or down-) typically requires a gussian smoothing to be run before, so that outliers are not introduced (esp. when downsampling). This PR adds this option, relying on scipy.

oesteban · 2020-08-13T14:37:07Z

@mgxd @effigies please double-check my changes make sense and address the problems identified. thanks!

effigies · 2020-08-13T14:46:12Z

Scaling factors are based on zooms, not dimensions. And I don't understand why you're using get_data_dtype to determine the type of the interpreted data.

oesteban · 2020-08-13T14:56:18Z

Scaling factors are based on zooms, not dimensions.

Right, thanks for catching this.

And I don't understand why you're using get_data_dtype to determine the type of the interpreted data.

TBH, because I thought that is the correct way of doing it. Please suggest what would be best.

effigies · 2020-08-13T15:04:32Z

get_data_dtype tells you what type things are on disk. If you want to know the type of the interpreted data, you can use np.asanyarray(img.dataobj).dtype or img.dataobj[0, 0, 0].dtype.

What's the goal of manually setting the dtype? Is there a case where you want to clip your data to integers?

oesteban · 2020-08-13T15:30:53Z

What's the goal of manually setting the dtype? Is there a case where you want to clip your data to integers?

The output image should be of the same dtype of the input - so I guess this is correct?

effigies · 2020-08-13T15:43:03Z

What should be correct? I'm having trouble following the conversation. Can you provide a concrete case where setting the dtype might help?

let ``map_coordinates`` do the trick.

oesteban · 2020-08-13T15:58:15Z

I have removed the casting in the latest commit. After thinking about it, there's no need to re-cast the data because map_coordinates will ensure the output data type matches that of the input image. I agree that introducing a rounding before map_coordinates just wasted accuracy for nothing. Is that what you meant?

Please note that we still need to interrogate the data type on disk because the output of this function should not modify the datatype. Is there any automated operation on the dtype by nibabel? (either casting the data array to the dtype stated on the header or vice-versa, changing the nifti header to reflect the datatype of the data array.

niworkflows/utils/images.py

effigies · 2020-08-13T16:42:09Z

Is there any automated operation on the dtype by nibabel? (either casting the data array to the dtype stated on the header or vice-versa, changing the nifti header to reflect the datatype of the data array.

The data will be cast to the dtype set in the header. If scale factors are needed to make it work, they will be set.

The only time you determine the on-disk dtype from the data object is if you don't pass in a header.

Co-authored-by: Chris Markiewicz <[email protected]>

oesteban · 2020-08-13T16:45:32Z

The data will be cast to the dtype set in the header.

Got it

If scale factors are needed to make it work, they will be set.

Does this mean the cast does not round? I'm not sure we don't want to round in this case.

effigies · 2020-08-13T17:48:33Z

Does this mean the cast does not round?

It will round to the precision of your output type. Basically what you get is:

def rescale_data(data, dtype):
    # assume dtype is int
    dinfo = np.iinfo(dtype)
    # Slope matches change in data range, so we get fixed precision
    slope = (data.max() - data.min()) / (dinfo.max - dinfo.min)
    # We want inter such that orig_min = new_min * slope + inter
    inter = data.min() - dinfo.min * slope
    # And generally, orig_data = new_data * slope + inter
    new_data = (data - data.min()) / slope + dinfo.min
    # Round and truncate to precision of the target data type
    return np.round(new_data).astype(dtype), slope, inter

def reconstruct_data(data, slope, inter):
    return data * slope + inter

data = np.random.default_rng().normal(1000, 200, size=(5,5,5))
new_data, slope, inter = rescale_data(data, 'int16')
reconstructed = reconstruct_data(new_data, slope, inter)
np.abs(reconstructed - data).max()

I'm not sure we don't want to round in this case.

Round to what?

oesteban · 2020-08-13T17:56:31Z

Round to what?

if issubclass(input_nii.header.get_data_dtype(), numbers.Integer) and not issubclass(data.dtype, numbers.Integer):
    data = np.round(data).astype(input_nii.header.get_data_dtype())

effigies · 2020-08-13T18:00:02Z

That's done implicitly.

First release in the 1.3.x series. This release includes enhancements and bug-fixes towards the release of the first LTS version of fMRIPrep. PyBIDS has been revised to use more recent versions, a series of ANTs' interfaces have been deemed ready to upstream into Nipype, and several improvements regarding multi-echo EPI are included. With thanks to Basile Pinsard for contributions. * FIX: Patch ``ApplyTransforms`` spec to permit identity in a chain (#554) * FIX: Add dots to extensions in PyBIDS' config file (#548) * FIX: Segmentation plots aligned with cardinal axes (#544) * FIX: Skip T1w file existence check if ``anat_derivatives`` are provided (#545) * FIX: Avoid diverting CIFTI dtype from original BOLD (#532) * ENH: Add ``smooth`` input to ``RegridToZooms`` (#549) * ENH: Enable ``DerivativesDataSink`` to take multiple source files to derive entities (#547) * ENH: Allow ``bold_reference_wf`` to accept multiple EPIs/SBRefs (#408) * ENH: Numerical stability of EPI brain-masks using probabilistic prior (#485) * ENH: Add a pure-Python interface to resample to specific resolutions (#511) * MAINT: Finalize upstreaming of ANTs' interfaces to Nipype (#550) * MAINT: Update to Python +3.6 (#541)

1.3.0rc3 First release in the 1.3.x series. This release includes enhancements and bug-fixes towards the release of the first LTS version of fMRIPrep. PyBIDS has been revised to use more recent versions, a series of ANTs' interfaces have been deemed ready to upstream into Nipype, and several improvements regarding multi-echo EPI are included. With thanks to Basile Pinsard for contributions. * FIX: Patch ``ApplyTransforms`` spec to permit identity in a chain (nipreps#554) * FIX: Add dots to extensions in PyBIDS' config file (nipreps#548) * FIX: Segmentation plots aligned with cardinal axes (nipreps#544) * FIX: Skip T1w file existence check if ``anat_derivatives`` are provided (nipreps#545) * FIX: Avoid diverting CIFTI dtype from original BOLD (nipreps#532) * ENH: Add ``smooth`` input to ``RegridToZooms`` (nipreps#549) * ENH: Enable ``DerivativesDataSink`` to take multiple source files to derive entities (nipreps#547) * ENH: Allow ``bold_reference_wf`` to accept multiple EPIs/SBRefs (nipreps#408) * ENH: Numerical stability of EPI brain-masks using probabilistic prior (nipreps#485) * ENH: Add a pure-Python interface to resample to specific resolutions (nipreps#511) * MAINT: Finalize upstreaming of ANTs' interfaces to Nipype (nipreps#550) * MAINT: Update to Python +3.6 (nipreps#541)

1.3.0 First release in the 1.3.x series. This release includes enhancements and bug-fixes towards the release of the first LTS (*long-term support*) version of *fMRIPrep*. *PyBIDS* has been revised to use more recent versions, a series of ANTs' interfaces have been deemed ready to upstream into *Nipype*, and several improvements regarding multi-echo EPI are included. With thanks to Basile Pinsard for contributions. * FIX: Patch ``ApplyTransforms`` spec to permit identity in a chain (nipreps#554) * FIX: Add dots to extensions in PyBIDS' config file (nipreps#548) * FIX: Segmentation plots aligned with cardinal axes (nipreps#544) * FIX: Skip T1w file existence check if ``anat_derivatives`` are provided (nipreps#545) * FIX: Avoid diverting CIFTI dtype from original BOLD (nipreps#532) * ENH: Add ``smooth`` input to ``RegridToZooms`` (nipreps#549) * ENH: Enable ``DerivativesDataSink`` to take multiple source files to derive entities (nipreps#547) * ENH: Allow ``bold_reference_wf`` to accept multiple EPIs/SBRefs (nipreps#408) * ENH: Numerical stability of EPI brain-masks using probabilistic prior (nipreps#485) * ENH: Add a pure-Python interface to resample to specific resolutions (nipreps#511) * MAINT: Upstream all bug-fixes in the 1.2.9 release * MAINT: Finalize upstreaming of ANTs' interfaces to Nipype (nipreps#550) * MAINT: Update to Python +3.6 (nipreps#541)

oesteban added the enhancement label Aug 12, 2020

oesteban added this to the 1.2.7 milestone Aug 12, 2020

oesteban requested a review from mgxd August 12, 2020 09:03

oesteban modified the milestones: 1.2.7, 1.3.0 Aug 12, 2020

effigies reviewed Aug 12, 2020

View reviewed changes

mgxd approved these changes Aug 12, 2020

View reviewed changes

oesteban added 2 commits August 13, 2020 05:49

ENH: Add smooth input to RegridToZooms

15d1aa0

Resampling (up- or down-) typically requires a gussian smoothing to be run before, so that outliers are not introduced (esp. when downsampling). This PR adds this option, relying on scipy.

Update niworkflows/interfaces/images.py

cf6d6a0

oesteban force-pushed the enh/smoothing-regridtozooms branch from 26b788b to e7fdd57 Compare August 13, 2020 13:01

oesteban requested a review from effigies August 13, 2020 14:36

fix: address code review comments

f803bdd

oesteban force-pushed the enh/smoothing-regridtozooms branch from e7fdd57 to f803bdd Compare August 13, 2020 14:56

enh: eliminate casting of smoothed data

5a2fdd1

let ``map_coordinates`` do the trick.

effigies reviewed Aug 13, 2020

View reviewed changes

niworkflows/utils/images.py Outdated Show resolved Hide resolved

Update niworkflows/utils/images.py

f52c7fb

Co-authored-by: Chris Markiewicz <[email protected]>

oesteban merged commit 587a15b into nipreps:master Aug 13, 2020

oesteban deleted the enh/smoothing-regridtozooms branch August 14, 2020 10:18

		data = gaussian_filter(in_file.get_fdata(),
		2 if smooth is True else smooth).astype(dtype)

ENH: Add smooth input to RegridToZooms #549

ENH: Add smooth input to RegridToZooms #549

Uh oh!

Conversation

oesteban commented Aug 12, 2020

Uh oh!

pull-assistant bot commented Aug 12, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

codecov bot commented Aug 12, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Codecov Report

Uh oh!

effigies Aug 12, 2020

Choose a reason for hiding this comment

Uh oh!

mgxd left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

mgxd Aug 12, 2020

Choose a reason for hiding this comment

Uh oh!

oesteban Aug 12, 2020

Choose a reason for hiding this comment

Uh oh!

effigies Aug 12, 2020

Choose a reason for hiding this comment

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

effigies commented Aug 13, 2020

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

effigies commented Aug 13, 2020

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

effigies commented Aug 13, 2020

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

Uh oh!

effigies commented Aug 13, 2020

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

effigies commented Aug 13, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

oesteban commented Aug 13, 2020

Uh oh!

effigies commented Aug 13, 2020

Uh oh!

Uh oh!

ENH: Add `smooth` input to `RegridToZooms` #549

ENH: Add `smooth` input to `RegridToZooms` #549

pull-assistant bot commented Aug 12, 2020 •

edited

Loading

codecov bot commented Aug 12, 2020 •

edited

Loading

effigies commented Aug 13, 2020 •

edited

Loading