Acquisition speed and spatial resolution in MR imaging are limited by the signal-to-noise ratio (SNR). Images acquired with modern multi-channel receiver systems inevitably suffer from a spatially varying noise distribution. For conventional noise-suppressing image filters working on single-slice magnitude data in integer format, and not having additional knowledge, it is difficult to address this. Filter parameters may have to be tuned manually for different acquisition and reconstruction parameters and SNR levels.

Recently proposed algorithms for iterative reconstruction and compressed sensing contain noise-suppressing steps, but they typically expand – at every iteration – the data back into multi-channel k-space, often non-Cartesian, which increases computation time. Besides, they often require pseudo-random sampling or multiple time points to play out their strength, which limits their application in routine practice.

Supplementary information about the noise in MR data can be acquired, but is typically only internal to the image reconstruction software. There are techniques for predicting^1,2 or measuring³ the k-space and/or image-space noise distribution. This information is crucial for high-quality noise suppression.

In the present work, an iterative denoising algorithm for routine Cartesian MRI is presented, which uses all available information of the acquisition and reconstruction process, does not require manual parameter tuning, and can ultimately be computed in reasonable time.

Routine imaging was performed on 1.5T and 3T clinical whole-body scanners (MAGNETOM Aera & Skyra, Siemens Healthcare, Erlangen, Germany), with a focus on 3D sequences, e.g. fat-saturated or 2pt-Dixon T1-weighted gradient-echo (VIBE) and T2-weighted variable-flip-angle spin-echo (SPACE) imaging of the abdomen, and magnetization-prepared 2-echo gradient-echo (MP2RAGE) imaging of the head. The acquisitions included parallel imaging acceleration with CAIPIRINHA⁴. Acquisition parameters were deliberately pushed for increased resolution and/or decreased acquisition and breath-hold times, leading to borderline SNR.

Complex-valued image data in float precision were exported from a prototypical modification of the image reconstruction pipeline of the scanner, after parallel imaging reconstruction and coil combination, but prior to interpolation and magnitude operation. In addition to the image data, supplementary information was exported, including noise calibration information, image normalization, and k-space filtering. Spatial noise enhancement from the parallel imaging reconstruction, as described by the g-factor, was calculated directly from the reconstruction coefficients². Figure 1 shows exemplary supplementary data.

Iterative denoising was implemented offline in Matlab (MathWorks, Natick, MA, USA). It consists of a bank of orthogonal wavelet transforms (Daubechies 1,2,3,4 with 3 levels), and at each iteration, the image is denoised by thresholding in all different wavelets with a garrote function⁵. The next iterate is computed as a weighted sum of a) all those thresholdings, b) the previous iterate, and c) the original image, see figure 2. The thresholds are spatially adapted to the local noise level, and the recombination weights are iteratively updated to minimize Stein's Unbiased Risk Estimator (SURE)⁶; which optimizes the performance for any image contrast. As a final step, a medium amount of edge enhancement was applied to the processed images to improve perceptual image quality. This was implemented as a prototype modification of a processing module already existing in the scanner image reconstruction system.

After denoising using a constant internal parameterization, the complex-valued image data were re-imported to the scanner reconstruction pipeline, and magnitude images were calculated in the usual manner. As an initial reference, the standard product image filter implementation was used for comparison.

Figures 3-5 show image examples from the abdomen and the head. The local effect of the iterative denoising varies with the local noise level as appropriate. No typical artifacts known from conventional noise suppression filters, like locally insufficient denoising, loss of fine structures, or “oil-painting” appearance, were observed in the initial data collection. No manual tuning of filter parameters depending on image SNR or contrast was needed.

Computation times in these examples for one 3D volume were approximately 7-15 minutes on a standard CPU architecture.

We have presented a new approach to denoising of MRI data, which appears to give better results than conventional noise suppression filters by exploiting supplementary information about the acquisition and reconstruction process and the noise, as well as by working on 3D complex-valued data.

Alternative to wavelet thresholding, the inner-core regularizing operations of the iterations could be chosen differently, e.g. as total variation, nonlocal means, or BM4D.

The results are promising, but need a more systematic investigation, including quantitative SNR, clinical evaluation and further comparisons with reference image filters. Current calculation times are still too long, but are based on a non-optimized implementation. A prototype inline implementation is currently being developed to run directly on the scanner.