DOI: 10.1148/rg.282075083
RadioGraphics 2008;28:345-356
© RSNA, 2008
Informatics in Radiology
Automatic and Adaptive Brain Morphometry on MR Images1
Qingmao Hu, PhD,
Guoyu Qian, MSc,
Michael Teistler, PhD, and
Su Huang, MSc
1 From the Shenzhen Institute of Advanced Integration Technology, Chinese Academy of Sciences, Chinese University of Hong Kong, 3A, Nanshan Medical Instruments Park, 1019 Nanhai Ave, Shenzhen 518067, China (Q.H.); and Biomedical Imaging Lab, Singapore Bio-imaging Consortium, Agency for Science, Technology and Research, Singapore (G.Q., M.T., S.H.). Presented as an Informatics exhibit at the 2006 RSNA Annual Meeting. Received April 25, 2007; revision requested June 22 and received July 24; accepted August 14. All authors have no financial relationships to disclose.
Address correspondence to Q.H. (e-mail: qm.hu{at}siat.ac.cn).
 |
Abstract
|
|---|
Automatic segmentation of brain tissue on magnetic resonance images remains a challenge due to variations in brain shape and size, use of different pulse sequences, overlapping signal intensities, and imaging artifacts. An image analysis system that combines robust image processing techniques with anatomic knowledge was developed to meet this challenge. The system is fast, accurate, and robust to the variability of brain anatomy and imaging conditions and is useful for studying the brain in healthy adults, patients with a shrunken brain due to brain atrophy, and children. With this new thresholding method, the range of the proportion of brain tissue can be determined, thereby making good segmentation possible even in the presence of intrasectional inhomogeneity and noise. The system can adaptively adjust the morphologic processing to break the connection between brain and nonbrain tissue while preserving small brain fragments. It can also segment the white matter and gray matter of the two hemispheres separated by the midsagittal plane. The segmentation results can be visualized in either two or three dimensions. The system has been validated against 53 public data sets and qualitatively tested on 47 clinical data sets, yielding a better accuracy than that of the four most popular methods of brain segmentation.
© RSNA, 2008
|
Introduction
|
|---|
The extraction of brain tissue from three-dimensional (3D) magnetic resonance (MR) imaging volumes is critical for many applications, including brain morphometry (1); registration between conventional and functional MR imaging data (2); visualization of activation in the cortex at functional MR imaging (3); visualization and quantification of the shape of the cortex (4); analysis of the spatial distribution of gray matter (5); localization of functional activation at magneto-encephalography and electroencephalography (6); determination of cortical landmarks for Talairach transformation (7); and characterization of neurologic disorders such as multiple sclerosis and stroke (8), Alzheimer disease (9), Parkinson disease (10), and Klinefelter syndrome (1).
Although several automated and semiautomated brain extraction algorithms (BEAs) are available, a well-accepted system has yet to be developed due to the wide variations in brain shape and size; factors that inhere in head MR imaging (noise, gray-scale inhomogeneity, partial volume effects, artifacts); and the proximity of brain tissue to nonbrain tissue (eg, skull, orbits, skin, optic nerves, meninges, sinuses), in terms of both anatomic location and signal intensity (7). In this article, we discuss and illustrate a brain morphometry system that we developed for use with clinical MR imaging data.
|
Brain Extraction Algorithms
|
|---|
Algorithms for brain extraction (also called skull stripping or head peeling) can be characterized in many ways. Like Rehm et al (11), we classify BEAs operating on single T1-weighted MR imaging volumes according to their dominant operations.
Threshold-with-Morphology–based Algorithms
Threshold-with-morphology–based algorithms typically compute or accept as input one or more signal intensity thresholds to help distinguish between gray matter and white matter. Mathematic morphologic operations, connected component analysis, and thresholding are typically used. Erosion and dilation are the two fundamental mathematic morphologic operations; the former is used to break the connection between brain and nonbrain tissue, and the latter is used to recover the pixels that have been eroded away. Connected component analysis is used to group similar tissues that are spatially connected, such as gray matter and white matter. Thresholding is a simple way to distinguish between different types of tissue, since different tissues have different signal intensity ranges, although these ranges may overlap. Examples of threshold-with-morphology–based algorithms are described in references 8 and 12–15.
Edge-based Algorithms
With edge-based algorithms, signal intensity gradients are used to define edges of variable strength. Closed contours are defined, and morphologic operations are applied to the enclosed regions to extract the brain. The method described by Shattuck et al (16) is an example of this approach.
Algorithms That Make Use of A Priori Information
A priori information can take the form of a brain atlas or models of the brain. Examples of algorithms that make use of such information are described by Ashburner and Friston (17) and by Atkins and Mackiewich (18).
Deformable-Surface Algorithms
Algorithms that make use of deformable surfaces create an initial shape or contour, either automatically or with operator input. Forces derived from signal intensities then iteratively drive the shape to either a cortical envelope or a convoluted gray matter surface. The method described by Smith (19) is an example of this approach.
Hybrid Models
Hybrid models combine existing BEAs. The method described by Rehm et al (11) relies on warping to a template, threshold-with-morphology operations, and edge detection. The method described by Rex et al (20) combines the techniques described by Shattuck et al (16) and Smith (19), 3DIntracranial from AFNI software (21), and watershed data obtained with FreeSurfer software (22). The method described by Segonne et al (23) combines watershed algorithms and deformable-surface models.
Thresholds may be determined for BEAs in one of two ways: gaussian curve fitting (8,13,18, 24) and empirical formulas (8,14). However, experience has shown that the mixture-gaussian model with a fixed number of mixtures may be inapplicable due to the complexity of real MR imaging data (25). Moreover, empirical formulas cannot be adapted to the substantial variations in imaging conditions and the wide anatomic variability among brains (14).
Thus, a method for determining thresholds that are robust to artifacts, adaptable to imaging conditions, and applicable in different individuals is highly desirable. As for the morphologic operations, a larger structuring element (SE) is helpful in breaking connections between brain and non-brain tissue, but small brain fragments may be discarded. On the other hand, a smaller SE is helpful in preserving small brain fragments, but undesirable connections may remain unbroken. Thus, an adaptable-size SE is desirable. Existing threshold-with-morphology–based BEAs that make use of a fixed-size SE and global thresholding have serious limitations: (a) they cannot handle data with significant signal intensity inhomogeneity due to the global thresholding; and (b) because of the variability in data, the SE may be too small to break the connection between brain and nonbrain tissue or too large to preserve small brain fragments.
Our goal has been to provide an image analysis system for brain segmentation that could process clinical MR images. Our system is of the threshold-with-morphology type and combines local thresholding, morphologic operations with an adaptable-size SE, and anatomic knowledge. The signal intensity inhomogeneity over the entire 3D data set is divided into intrasectional and intersectional inhomogeneity. The former type is addressed with knowledge-based thresholding, whereas the latter type is addressed with local thresholding (ie, thresholds that are pixel position dependent). The smallest SE that can break the connection between brain and nonbrain tissue is automatically determined for each data set.
In image processing, the picture element is called a pixel in two dimensions (2D) and a voxel in 3D. In brain segmentation of 3D data, sometimes the distinction between 2D and 3D is not clear, so we use the term pixel for both 2D and 3D picture elements for the sake of simplicity.
|
Materials and Methods
|
|---|
The performance of BEAs is highly dependent on the quality of the chosen data; thus, it is desirable to test BEAs with ground truth on publicly available data that have been widely used. In our study, 53 public data sets were used for validation and quantification, including 18 BrainWeb phantom (BWP) T1-weighted data sets (181 x 217 x 181 1-mm isotropic pixels) with 0%, 1%, 3%, 5%, 7%, or 9% noise and 0%, 20%, or 40% inhomogeneity (http://www.bic.mni.mcgill.ca/brainweb); 20 normal T1-weighted data sets (1.5 T, 256 x 256 section dimension, 60–65 coronal sections, 1 x 1 x 3-mm3 pixels) from the Internet Brain Segmentation Repository (IBSR) (http://www.cma.mgh.harvard.edu/ibsr); and 15 clinical T1-weighted data sets (1.5 T, 256 x 256 section dimension, 180 axial sections, 0.86 x 0.86 x 1-mm3 pixels) from Brain Extraction Evaluation Service (BEES) (http://www.neurovia.umn.edu/webservice/BEE_service.html). Forty-seven clinical spoiled gradient-echo recovery data sets with 3-T field strength (Singapore General Hospital, Singapore) or 1.5-T field strength (Nagasaki University School of Medicine, Nagasaki, Japan) (pixel size, 0.86–3 mm) were also used for qualitative evaluation.
The brain (the sum of white matter and gray matter) is segmented by (a) determining signal intensity thresholds for all axial sections, (b) determining the optimum-sized SE for morphologic operations, and (c) refining the study. Our assumptions are that (a) the MR imaging data may have serious inhomogeneity, which demands local thresholding; and (b) the connection between brain and nonbrain tissue is variable, so that the morphologic operation should be adaptable to allow the connection to be broken while preserving small fragments of brain tissue.
After the brain is extracted, white and gray matter are differentiated with use of fuzzy c-means clustering (26). The extracted brain or the gray and white matter in the two hemispheres can then be displayed as 2D sections or as a 3D model with volume rendering (VR).
Determination of Signal Intensity Thresholds
Thresholding consists of (a) finding a middle axial section (the reference section) to determine the threshold 
ref by exploring the proportion of brain tissue within the skull; (b) finding the superior and inferior axial sections and calculating the threshold sup/inf by exploring the proportion of brain tissue within the skull; and (c) determining the thresholds for all axial sections on the basis of these three thresholds by means of linear interpolation.
The image is usually brightest in the middle of the imager; thus, the reference section is taken to be the axial section with the maximum average signal intensity. The threshold ref is determined with supervised range-constrained thresholding (25). Specifically, the region of interest (ROI) includes all the pixels enclosed within the skull (Fig 1b). Within the ROI, the proportion of non-brain tissue varies between 13% and 30%. For all pixels within the ROI, the minimum/maximum signal intensity i1/i2 is determined such that the proportion of pixels with signal intensities within the range i1–i2 is within [13%,30%]. The threshold ref is determined on the basis of maximized between-class variance; in other words, ref divides all pixels within the ROI with signal intensities in the range of i1–i2 into two classes to maximize between-class variance. In this way, the brain anatomy (the range of foreground proportion) and the image processing technique (maximized between-class variance) are combined to yield a threshold that is robust to intrasectional inhomogeneity and noise (7). Figure 1c shows the brain tissue as well as skull pixels with signal intensities greater than ref. The remaining non-brain tissue with a signal intensity similar to that of the brain but in a different anatomic location is removed by means of morphologic operations and connected component analysis.
View larger version (166K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 1a. Reference section of the BrainWeb phantom 1_9_40: the reference image (a), the head mask or ROI (b), and the image as it appears with the pixels whose signal intensity is below the threshold ref set to 0(c).
|
|
View larger version (18K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 1b. Reference section of the BrainWeb phantom 1_9_40: the reference image (a), the head mask or ROI (b), and the image as it appears with the pixels whose signal intensity is below the threshold ref set to 0(c).
|
|
View larger version (131K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 1c. Reference section of the BrainWeb phantom 1_9_40: the reference image (a), the head mask or ROI (b), and the image as it appears with the pixels whose signal intensity is below the threshold ref set to 0(c).
|
|
Finding the superior or inferior axial section requires an approximation of the brain, which we perform by (a) binarizing the data by setting any pixel with a signal intensity greater than ref as foreground (nonblack) and as background otherwise (black); (b) eroding the foreground components with a cuboid SE of 6 mm (equivalent to removing the outermost 3 mm of foreground pixels), the concern being to break the connection rather than preserve small brain fragments; (c) finding the largest foreground components; and (d) dilating the largest foreground component (ie, the foreground component with the largest number of foreground pixels) with the same SE (equivalent to adding 3 mm of pixels to the boundary of the foreground component). The superior or inferior axial section is then found to be the most superior (or inferior) axial section with at least 100 (or 1000) mm2 foreground pixels. The threshold sup/inf can be determined similar to the threshold ref by specifying the proportion of nonbrain pixels as ranging from 20% to 55% for sup and from 28% to 58% for inf. The threshold (z) for axial section z can then be determined with linear interpolation (Fig 2); that is, axial sections between the middle and superior axial sections will have thresholds between ref and sup, whereas those between the middle and inferior axial sections will have thresholds between ref and inf.
View larger version (4K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 2. Graph illustrates that the signal intensity threshold (z) at axial section z is determined by linear interpolation from the thresholds at the superior (left), middle, and inferior (right) axial sections.
|
|
Determination of the Optimum-sized SE for Morphologic Operations
The optimum-sized SE is the smallest SE that could break the connection between brain and nonbrain tissue. For ease of implementation, we use only a cuboid SE. In our implementation, each side of the SE can be 6 mm (SE-6), 4 mm (SE-4), or 2 mm (SE-2). For each SE, the corresponding brain tissue can be derived in four steps: (a) set pixels at section z to foreground if their signal intensities are greater than the threshold (z), and to background otherwise; (b) erode all the foreground components with the SE; (c) find the maximum foreground component (in terms of number of pixels) of all the eroded foreground components; and (d) dilate the maximum foreground component with the same SE. To determine if the foreground output contains substantial nonbrain tissue, the proportion of pixels less than 10 mm from the boundary pixels of the skull (Fig 1b) is calculated. The smallest SE with a calculated proportion of less than 5% is the optimum SE. SE-2 is optimal for preserving small brain fragments (Fig 3), whereas SE-6 is optimal for breaking the connection between brain and nonbrain tissue (Fig 4). Next, the dilated maximum foreground component bin0 of the optimum SE is fine tuned to obtain the final brain segmentation.
View larger version (169K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 3a. Axial section of the BrainWeb phantom 1_9_40: the original image (a), the segmented brain with SE-2 (b), and the segmented brain with SE-4 (c). SE-4 is too large to preserve small brain fragments.
|
|
View larger version (94K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 3b. Axial section of the BrainWeb phantom 1_9_40: the original image (a), the segmented brain with SE-2 (b), and the segmented brain with SE-4 (c). SE-4 is too large to preserve small brain fragments.
|
|
View larger version (87K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 3c. Axial section of the BrainWeb phantom 1_9_40: the original image (a), the segmented brain with SE-2 (b), and the segmented brain with SE-4 (c). SE-4 is too large to preserve small brain fragments.
|
|
View larger version (129K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 4a. Axial section of IBSR data: the original image (a), the segmented brain with SE-4 (b), and the segmented brain with SE-6 (c). SE-4 is too small to break the connection between brain and nonbrain tissue.
|
|
View larger version (100K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 4b. Axial section of IBSR data: the original image (a), the segmented brain with SE-4 (b), and the segmented brain with SE-6 (c). SE-4 is too small to break the connection between brain and nonbrain tissue.
|
|
View larger version (97K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 4c. Axial section of IBSR data: the original image (a), the segmented brain with SE-4 (b), and the segmented brain with SE-6 (c). SE-4 is too small to break the connection between brain and nonbrain tissue.
|
|
Refinement Stage
The refinement stage includes recovering small brain fragments due to morphologic erosion with the optimum SE and removing the sagittal sinus.
If the optimum SE is not SE-2, bin0 will undergo a conditional dilation. For each foreground pixel of bin0, change all its neighboring pixels at axial section z with a signal intensity greater than the threshold (z) from background to foreground to obtain a new binary volume bin1.
The sagittal sinus is situated near the midsagittal plane (27) and has signal intensities slightly lower than those of gray matter. Thus, the sagittal sinus is removed in three steps: (a) find the axial section range for removing the sagittal sinus, starting from the most superior axial section of bin1 (containing at least one foreground pixel) and ending at an axial section 15 mm inferior to the most superior axial section; (b) calculate the sinus threshold for each axial section (average signal intensity plus the standard deviation of all foreground pixels of the axial section in the midsagittal plane); and (c) remove the sagittal sinus pixels if their distance from the midsagittal plane is less than 3 mm and their signal intensities are below the sinus threshold to yield the final brain binf. Figure 5 shows clinical data with most of the sagittal sinus removed.
Differentiation between White Matter and Gray Matter
Differentiation between white matter and gray matter is important in the quantitative analysis of brain images and has sparked much interest among researchers in recent years. Many approaches have been studied and developed. These approaches can be broadly divided into three categories: classification-based, region-based, and contour-based methods (28). Classification-based fuzzy c-means clustering is used because it is one of the unsupervised clustering methods popularly used for tissue segmentation (29–32). Although fuzzy c-means clustering is a reliable means of white matter–gray matter classification, it does not address the problem of signal intensity inhomogeneity artifact. To circumvent this deficiency, an automatic inhomogeneity correction method (33) was used prior to segmentation of white matter and gray matter. This correction method automatically selects reference points for the construction of a low-order polynomial function to correct for image inhomogeneity.
Sequential application of the inhomogeneity correction operation and fuzzy c-means clustering allows the segmentation of white matter and gray matter (Fig 6).
View larger version (57K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 6a. Segmentation of a 3-T clinical data set. (a) Axial section before (left) and after (right) the segmentation of white matter (blue) and gray matter (green). (b) Three-dimensional models of the brain show brain tissue in the left (green) and right (purple) hemispheres.
|
|
View larger version (62K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 6b. Segmentation of a 3-T clinical data set. (a) Axial section before (left) and after (right) the segmentation of white matter (blue) and gray matter (green). (b) Three-dimensional models of the brain show brain tissue in the left (green) and right (purple) hemispheres.
|
|
Rendering of the Extracted Brain, Gray Matter, and White Matter
For the presentation of the segmentation results, a software module was developed for interactive 2D and 3D visualization of segmented volumetric data sets, using and extending the software framework of an in-house application (VolumeExplorer [http://www.virtusmed.info/VE]) (34). The software allows section-based visualization and VR as shown in Figure 7. Each segmented structure can be interactively assigned an arbitrary color intensity and opacity and displayed on a VR image or overlaid on a 2D image. Lower original MR imaging signal intensities are displayed with lower color intensity and opacity in a user-definable way. In addition, the software allows interactive visualization of subregions of a segmented structure by enabling the user to define a point of interest and a threshold for growing a local region (ie, grouping pixels into a region) based on the original MR imaging signal intensities. Both 2D and 3D visualization are based on OpenGL technology (35), making use of the Open Inventor graphics software development kit (36). In particular, VR is performed on the basis of 3D textures (37), thereby allowing interactive rendering through the use of standard OpenGL features that are supported by many graphics cards nowadays. Each pixel of a segmented volume data set comprises 16 bits. Eight bits are used for the original MR imaging signal intensities, and the remaining bits are used for binary segmentation information (eg, differentiation between brain and nonbrain tissue or between white and gray matter). For the VR process, a color intensity and an opacity value are assigned to each combination of original signal intensity bits and segmentation bits with use of a hardware-supported look-up table. The software has been used on a Pentium IV personal computer (Intel, Santa Clara, Calif) with 2.6 GHz, 1024 MB of main memory, and an Nvidia Geforce 7900 graphics card with 256 MB of graphics memory.
View larger version (79K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 7a. VR image (a) and views obtained in different planes with overlays (b) depict segmentation results in 3D and 2D, respectively. Various visual properties can be interactively assigned to each segmented structure.
|
|
View larger version (80K):
[in this window]
[in a new window]
[Download PPT slide]
|
Figure 7b. VR image (a) and views obtained in different planes with overlays (b) depict segmentation results in 3D and 2D, respectively. Various visual properties can be interactively assigned to each segmented structure.
|
|
|
Results
|
|---|
Brain extraction and white matter–gray matter classification were implemented in C++ and have been quantitatively validated (only the brain segmentation) against the 53 public data sets and qualitatively tested against 47 clinical data sets obtained in healthy individuals. It took less than 2 minutes to process any one data set with the Pentium IV personal computer (Intel). Figure 7 shows segmentation results in both 2D and 3D. Figure 8 shows the segmentation results from a 1.5-T data set obtained in an aged individual with a shrunken brain. Figure 9 shows axial sections (and the corresponding postsegmentation images) from a 3-T clinical data set with serious inhomogeneity.
To quantify BEAs, the Jaccard similarity index (Tamimoto coefficient) and the 
index (Dice coefficient) were calculated in the manner frequently reported in the literature (8,11,38). In addition, the segmentation errors (false-positive or false-negative rate) were estimated by dividing the number of incorrectly labeled nonbrain (or brain) pixels by the total number of brain pixels in the ground truth volume (38).
For the 18 BWP, 20 IBSR, and 15 BEES data sets, our method yields average indexes of 0.976 ± 0.003, 0.964 ± 0.013, and 0.972 ± 0.011, respectively; average false-positive rates of 3.61% ± 0.19%, 3.62% ± 2.45%, and 0.63% ± 0.34%, respectively; and average false-negative rates of 1.18% ± 0.55%, 2.85% ± 1.67%, and 4.18% ± 2.01%, respectively.
We compared our BEA (adaptive brain segmentation [ABS]) with up-to-date freely available BEAs: statistical parametric mapping (SPM2) (17), hybrid watershed algorithm (HWA) (23), brain surface extraction (BSE) (16), and brain extraction tool (BET) (19). The average accuracies of these methods against the 53 public data sets are listed in the Table. BSE is interactive in that parameters need to be changed for good segmentation.
|
Discussion
|
|---|
For the 53 public data sets, the proposed algorithm, adaptive brain segmentation, is consistently more accurate than the freely available BEAs (Table). In a study by Lemieux et al (24), the BWP data set with 3% noise and 20% nonuniformity was tested and yielded a index of 0.977, which was almost identical to our results (0.978). The McStrip (11) has been tested against the 15 nonanonymized BEES data sets, yielding an average index of 0.988. According to http://www.neurovia.umn.edu/webservice/BEE_interp.html, the average index for the Mc-Strip against the anonymized data should be around 0.968, which is slightly lower than our results (0.972). The level set–based method (39) yields an average index of 0.96 against the 20 IBSR data sets, which is similar to our results (0.964). However, this method is much more time consuming (5.3 seconds per section, or about 22 minutes for 256 sections). Direct comparison with BEMA (brain extraction meta-algorithm) (20) was not made because we were unable to implement BEMA and could not determine its accuracy against any of the 53 public data sets.
We have not carried out more extensive quantitative validation due to the huge amount of work and difficulty involved in generating high-quality ground truth.
There are several parameters associated with the algorithm that are based on training as well as our experience. In training, to account for possible tilts during imaging, we rotate clinical data (from Japan and Singapore, both 1.5 T and 3 T, not included in the 47 clinical data sets used for qualitative validation) by 15° in the coronal, axial, and sagittal directions. The background proportions for the three thresholds ref, sup, and inf, and the two constants 100 mm2 and 1000 mm2, are based on extensive validation of clinical data. On the basis of our experience and tests, we observed that (a) the SE should not be larger than 6 mm per side (SE-6) to preserve brain fragments, (b) pixels less than 10 mm from the boundary pixels of the skull can be considered possible nonbrain pixels for calculating the non-brain proportion of a brain segmentation, and (c) a calculated proportion of more than 5% of nonbrain tissue is likely due to unexpected data such as some sections being significantly brighter. We also successfully tested our algorithm qualitatively against other clinical data (including data obtained in both young and old individuals, spoiled gradient-echo recovery data, and T1-weighted data). Because these parameters have been tested on data sets obtained at different centers with imagers of various ages (range, 13–65 years), obtained in individuals with different ethnic backgrounds (Asian, Caucasian), and containing various imaging artifacts, we believe that they may be generally applicable.
The main advantage of our system is that it can extract the brain from clinical images in the presence of signal intensity inhomogeneity and break the connection between brain and nonbrain tissue with an optimized morphologic operation to preserve small brain fragments. These functions are realized through (a) thresholding that is based on the variable proportion of background in the head ROI on the superior, middle, and inferior axial sections and is used to determine the local thresholds to account for individual anatomy and imaging conditions; (b) changing the size of the SE such that the proportion of pixels less than 10 mm from the head mask is less than 5%, which makes the SE adaptable to connections between brain and nonbrain tissue at varied connection strengths; and (c) refinement based on anatomic knowledge. Another exclusive feature is the capacity to warn of the need for human intervention due to excessive inclusion (>5%) of nonbrain tissue in the segmented brain. The system also provides editing functionality (including zoom, pan, and addition and removal of points and regions) for cases in which the user is not satisfied with the segmentation results.
It should be stressed that the system is not without limitations. The tested data were all obtained in healthy individuals, and only T1-weighted and spoiled gradient-echo recovery MR imaging volumes have been tested. In addition, we have not tested any data with a pixel size greater than 3 mm.
Actually, testing on "normal" data sets does not necessarily limit the applicability of the system in unhealthy individuals. As pointed out by Lemieux et al (14), the main factors that may limit applicability are (a) the presence of brain tissue with a signal intensity significantly lower than the threshold (z), and (b) the presence of brain structures smaller than a side of the optimum SE. Therefore, if lesions or abnormalities within the brain are hypointense relative to cerebrospinal fluid, they are difficult to identify and segment; in other words, any abnormalities with a signal intensity greater than that of cerebrospinal fluid can be included, and if combined with the anatomic model, these abnormalities could be segmented (eg, in cases of brain hemorrhage or infarction). Because the accuracy is around 3% in a typical clinical setting, a variation in the brain of more than 3% is expected, which could clinically be used for skull peeling and brain morphometry as well as for tracing the morphologic features of the brain after surgery or treatment. Moreover, the warning capability of the system can remind the imaging expert to enhance quality or intervene when the segmented brain contains substantial skull and is ideal for clinical use.
Compared with methods based on more advanced processing techniques such as level set–based methods (39) and hybrid methods (20), the proposed method is faster and simpler to use (parameters are automatically determined and validated) and has comparable accuracy. We believe that the proposed system can be a useful tool for both research in and the clinical practice of brain morphometry.
|
Conclusions
|
|---|
We have developed an adaptive brain segmentation system that has been validated against 53 public data sets to yield a better accuracy ( = 0.970) than the four most popular methods. The system is fast (<2 minutes per data set on a standard personal computer), accurate, and robust to the variability of brain anatomy and imaging conditions. It can be used for brain morphometry in healthy adults as well as in individuals with a shrunken brain due to brain atrophy (eg, from Alzheimer disease) and in children. The proposed ideas of adaptive thresholding based on local thresholds of the superior, middle, and inferior axial sections and adaptive determination of the optimum-sized SE may be applicable in brain segmentation performed with other pulse sequences and in improving other existing methods such as brain surface extraction (16) in which signal intensity thresholding or morphologic processing is used.
|
TAKE-HOME POINTS
|
|---|
Erosion and dilation are the two fundamental mathematic morphologic operations; the former is used to break the connection between brain and nonbrain tissue, and the latter is used to recover the pixels that have been eroded away.
The signal intensity inhomogeneity over the entire 3D data set is divided into intrasectional and intersectional inhomogeneity. The former type is addressed with knowledge-based thresholding, whereas the latter type is addressed with local thresholding.
It is assumed that the connection between brain and nonbrain tissue is variable, so that the morphologic operation should be adaptable to allow the connection to be broken while preserving small fragments of brain tissue.
The optimum-sized SE is the smallest SE that could break the connection between brain and nonbrain tissue.
The main advantage of our system is that it can extract the brain from clinical images in the presence of signal intensity inhomogeneity and break the connection between brain and nonbrain tissue with an optimized morphologic operation to preserve small brain fragments. These functions are realized in part through thresholding that is based on the variable proportion of background in the head ROI on the superior, middle, and inferior axial sections and is used to determine the local thresholds to account for individual anatomy and imaging conditions.
|
Acknowledgments
|
|---|
The authors would like to thank Sauw Ming Liong, Liang Yoong Ho, and Chee Kheong Wong for their work in visualization of the segmented brain.
|
Footnotes
|
|---|
Abbreviations: BEA = brain extraction algorithm, BEES = Brain Extraction Evaluation Service, BWP = BrainWeb phantom, IBSR = Internet Brain Segmentation Repository, ROI = region of interest, SE = structuring element, VR = volume rendering, 2D = two-dimensional, 3D = three-dimensional
|
References
|
|---|
- Shen D, Liu D, Liu H, Clasen L, Giedd J, Davatzikos C. Automated morphometric study of brain variation in XXY males. Neuroimage 2004; 23(2):648–653.[CrossRef][Medline]
- Friston KJ, Ashburner J, Frith CD, Poline JB, Heather JD, Frackowiak RS. Spatial registration and normalization of images. Hum Brain Mapp 1995;3(3):165–189.[CrossRef]
- Desai R, Liebenthal E, Possing ET, Waldron E, Binder JR. Volumetric vs. surface-based alignment for localization of auditory cortex activation. Neuroimage 2005;26(4):1019–1029.[CrossRef][Medline]
- Thompson PM, MacDonald D, Mega MS, Holmes CJ, Evans AC, Toga AW. Detection and mapping of abnormal brain structure with a probabilistic atlas of cortical surfaces. J Comput Assist Tomogr 1997;21(4):567–581.[CrossRef][Medline]
- Woermann FG, Sisodiya SM, Free SL, Duncan JS. Quantitative MRI in patients with idiopathic generalized epilepsy. Brain 1998;121(pt 9):1661–1667.[Abstract/Free Full Text]
- Dale AM, Sereno MI. Improved localization of cortical activity by combining EEG and MEG with MRI cortical surface reconstruction: a linear approach. J Cognit Neurosci 1993;5(2):162–176.[CrossRef]
- Hu QM, Qian GY, Nowinski WL. Fast, accurate, and automatic extraction of the modified Talairach cortical landmarks from magnetic resonance images. Magn Reson Med 2005;53(4):970–976.[CrossRef][Medline]
- Shan ZY, Yue GH, Liu JZ. Automated histogram-based brain segmentation in T1-weighted three-dimensional magnetic resonance head images. Neuroimage 2002;17(3):1587–1598.[CrossRef][Medline]
- Juottonen K, Laakso MP, Partanen K, Soininen H. Comparative MRI analysis of the entorhinal cortex and hippocampus in diagnosing Alzheimer’s disease. AJNR Am J Neuroradiol 1999; 20(1):139–144.[Abstract/Free Full Text]
- Price S, Paviour D, Scahill R, et al. Voxel-based morphometry detects patterns of atrophy that help differentiate progressive supranuclear palsy and Parkinson’s disease. Neuroimage 2004;23(2):663–669.[CrossRef][Medline]
- Rehm K, Schaper K, Anderson J, Woods R, Stoltzner S, Rottenberg D. Putting our heads together: a consensus approach to brain/non-brain segmentation in T1-weighted MR volumes. Neuroimage 2004;22(3):1262–1270.[CrossRef][Medline]
- Hohne KH, Hanson WA. Interactive 3D segmentation of MRI and CT volumes using morphological operations. J Comput Assist Tomogr 1992; 16(2):285–294.[Medline]
- Brummer ME, Mersereau RM, Eisner RL, Lewine RR. Automatic detection of brain contours in MRI dataset. IEEE Trans Med Imaging 1993; 12(2):153–166.[CrossRef][Medline]
- Lemieux L, Hagemann G, Krakow K, Woermann FG. Fast, accurate, and reproducible automatic segmentation of the brain in T1-weighted volume MRI data. Magn Reson Med 1999;42(1):127–135.[CrossRef][Medline]
- Stokking R, Vincken KL, Viergever MA. Automatic morphology-based brain segmentation (MBRASE) from MRI-T1 data. Neuroimage 2000;12(6):726–738.[CrossRef][Medline]
- Shattuck DW, Sandor-Leahy SR, Schaper KA, Rottenberg DA, Leahy RM. Magnetic resonance image tissue classification using a partial volume model. Neuroimage 2001;13(5):856–876.[Medline]
- Ashburner J, Friston KJ. Voxel-based morphometry: the methods. Neuroimage 2000;11(6 pt 1): 805–821.[CrossRef][Medline]
- Atkins MS, Mackiewich BT. Fully automatic segmentation of the brain in MRI. IEEE Trans Med Imaging 1998;17(1):98–107.[CrossRef][Medline]
- Smith S. Fast robust automated brain extraction. Hum Brain Mapp 2002;17(3):143–155.[CrossRef][Medline]
- Rex DE, Shattuck DW, Woods RP, et al. A meta-algorithm for brain extraction in MRI. Neuroimage 2004;23(2):625–637.[CrossRef][Medline]
- Cox RW. AFNI: software for analysis and visualization of functional magnetic resonance neuroimages. Comput Biomed Res 1996;29(3):162–173.[CrossRef][Medline]
- Dale AM, Fischl B, Sereno MI. Cortical surface-based analysis. I. Segmentation and surface reconstruction. Neuroimage 1999;9(2):179–194.[CrossRef][Medline]
- Segonne F, Dale AM, Busa E, et al. A hybrid approach to the skull stripping in MRI. Neuroimage 2004;22(3):1060–1075.[CrossRef][Medline]
- Lemieux L, Hammers A, Mackinnon T, Liu RS. Automatic segmentation of the brain and intracranial cerebrospinal fluid in T1-weighted volume MRI scans of the head, and its application to serial cerebral and intracranial volumetry. Magn Reson Med 2003;49(5):872–884.[CrossRef][Medline]
- Hu QM, Hou ZJ, Nowinski WL. Supervised range-constrained thresholding. IEEE Trans Image Proc 2006;15(1):228–240.[CrossRef][Medline]
- Bezdek JC, Hathaway RJ, Sabin MJ, Tucker WT. Convergence theory for fuzzy C-means: counterexamples and repairs. IEEE Trans Syst Man Cybern 1987;17(5):873–877.
- Hu QM, Nowinski WL. A rapid algorithm for robust and automatic extraction of the midsagittal plane of the human cerebrum from neuroimages based on local symmetry and outlier removal. Neuroimage 2003;20(4):2153–2165.[CrossRef][Medline]
- Liew AW, Yan H. Current methods in the automatic tissue segmentation of 3D magnetic resonance brain images. Curr Med Imaging Rev 2006; 2(1):91–103.[CrossRef]
- Brandt ME, Bohan TP, Kramer LA, Fletcher JM. Estimation of CSF, white and gray matter volumes in hydrocephalic children using fuzzy clustering of MR images. Comput Med Imaging Graph 1994; 18(1):25–34.[CrossRef][Medline]
- Zhu C, Jiang T. Multicontext fuzzy clustering for separation of brain tissues in magnetic resonance images. Neuroimage 2003;18(3):685–696.[CrossRef][Medline]
- Selvathi D, Arulmurgan A, Thamarai Seivi S, Alagappan S. MRI image segmentation using unsupervised clustering techniques. Presented at the Sixth International Conference on Computational Intelligence and Multimedia Applications, Las Vegas, Nev, August 16–18, 2005.
- Pham DL, Prince JL. An adaptive fuzzy c-means algorithm for image segmentation in the presence of intensity inhomogeneities. Pattern Recognit Lett 1999;20(1):57–68.[CrossRef]
- Hou ZJ, Huang S, Hu QM, Nowinski WL. A fast and automatic method to correct intensity inhomogeneity in MR brain images. Lect Notes Comput Sci 2006;4191(2):324–331.
- Teistler M, Breiman RS, Liong SM, Ho LY, Shahab A, Nowinski WL. Interactive definition of transfer functions in volume rendering based on image markers. Int J Comput Assist Radiol Surg 2007;2(1):52–64.
- Neider J, Woo M. The official guide to learning OpenGL, Version 1.2. OpenGL Architecture Review Board. Amsterdam, the Netherlands: Addison-Wesley Longman, 1999.
- Wernecke J. The inventor mentor: programming object-oriented 3D graphics with Open Inventor, Release 2. Amsterdam, the Netherlands: Addison Wesley Longman, 1994.
- Wilson O, Gelder AV, Wilhelms J. Direct volume rendering via 3D textures. Technical Report UCSC-CRL-94–19. Santa Cruz, Calif: University of California, Santa Cruz, 1994.
- Lee JM, Yoon U, Nam SH, Kim JH, Kim IY, Kim SI. Evaluation of automated and semi-automated skull-peeling algorithms using similarity index and segmentation error. Comput Biol Med 2003; 33(6):495–507.[CrossRef][Medline]
- Zhuang AH, Valentino DJ, Toga AW. Skull-stripping magnetic resonance brain images using a model-based level set. Neuroimage 2006;32(1): 79–92.[CrossRef][Medline]
Top
Abstract
Introduction
