Volumetric T 1 and T 2 magnetic resonance brain toolkit for relaxometry mapping simulation

Introduction: Relaxometry images are an important magnetic resonance imaging (MRI) technique in the clinical routine. Many diagnoses are based on the relaxometry maps to infer abnormal state in the tissue characteristic relaxation constant. In order to study the performance of these image processing approaches, a controlled simulated environment is necessary. However, a simulated relaxometry image tool is still lacking. This study proposes a computational anatomical brain phantom for MRI relaxometry images, which aims to offer an easy and flexible toolkit to test different image processing techniques, applied to MRI relaxometry maps in a controlled simulated environment. Methods: A pipeline of image processing techniques such as brain extraction, image segmentation, normalization to a common space and signal relaxation decay simulation, were applied to a brain structural ICBM brain template, on both T1 and T2 weighted images, in order to simulate a volumetric brain relaxometry phantom. The FMRIB Software Library (FSL) toolkits were used here as the base image processing needed to all the relaxometry reconstruction. Results: All the image processing procedures are performed using automatic algorithms. In addition, different artefact levels can be set from different sources such as Rician noise and radiofrequency inhomogeneity noises. Conclusion: The main goal of this project is to help researchers in their future image processing analysis involving MRI relaxometry images, offering reliable and robust brain relaxometry simulation modelling. Furthermore, the entire pipeline is open-source, which provides a wide collaboration between researchers who may want to improve the software and its functionality.


Introduction
Magnetic resonance imaging (MRI) is crucial for many applications in the clinic routine.Due to its clinical importance, several MRI techniques are being intensely researched for better image acquisitions and processing.One classical MRI technique is the relaxometry, which measures the characteristic energy decay after a radio-frequency pulse, also known as the relaxometry constants (T2 or T1, depending the weighted imaging technique used) (Carneiro et al., 2006;Deoni, 2010;Haacke et al., 1999).This process is widely used for different studies, especially for brain research in clinical applications (Cheng et al., 2012;Ellingson et al., 2012;Hasan et al., 2012;Kosior et al., 2011), which plays an important role in the diagnosis of several brain diseases, e.g. in Parkinson (Barbosa et al., 2015), Alzheimer (House et al., 2006) and Multiple Sclerosis (Burgetova et al., 2010;Ellingson et al., 2012).For this reason, this classical imaging technique still has an intensive investment.
Following the MRI acquisition and hardware improvement, investigation into image processing related to the relaxometry acquisitions also increased in the scientific literature.Signal decay modelling (Lebel and Wilman, 2010), image denoising algorithms (Feng et al., 2014;Senra et al., 2014) and multi-echo spatial acquisition (Kumar et al., 2012) are some examples of the different image acquisition and processing applied to relaxometry images in order to obtain a more precise relaxation estimate.However, a simulated environment is also important for the development of the research, giving a standard and controlled environment to test new image acquisitions and image processing techniques.In contrast, some computational brain simulations for structural and functional MRI images are well known and widely used in many studies (Chau and McIntosh, 2005;Cocosco et al., 1997;Drobnjak et al., 2006;Rykhlevskaia et al., 2008), but, unfortunately, the relaxometry image simulation is still lacking for brain studies.
In this study, a computational tool for brain MRI relaxometry simulation is proposed.The method described here is based on structural brain templates and results in three-dimensional relaxometry simulations, with controlled image parameter settings.The purpose of this software is to offer a simple, flexible and reliable Volume 32, Number 3, p. 301-305, 2016 image processing toolkit to generate brain relaxometry images for research development.

Methods
The image used as ground-truth tissues to build the simulated relaxometry dataset here is based on the ICBM MRI template (MNI152 brain templates), providing both T1 and T2 weighted MRI images (Chau and McIntosh, 2005;Grabner et al., 2006).The ICBM template is well-known in the scientific literature in several brain studies, which provides a standard brain anatomical template from a sample of registered and normalized group of healthy individuals (Chau and McIntosh, 2005).Furthermore, with the ICBM template, the brain could be divided in different tissues, namely the white matter (WM), gray matter (GM) and cerebrospinal fluid (CSF), among others (Cocosco et al., 1997;Collins et al., 1998).For simplicity, this version of relaxometry simulation only uses the WM, GM and CSF tissues in the simulation process.
In summary, the brain tissues classification and further signal decay simulation are based on a sequential image processing pipeline (Figure 1).Firstly, the brain extraction procedure follows with the neck removal option in the FSL-BET toolkit (Jenkinson et al., 2012;Popescu et al., 2012), which seems to be more appropriate for structural MRI brain extraction (Popescu et al., 2012) and results in an image without voxels belonging on non-brain tissue.Secondly, the image segmentation process follows to extract each brain tissue, i.e.WM, GM, and CSF.In the tissue segmentation step, the expectation maximization algorithm, based on the image histogram classification by the Gaussian mixture models (Jenkinson et al., 2012) and a Markov Random Field classifier are used to separate each brain tissue (Woolrich et al., 2009).Finally, the signal decay is simulated using a first order exponential decay, providing the T1 and T2 characteristic tissue signal behaviour for each voxel in the whole brain.The mathematical model (Equation 1) used here is a first order exponential, which is expected to simulate a first approximation of the natural magnetization relaxation seen in the MRI phenomena (Haacke et al., 1999).The T1 and T2 tissue decay constants are set using the normal brain tissue characteristics, reported in the literature (Wansapura et al., 1999), but the user can set other brain tissue relaxation constants according to their own criteria.(1) Where S 0 is the first signal acquired from the first echo, the S(t) is the signal acquired after a time t, T is the tissue characteristic relaxation constant (T1 or T2) and C is an offset for exponential adjustment.The entire image processing pipeline used here uses the FMRIB Software Library (FSL) tools (Jenkinson et al., 2012).
In order to add more complexity in the brain simulations, two different image artefact possibilities are available: the Rician noise background and the radio-frequency (rf) coil inhomogeneity.Both types of artefacts are commonly assumed in anatomical MRI acquisition techniques and have an important impact for the majority of image processing approaches (Cárdenas-Blanco et al., 2008;Haacke et al., 1999).The Rician noise intensity is based on the brightest tissue for each image acquisition, i.e. the WM signal intensity for T1 images and the CSF signal intensity on T2 images, and the rf inhomogeneity is based on low-frequency rf coils templates (Cocosco et al., 1997).Therefore, three different artefact levels for both the Rician noise (0%, 3% and 9%) and rf inhomogeneity artefacts (0%, 20% and 40%) are available.Here, for the image noise estimate, the same procedure was applied as that described in the BrainWeb brain simulation toolkit (Cocosco et al., 1997).
Furthermore, to provide a comparison between the simulation and real relaxometry estimative, a statistical analysis was conducted using real relaxometry data.A set of 10 real relaxometry images from control asymptomatic subjects, with the same imaging acquisition protocol, was used.The imaging parameters were set as following: GRASE sequence, TR = 3000 ms, EPI factor 3, FOV = 230 × 230 mm 2 , in-plane resolution = 0.479 × 0.479 mm 2 , 20 slices, slice thickness = 2 mm, eleven time echos (TE = 24,36,48,60,72,84,96,108,120,132 and 144 ms), and acquisition time = 6 min.A two-tailed statistical t-test (significance level of 95%) was used to show whether the simulation images provide the same T2 estimate as given in the real data.Furthermore, the absolute error (AE) analysis, between the simulated and real T2 estimate, was also calculated in the same dataset.Only a T2 relaxometry was conducted, but it is assumed that these results could reflect a similar performance in a T1 weighted relaxometry case.

Results
After the image preparation, i.e. brain extraction, brain tissue segmentation and noise simulation procedures, the relaxometry reconstruction is conducted.These initial image preparation steps are illustrated in Figure 1, where the image voxel-wise relaxometry signals are modelled using Equation 1.Each brain tissue has a different relaxation constant value and for this reason, the WM, GM, and CSF tissues must be segmented before the signal decay simulation.
It were reconstructed a whole-brain T2 relaxometry map for each real relaxometry subject.The same first order exponential was used here.The mean T2 value, for each brain tissue, found in the real data was used as initial parameter to create the T2 simulation.All segmented brain tissues (GM, WM, and CSF) were used as samples to calculate the mean T2 values.Figure 2A, B show an example for each T2 image, simulated and real, respectively.The noise and rf inhomogeneity intensity for the simulated comparison were set as 3% and 20%, respectively, which gives a reasonable approximation for the artefact found in the real data.In addition to the T2 estimate, it is possible to compare the Absolute Error (AE) given by the absolute difference between the simulated and real relaxometry maps, which infer the accuracy of the simulation procedure.Figure 2C, D illustrate the results for the T2 relaxometry estimate and AE analysis.

Discussion
The brain relaxometry simulation models created here show a reliable brain structural dependence and also offer a better alternative when compared with the geometrical brain phantom provided in other image simulation tools (Hellerbach et al., 2013;Koay et al., 2007;Van De Walle et al., 2000).Brain structural and functional MRI image simulations are already available for different open-source toolkits, but only for limited options of MRI acquisitions (Jenkinson et al., 2012).The relaxometry mapping simulation still needs further developments and the toolkit presented here could be a further step in this necessary development.
Figure 2C shows the T2 values for each subject for GM, WM, and CSF tissues.Firstly, the T2 value variability is expected, due to image artefacts, where for both the real and simulated dataset offered a T2 variation of around 10%.For this reason, it is important to simulate those Rician and rf coil inhomogeneities in order to obtain a realistic signal variability.Secondly, the T2 simulated values (GM, WM, and CSF) also showed a statistical equivalence compared with the real data.Here the statistical test was evaluated, for all brain tissues, and maintained the H 0 hypothesis, i.e. the simulated and real T2 estimate are not statistically different from each other.Finally, the AE analysis reveals that the simulated dataset offers T2 variability similar to what is seen in the real situation, roughly 10% of the mean value, where its signal fluctuation is directly referred with the noise intensity given in the initial simulation parameters.
Additionally, this entire project is freely available for download on the internet (http://acsenrafilho.github.io/relaxophantom) and it has an open-source license, which provides the liberty to use, adapt and share.Suggestions and collaborations of other researchers are also welcome for the improvement of this tool.At this stage, the toolkit only offers a limited range of noise and simulation parameters, which restricts a more complex study approach.However, in the near future, there is a plan for the expansion of new features such as other non-linear exponential decay simulations, magnetic susceptibility field inhomogeneity, partial segmentation to simulate the partial volume effect in near structural boundaries (WM/GM interface and deep gray matter structures) as well as adding other brain tissues such as fat, blood and others (Cocosco et al., 1997).Moreover, a graphic user interface could be implemented to offer a friendly interaction between the user and the software, which is at the moment, only via command terminal lines.
In conclusion, the tool described here is a suitable alternative for recreating a reasonable brain relaxometry acquisition with controlled imaging artefacts, bringing the convenience of a robust neuroscience toolkit for future research and brain image processing studies.

Figure 1 .
Figure 1.Computational pipeline to reconstruct the tridimensional magnetic resonance imaging relaxometry maps.The entire image processing (brain extraction, segmentation and signal decay simulation) is made using FMRIB Software Library (FSL).It is possible to change the noise intensity and the type of relaxometry image, i.e.T1 or T2 mapping.

Figure 2 .
Figure 2. Comparison between the real and simulated T2 relaxometry maps.A middle brain axial slice is given for (A) simulated and (B) real T2 relaxometry, which the first echo is illustrated.The cerebrospinal fluid (CSF), gray (GM) and white matter (WM) are used here with noise and rf non-uniformity set as 3% and 20%, respectively.(C) Mean T2 value is given by each subject in the real (GM-R, WM-R, and CSF-R) and simulated (GM-S, WM-S, and CSF-S) dataset, which showed that both are equivalent (p<0.05).(D) Absolute error (AE) between real and simulated T2 estimative infer the precision of the simulated reconstruction, recovering the expected variability given by real estimative.