Display content:


In [1]:
# PYTHON CODE
# Module imports

import matplotlib.pyplot as plt
import plotly.plotly as py
import plotly.graph_objs as go
import numpy as np
from plotly import __version__
from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot
config={'showLink': False, 'displayModeBar': False}

init_notebook_mode(connected=True)

from IPython.core.display import display, HTML

New B1 filtering module

Most studies that use B1 maps do not use the raw maps; instead, they are typically filtered during a pre-processing step to reduce effects like noise and minor artifact, as the B1 distribution in tissue is expected to be a smoothly varying function. However, there appears to be no clear community consensus on which filtering techniques, nor which settings, to use. As mentioned in one of my recent publications (Boudreau et al. (2017), B1 mapping for bias‐correction in quantitative T1 imaging of the brain at 3T using standard pulse sequences), there is a wide range of filtering techniques and setting values reported in the literature, and none refer to publications to explain why those were chosen. It's very likely that different techniques will behave better in certain situations over others; for example, it's possible that one technique works best at filtering large anatomical contrast. For example, some techniques may be more sensitive to small local artifacts with large values, and some filtering implementations may handle the interface at the edge of the brain better than others (where there is lack of data outside the brain).

Beginning with version 2.3.0 of qMRLab, a filtering module has been added, and our double angle B1 mapping module has been extended to include this feature in the options tab. This new feature was implemented by Ilana Leppert, a researcher at the McConnell Brain Imaging Center at the Montreal Neurological Institute, and one of the core qMRLab developers. This qMRLab module also includes contributions by Tim Wu, a graduate student at the Montreal Neurological Institute, who implemented the 3D polynomial filtering algorithm.

Figure 1. Screenshot of filtering module options in qMRLab.

Screenshot of filtering module options in qMRLab. This new module contains a broad range of options that are customizable for the users' needs. Four filtering methods are available: Gaussian kernel, median kernel, polynomial fitting, and spline fitting. Filters can be applied in 2D or 3D mode, meaning either they are applied for single slices (in-plane) or across all three spatial directions. The size of the filters can also be customized. For Gaussian and median filters, the size (in terms of voxels) along each direction can be set. For polynomial and spline filters, the filter order can be set.

Here's an example B1 map dataset compared against the filtered output using each of the four filtering modes.

In [2]:
%% MATLAB/OCTAVE CODE
% Adds qMRLab to the path of the environment

cd ../qMRLab
startup
In [3]:
%% MATLAB/OCTAVE CODE
% Download brain MRI data

cmd = ['curl -L -o filter_map.zip https://osf.io/d8p4h/download/'];
[STATUS,MESSAGE] = unix(cmd);
unzip('filter_map.zip');
In [4]:
%% MATLAB/OCTAVE CODE
%

clear all

Model = filter_map; 

% Format data structure so that they may be fit by the model
data = struct();
data.Raw=double(load_nii_data('Raw.nii.gz'));
data.Mask=double(load_nii_data('Mask.nii.gz'));

Model.options.Smoothingfilter_Type = 'gaussian';
filterSize = 5;
Model.options.Smoothingfilter_sizex = filterSize;
Model.options.Smoothingfilter_sizey = filterSize;
Model.options.Smoothingfilter_sizez = filterSize;
FitResults_gaussian = FitData(data,Model,0);

Model.options.Smoothingfilter_Type = 'median';
filterSize = 5;
Model.options.Smoothingfilter_sizex = filterSize;
Model.options.Smoothingfilter_sizey = filterSize;
Model.options.Smoothingfilter_sizez = filterSize;
FitResults_median = FitData(data,Model,0);

Model.options.Smoothingfilter_Type = 'polynomial';
Model.options.Smoothingfilter_order = 6;
FitResults_polynomial = FitData(data,Model,0);

Model.options.Smoothingfilter_Type = 'spline';
Model.options.Smoothingfilter_order = 2;
FitResults_spline = FitData(data,Model,0);

...done
...done
...done
...done
In [5]:
%% MATLAB/OCTAVE CODE
% Code used to re-orient the images to make pretty figures, and to assign variables with the axis lengths.

sliceIndex = 20;

Mask = imrotate(squeeze(data.Mask(:,:,sliceIndex)), -90);

B1map = imrotate(squeeze(data.Raw(:,:,sliceIndex)), -90).*Mask;

B1map_gaussian = imrotate(squeeze(FitResults_gaussian.Filtered(:,:,sliceIndex)), -90).*Mask;
B1map_median = imrotate(squeeze(FitResults_median.Filtered(:,:,sliceIndex)), -90).*Mask;
B1map_polynomial = imrotate(squeeze(FitResults_polynomial.Filtered(:,:,sliceIndex)), -90).*Mask;
B1map_spline = imrotate(squeeze(FitResults_spline.Filtered(:,:,sliceIndex)), -90).*Mask;

xAxis = [0:size(B1map,2)-1];
yAxis = [0:size(B1map,1)-1];
In [6]:
In [7]:

Figure 2. Raw vs Filtered B1 maps using four filters: Gaussian, median, polynomial, and spline. Gaussian filter settings: 3D filter, filter size = 5 voxels (isotropic). Median filter settings: 3D filter, filter size = 5 voxels (isotropic). Polynomial fit settings: 3D filter, filter order = 6. Spline fit settings: 3D filter, filter order = 2.

Here's a comparison of the 2D and 3D options using the Median filtering mode.

In [8]:
%% MATLAB/OCTAVE CODE
%

clear all

Model = filter_map; 

% Format data structure so that they may be fit by the model
data = struct();
data.Raw=double(load_nii_data('Raw.nii.gz'));
data.Mask=double(load_nii_data('Mask.nii.gz'));

Model.options.Smoothingfilter_Type = 'median';
Model.options.Smoothingfilter_Dimension ='2D';
filterSize = 3;
Model.options.Smoothingfilter_sizex = filterSize;
Model.options.Smoothingfilter_sizey = filterSize;
Model.options.Smoothingfilter_sizez = filterSize;
FitResults_2D = FitData(data,Model,0);

Model.options.Smoothingfilter_Type = 'median';
Model.options.Smoothingfilter_Dimension ='3D';
filterSize = 3;
Model.options.Smoothingfilter_sizex = filterSize;
Model.options.Smoothingfilter_sizey = filterSize;
Model.options.Smoothingfilter_sizez = filterSize;
FitResults_3D = FitData(data,Model,0);

...done
...done
In [9]:
%% MATLAB/OCTAVE CODE
% Code used to re-orient the images to make pretty figures, and to assign variables with the axis lengths.

sliceIndex = 20;

Mask = imrotate(squeeze(data.Mask(60,:,:)), -90);

B1map = imrotate(squeeze(data.Raw(60,:,:)), -90).*Mask;

B1map_2D = imrotate(squeeze(FitResults_2D.Filtered(60,:,:)), -90).*Mask;
B1map_3D = imrotate(squeeze(FitResults_3D.Filtered(60,:,:)), -90).*Mask;

xAxis = [0:size(B1map,2)-1];
yAxis = [0:size(B1map,1)-1]*3;
In [10]:
In [11]: