dwidenoise: Enhancements #3029

Lestropie · 2024-11-04T04:52:20Z

Doing some work on #3023 (specifically #2742) for a Hackathon.

Draft PR; scope will change over the next couple of days.

Better separation of code responsible for fetching a batch of input data within a sliding spatial window from the code responsible for the denoising of the image data.

New default behaviour is to use an expanding spherical kernel with number of voxels at least 1.1 times the number of volumes. For voxels near the edge of the image FoV, the radius of the kernel will increase until the requisite number of voxels is obtained. Note that execution speed of this implementation seems to be reduced, even when using the cuboid kernel; this may be due to use of Eigen Blocks to denoise voxels with kernels smaller than the maximum processed.

- Added ability to define a spherical kernel of a fixed radius. This will result in voxels near the edges of the image FoV having fewer voxels within the PCA kernel. - Added option -voxels, to generate a spatial map of the number of voxels utilised in the PCA decomposition at each voxel location.

- Add option -max_dist, which exports a map encoding for each voel the maximal distance between that voxel and a voxel included in the PCA decomposition. - Refactor KernelData to keep track of the distance between each voxel within the kernel and the voxel being processed. - Bug fix to spherical kernels; kernel was erroneously being applied from the outside inwards rather than from the origin outwards. - Change type used for data encoding the spherical kernel, which should improve computational performance.

Change from recollection-from-memory ratio of 1.1 to that reported in the corresponding publication of 1.0/0.85.

- Default behaviour is now to use optimal shrinkage based on minimisation of the Frobenius norm. - Prior behaviour can be accessed using "-filter truncate". Closes #3022.

Closes #2591.

Default behaviour is to use a Gaussian kernel (as used in Cordero-Grande et al. 2019) with FWHM equal to twice the voxel spacing. Closes #3024.

An estimated PCA rank that is the weighted average of the patches that contributed to the rank estimation via overcomplete local PCA reconstruction is likely a superior measure of signal rank for downstream use.

Closes #3035.

Lestropie · 2024-11-17T11:53:32Z

@Arshiyasan: The code base has changed a whole lot since you last looked at it. If you want to contribute code for dwidenoise taking an input noise map, it'll need to be refactored, achieving it with a new class derived from Estimator::Base. If you still want to try it, I'll hold out, otherwise LMK and I'll get it done, I've kind of been scratching an itch here.

Closes #3034.

- Add -estimator med, which provides an estimate of the noise level based on the median eigenvalue as per Gavish and Donoho 2014. - Perform linear demodulation of complex input data prior to PCA; and in the case of dwidenoise, re-introduce this phase ramp after denoising has completed. A linear phase ramp is determined per k-space group (2D slice / 3D volume for 3D acquisition) based on the maximum of the interpolated k-space data. - dwidenoise can now filter the eigenspectrum based on noise level truncation (3.0.x behaviour), optimal shrinkage (as implemented previously), or additionally now based on optimal truncation as per Gavish and Donoho 2014. - mrfilter: New mode "demodulate", which provides a direct interface to the same phase demodulation functionality as is used in dwi*noise.

- Estimators now explicitly yield an estimate of the upper end of the noise distribution in addition to the estimated noise level. - Fix to -estimator med (had gross incorrect scaling). - Fix centering of searchlight for even subsample factors. - Fix alignment of voxel grids between input image and output images that are defined with respect to subsampling. - Fix creation of output data where estimator fails completely. - Revert optimal shrinkage to use the upper bound of the MP distribution as the unity reference rather than the noise level. - Fix optimal thresholding failing to construct the local dictionary mapping from aspect ratio to lambda_star parameter to reduce redundant calculations.

Functionality is incompatible with latest enhancements to / capabilities of denoising.

Option gives the ability to compensate for pre-estimated non-fixed noise level within each spatial patch.

File missing from bc33d4d.

- dwi2noise and dwidenoise: Ability to estimate and remove non-linear phase across images prior to denoising. In the case of dwidenoise this phase is re-introduced into the data after PCA denoising. This is now the default over linear phase ramp removal. - mrfilter: New filter "kspace"; this currently only supports a Tukey window, therefore demonstrating the operation of the demodulation that is in use for DWI denoising, but can be used as a standalone operation, and is intended to be extensible to incorporate other k-space filter designs.

Results in the output being reconstructed from a fixed number of eigenvectors, without making any attempt at estimating the signal rank / noise level from the data. Closes 3046.

- Fix erroneous spatial localisation of voxels within patch for reading noise level estimate values. - New option -noise_cov, for investigating the variance in within-patch noise level.

Defining this operation as a variance-stabilising transformation is more accurate, and also flags the fact that that future augmentation accounting for the Rician bias in magnitude data will fall under the same mechanism.

Erroneously omitted from bba4fe1.

By default, if the input image series contains a gradient table that appears to be arranged into shells, then the demeaning will be computed per shell per voxel; otherwise a single mean will be computed per voxel. Demeaning can also be disabled using the -demean option. The code responsible for the demeaning has been integrated with the code for phase demodulation and the variance-stabilising transform as "preconditioning". The presence of demeaning does not affect the PCA or the noise level estimation in any way; it only influences the _reported_ signal rank. Supersedes #2363.

dwidenoise: Modularise kernel

4424cb4

Better separation of code responsible for fetching a batch of input data within a sliding spatial window from the code responsible for the denoising of the image data.

Lestropie added the enhancement label Nov 4, 2024

Lestropie self-assigned this Nov 4, 2024