Hi,
I’m a computer science student and I do my internship in a medical imaging laboratory where I work on machine learning methods for predicting brain age from T1 MRI. I learn every day different things about medical imaging, but I’m still very novice in this field.
In a paper, authors talk about “smoothing MRI with a FWHM kernel of 8mm and resampling with spatial resolution of 8mm”. I understand that they apply a gaussian filter and they reduce MRI resolution, but I have two unanswered questions :
- How do we reduce MRI resolution ? Is it correct to do this simply by using mean operator among neighbors ? That what they mean by saying “resampling” ?
- Is it better to apply gaussian kernel before or after the resampling/reduction ?
Thanks in advance.
You may want to look at my spatial processing demos (executables linked at bottom of the page). This also links to a nice Python demo by @satra. Further, you may want to skim the spatial processing lecture from my PowerPoint files.
For resampling an image to a different resolution or a differnet shape, we can choose our interpolation method. Common ones in our field include:
- Box filter, aka nearest neighbor. Typically only for indexed atlases (e.g. consider a voxel in the output image that is a mixture of area 17 and area 19, it is incorrect to refer to this voxel as ‘area 18’).
- Tent filter, aka trilinear interpolation. This does inherently lose high frequency information. For most neuroimaging this is not an issue: we typically remove high frequencies by blurring our data for our statistics (using a Gaussian Kernel where we specify a wide full-width half maximum). Unless a neuroimager explicitly mentions the interpolation used, it is likely a tent filter.
- Higher order filters, which emulate a windowed-Sinc function. These are often referred to as
spline, bicubic or sinc interpolation depending on the algorithm used. These are better at preserving high frequency information, but this comes at the expense of ringing effects.
There are two issues:
- Different interpolation methods trade off looking blurry versus showing ringing artifacts. The most subjectively pleasing solutions tend to mix interpolation with a bit of blurring. This explains the popularity of the Catmull-Rom, Mitchell-Netravali and the Keys Cubic Filter Family.
- You need to be especially careful when reducing the spatial resolution. Specifically, one needs to be aware of aliasing artifacts. See “General Filtered Image Rescaling” by Dale Schumacher from Graphics Gems III for a seminal description of the issue and elegant code.
In terms of both accuracy and speed, it is ideal to compute the blur and interpolation together as an integrated step. In particular if you will be saving your image data as low-precision integer values. This explains why SPM’s DARTEL allows you to specify the blur during the spatial resampling. However, most neuroimaging pipelines process data as sequential stages. Perhaps @satra can respond regarding the preferred order. I think there is a good case for applying the blur first, as it will remove the high frequencies that can induce higher order interpolation to generate ringing artifacts. On the other hand, if you interpolate first the subsequent blur will at least attenuate the ringing artifacts. Perhaps a good solution is to apply the blur on the highest resolution, e.g. blur first if you the interpolation downsamples the data, blur after if you are upsampling.
2 Likes