Algorithm for brain extraction on Magnetic Resonance Images T1 using Morphological 3D Transformations

Algoritmo para la extracción de imágenes de resonancia magnética T1 usando transformaciones morfológicas en 3D

J. D. Mendiola-Santibañez* I. M Santillán Méndez ** C. Paredes Orta*** I. R Terol Villalobos ^****

*,***,**** CIDETEQ Centro de Investigación y Desarrollo Tecnológico en Electroquímica.

Correspondencia:
J. D Mendiola Santibañez
Parque Tecnológico S/N, San
Fandila-Pedro Escobedo, CP
76703, Querétaro, México. Tel.
+52 442 2116027.
Correo electrónico: mendijor@uaq.mx

ABSTRACT

In this paper a 3D morphological composition of transformations for brain extraction on brain Magnetic Resonance Images T1 (MRI T1) is presented. The proposal makes use of two morphological connected transformations, the lower leveling and a family of the viscous alternating sequential filters (VASFs). The properties of these operators -which consist in the control of the reconstruction process of a marker into the original image-, are exploited to segment the brain in 20 volumes of MRI T1. The segmented brains are compared with respect to: i) the segmentations obtained from BET which is popular among the scientific community for segmenting the brain; and ii) manual segmentations. The computed indices indicate that the proposed transformation produces good results during its performance. The consumed time for the algorithm during the execution is acceptable and it can be implemented in Matlab.

Keywords: brain segmentation, connected transformations, filtering, viscous transformations.

RESUMEN

En este artículo se presenta una composición de transformaciones morfológicas para la extracción del cerebro en imágenes de resonancia magnética T1 (MRI T1) en 3D. La propuesta hace uso de dos transformaciones morfológicas conexas, la nivelación inferior y una familia de filtros viscosos alternados secuenciales (VASFs). Las propiedades de estos operadores- las cuales consisten en el control del proceso de reconstrucción de un marcador dentro de la imagen original-, se explotan para segmentar el cerebro de 20 volúmenes de MRI T1. Los cerebros segmentados se comparan respecto a: i) segmentaciones obtenidas a partir del algoritmo BET, el cual es popular en la comunidad científica para la segmentación del cerebro; y ii) segmentaciones manuales. Los índices calculados indican que la transformación propuesta produce buenos resultados en su desempeño. El tiempo empleado por el operador durante su ejecución es aceptable además de que la propuesta puede ser implementada en Matlab.

Palabras clave: segmentación del cerebro, transformaciones conexas, filtrado, transformaciones viscosas.

INTRODUCTION

Mathematical Morphology (MM) is a technique where the images can be analyzed from a structural point of view. A special field of application of the MM is the segmentation of medical imaging. In particular, in this paper a composition of two morphological transformations defined for the 3D case and capable of separating brain tissue from skull is presented as proposal. The composed transformation takes into consideration the lower leveling[1] and the viscous alternated sequential filters[2]. The lower leveling propagates a marker (the marker is a part of the original image obtained from an operator) at the interior of the original volume in a controlled manner. Hence, the output image is contained into the original one. While the viscous alternated sequential filter allows disconnecting and smoothing chained components[2]. The viscous name appears because when chained components are separated, the output image visually resembles with a viscous effect.

The proposal is tested on 20 volumes of MRI T1 taken from the IBSR database[3] to separate brain. In order to compare our results, the Brain Extraction Tool (BET) -which is a popular skull-stripping algorithm- will be used. An important characteristic of BET is its quickness to separate the skull from MRI volumes of the brain. Other algorithms reported in the current literature used to carry out the same task can be found in references [4,5,6,7,8,9,10,11,12]. The composed operator provided in this paper is similar to that presented in [12]. However, there are two significant differences, i) here, the composition is obtained from the application of the lower leveling and the alternated sequential filters; while in [12], the lower leveling and the viscous opening are utilized; and ii) the execution time to segment the brain is improved in this paper.

The computed segmentations will be compared using the Jaccard[13] and Dice[14] indices. It is noteworthy to mention that, in the literature, two problems are associated during the deskulling process, 1) time used to separate brain, and 2) the quality of the segmentation. For example, in [12], the Jaccard and Dice indices have values around 0.93 and the proposal given to separate brain consumes around 15 minutes. In this paper it is introduced a morphological transformation capable of processing a MRI T1 in 3D with 60 slices in 83 seconds with Jaccard indices around 0.93.

This paper is organized as follows, in section Methodology some basic morphological transformations are defined which will be useful to introduce our proposal. In section Results and Discussion, a comparison is made between the segmented brains using our proposal and those obtained from: i) manual segmentations provided by the IBSR database[3]; and ii) BET algorithm implemented in the MRicro software[15]. In Discussion section, the advantages and disadvantages of our proposal are presented. The conclusions correspond to the last section.

METHODOLOGY

Opening by reconstruction

In MM, the basic transformations are the erosion ε_μB(f) and dilation δ_μB(f) [16], where B represents the 3D structuring element which contains its origin in the center. Figure 1 illustrates the shape of the structuring element used in this paper. denotes the transposed set of B respect its origin, = {—x : x ∈ B}, μ is a size parameter, f is the input image defined on Z³, Z represents the integer set and x is a point on the definition domain.

The morphological erosion and dilation are expressed as follows:

where ∧ and ∨ are the inf and sup operators. The morphological erosion and dilation permit to build other type of transformations; these include the morphological opening_μB(f) and closing _μB(f):

In order to reduce the notation, consider that μB = μ. The transformations having the characteristics of modifying the regional maxima or minima, without considerably affecting the remaining components, are the opening and closing by reconstruction. These operators use the geodesic transformations[17]. The geodesic dilation and erosion are expressed as = f ∧ δ_B (g) with g ≤ f, and = f ∨ ε_B (g) considering g ≥ f, respectively. When function g is equal to the morphological erosion or dilation, and the basic geodesic transformations are iterated until stability, the opening or closing by reconstruction are obtained,

Viscous opening and closing

Viscous opening and closing utilize the parameters λ and μ. The ε_λ(f) and δ_λ(f) are used to detect those components that can be separated. To illustrate this situation consider the image in Fig. 2. Original image is presented in Fig. 2(a). In Fig. 2(b), the erosion size λ = 5 disconnects chained components; posteriorly in Fig. 2(c) the opening by reconstruction eliminates the components of size minor than μ — λ = 3 while the remaining components are maintained unchanged. Finally, in Fig. 2(d) the dilated image with λ = 5 recover the size of the original components but now they are disconnected. The components that do not support the erosion and the opening by reconstruction are merged with the background.

Viscous Alternating Sequential Filters (VASFs)

Serra[18] defined and characterized four operators where the size of the structuring element μ is indexed over a size distribution with 1 ≤ μ ≤ k.

Let us take one of these operators defined as follows in terms of connected viscous openings and closings[2]:

For fixed parameters λ₁, λ₂, with λ₁≤ λ₂ ≤ μ,

From last equation and for a family {λ_i}, with λ_j≤ λ_k if j < k, the following VASF is defined[2]:

with the condition λ_n ≤ μ.

The VASFs will be used to segment some regions of interest. The segmentation process consist in flat the maxima and minima iteratively through the viscous opening and closing. At the end with a threshold is posible to recover the flatened zones under study. Fig. 3 illustrates the performance of N_λn,μ. Notice that, the brain has been completely separated.

Nevertheless, the problem of Eq. 3 is that not always produces the same result and parts of the skull can be obtained in the final segmentation. Due to this situation, the VASFs will be used with the lower leveling1 to improve the final result. Lower leveling transformation is presented below.

Lower Leveling

The extended lower leveling is expressed as follows[1]:

where f is the reference image, g is a marker, α ∈ [0, 255] is a positive scalar called slope, μ and i is the size of the structuring element. The marker g used in this paper corresponds to,

The Eq. 5 is iterated until stability with the purpose of reconstructing the marker g at the interior of the original mask f, i.e.:

The image in Fig. 4, ilustrates the performance of Eq. 6. In Fig. 4(b) a marker obtained from the internal gradient is shown. The internal gradient is defined as the arithmetic difference between the original image and the eroded one, i.e., gradi_μ(f) = f — ε_μ(f). The internal gradient uniquely detects regions on the skull. The idea of applying Eq. 6 using as marker the internal gradient is to reconstruct precisely the skull. The output image after applying Eq. 6 with α = 3 is displayed in Fig. 4(c). In this image the skull has been separated and reconstructed. Hence the arithmetic difference f — (f,g), allows us to recover the remaining components, i.e., brain, dura matter, among others, as is illustrated in Fig. 4(d). The problem of separating brain in this way is that the intensity levels of the white matter are similar to the skull bringing as consequence that the marker cannot be obtained uniquely on the bone.

On the other hand, the appropriated size of the structuring element to detect the marker utilized in this paper is found through a granulometric study. The granulometry definition is introduced below.

Granulometry

Granulometry is the distribution by sizes of particles that constitute an aggregate; it is employed in diverse areas to describe the qualities of size and shape of individual grains within a product. The concept of granulometry was introduced by G. Matheron at the end of the sixties and is presented as follows[19].

Definition (Granulometry). Let (ψ_λ≥0) _λ∈Z+ be a family of transformations depending on an unique positive parameter λ. This family constitutes a granulometry ifand only if the next three properties are verified: (i)∀ positive λ, ψ_λ is increasing; (ii) ∀ positive λ, ψ_λ is antiextensive, and (iii) ∀ positive λ and μ, ψ_μ, = ψ_μψ_λ = ψ_max(λ,μ).

The family of morphological openings and closings for the numerical case {_μ}, {_μ} with μ = {1,...,n} fulfills this last definition. The volume distribution x is computed with Eq. 7:

where vol represents the volume of f, i.e., the sum of all intensity levels into the image. Figure 5 presents a granulometric study applied on a volume of MRI T1 segmented in this paper. Notice in this figure that, for sizes of μ between 1-5 the openings detect a component corresponding to the skull while openings between sizes μ in the interval 5-21 detect an important presence of the brain. From this analysis, appropiated sizes of μ to detect the marker of the brain using the morphological opening are within the interval μ ∈ [6, 22].

The operators presented so far will be used to introduce a methodology to separate brain in volumes of MRI T1 for the 3D case.

Methodology to segment brain

The idea to separate brain and skull consists in combining Eq. 3 and Eq. 6 using as marker the Eq. 5. The separation of the brain consists in: i) to compute a marker of the brain using the morphological opening; ii) posteriorly, the lower leveling will reconstruct the portion of the brain into the head. However, the slope α will avoid the complete reconstruction of the marker at the interior of the volume; and iii) to clean those regions adhered to the reconstructed brain, Eq. 3 will be applied together with a threshold operation.

Formally, next filter expresses the aforementioned in steps i, ii and iii,

where _μy (f) indicates that a portion of the brain will be reconstruted into the head f through the transformation , and lately all regions adhered to the brain are smoothed by applying the VASFs, i.e., .

Equation 8 is in term of the parameters λ_n, μ_z, μ_x, α and μ_y. Now, an explanation of how selecting the values for these parameters is provided.

]]> μ_y: This parameter represents the size of the structuring element of the morphological opening _μy (f) used to obtain an initial marker of the volume f. According to the granulometric study presented in Fig. 5 (see section Granulometry), adequate sizes for μ_y are contained in the interval μ_y ∈ [6, 22]. In particular, in this paper, μ_y= 12 is utilized.

μ_x: In this paper, μ_x takes the value of 1 to reconstruct the marker_μy(f) as close to the original volume f.

α: This parameter is used in Eq. 6 and its purpose is to control the propagation of the marker at the interior of the original volume f. When α increases, less volume is reconstructed. In Fig. 6 this effect is illustrated. Adequated values for a according to the image in Fig 6 are in the interval 6 to 12. With a in this interval, the skull component is almost eliminated. In our case α = 10 is considered.

λ_n, μ_z: These parameters correspond to the VASFs. The task of the VASFs is smoothing the remaining components of the skull after the application of the lower leveling.

An example of this situation can be seen in Fig. 6(g) where the lower leveling is applied considering α = 9. This example is important because it illustrates the approximate size of the remaining components of the skull that should be eliminated. As described before when the VASFs were introduced, once has been considered the dimension of the bone to be removed, the complete skull is supressed by applying the VASFs with λ_n= 2 and μ_z= 5 (see Fig. 3). However, due to the volumes processed with the VASFs will contain less information after applying the lower leveling, then the μ_z parameter can be reduced to 3 instead of 5 with the purpose of reducing the execution time. If λ_n = 2, the minimum value taken by μ_z is 3 because of next relation must be fulfilled, λ₁ ≤ λ₂≤ μ_z, see section of Viscous Alternating Sequential Filters.

 

RESULTS AND DISCUSSION

In this paper, 20 MRIs T1 (denoted as IBSR) of the head taken from the Centre for Morphometric Analysis (CMA) at the Massachusetts General Hospital (http://www.nitrc.org)[3] are processed with Eq. 8. The parameters utilized to separate the brain are shown in Table 1. From Table 1 next comments are presented:

i. The parameters λ_n, μ_z, μ_x, α, μ_y associated to Eq. 8 are similar for all segmented brains.

ii. The parameter Th₁ represents a threshold used to eliminate cerebrospinal fluid with the purpose of disconnecting several regions between the brain and skull. According to the results presented in 20, the cerebrospinal fluid is eliminated with threshold values among 20-255 sections. Fig. 7 illustrates this situation. Original image is located in Fig. 7(a), whilst the threshold among 20-255 sections is presented in Fig. 7(b). Notice that, several components among the brain and the skull have been merged with the background, while the remaining regions are maintained unchanged. In Fig. 7(c) it is presented a mask between the original image (7(a)) and the image in Fig. 7(b). The image in Fig. 7(c) is convenient during the processing because when the leveling transformation reaches the regions where cerebrospinal fluid was located (dark regions); the propagation of the marker will be slower, obtaining in this way almost the separated brain with the leveling transformation.

iii. The selected size of the opening to obtain a marker is μ_y= 12. This value belongs to the interval 6 to 22 according to the graph in Fig. 5.

iv. The parameter λ_n= 2 indicates that the following filter will be applied:

with μ fixed to 3.

v. The parameter Th₂ allows recovering the brain region. The threshold value is according to the information provided in [20], to eliminate dark regions.

a. Load the 3D image and assign to f

b. threshold1 =1*( f >Th₁)

c. h = f .* threshold11

d. μ_y= 12; α = 10, λ_n= 2; μ_z= 3

e. Compute the marker g = _μy (h)

f. Iterate Eq. 6 using g = _μy (h), the result is kept in J that represents α

g. Apply the filter N_λn,μ to the image J and the resultant image is kept in R that represents ξ_{λn, μz, μx, α, μy}

h. threshold 2 =1*(R >Th₂)

]]> i. brain=h.* threshold 2

In Fig. 8 several of the steps previously mentioned are illustrated. For implementing the step f of the algorithm, the condition used as stop criterion is the volume, i.e., when the volume of the processed image does not change then the idempotence or stability has been reached.

In Fig. 9 some output images belonging to the volume IBSR_7 (see Fig. 9(a)) and processed with the algorithm aforementioned are exhibited, see Fig. 9(c). In order to compare our results, the Jaccard and Dice indices are computed on: i) manual segmentations obtained from the IBSR volumes, and ii) segmentations obtained from the BET algorithm implemented in the MRicro software with default parameters.

Some slices of the manual segmentation corresponding to the volume IBSR_7 are presented in Fig. 9(b), while a set of slices belonging to the volume segmented with the MRicro software are displayed in Fig. 9(d). The Jaccard and Dice indices associated to the aforementioned segmentations are presented in Table 2.

According to the computed indices, the Eq. 8 outperforms the BET algorithm. However there are some points to mention:

1. BET algorithm is faster and only requires 2 seconds to segment the volume. Nevertheless, during the segmentation of the brain, this is cropped like a circle, due to this, the Jaccard and Dice indices have minor values than those computed with Eq. 8.

2. The proposed algorithm utilizes 83 seconds to segment a brain with 60 slices, and although the time is superior to the BET algorithm, the segmentations are better. In 12, the maximum time reported for other methodologies to separate brain is 35 minutes. Also in [12], two expressions to segment brain are proposed, the Eq. 9 and Eq. 12. The results obtained in this paper outperform those obtained with Eq. 12 and the reduction of the execution time is considerable. Our algorithm utilizes 83 seconds while Eq. 9 and Eq. 12 employ 174 and 354 seconds to process a volume of 60 slices respectively.

]]>  

CONCLUSIONS

In this paper a morphological transformation is presented to segment brain from MRI T1 of the head. The Jaccard and Dice indices indicate that our proposal provides better segmentations than those obtained from the BET algorithm which is very popular software among the scientific community. The time used to separate the brain with our proposal is adequate (83 seconds for a volume of 60 slices), in addition the algorithm can be implemented with Matlab because several of the required transformations are already implemented in such software.

REFERENCES

1. Mendiola-Santibañez JD, Terol-Villalobos IR, Jimenez Sanchez AR, Gallegos-Duarte M, Rodriguez-Reséndiz J, Santillán I. Application of morphological connected openings and levelings on magnetic resonance images of the brain. Int J Imaging Syst Technol, 2011; 21(4): 336-348.         [ Links ]

2. Santillán I, Herrera-Navarro AM, Mendiola-Santibáñez JD, Terol-Villalobos IR. Morphological Connected Filtering on Viscous Lattices. Journal of Mathematical Imaging and Vision, 2010; 36(3): 254-269.         [ Links ]

3. Center for Morphometric Analysis. Massachu setts General Hospital, The Internet Brain Segmentation Repository (http://www.nitrc.org/projects/ibsr), IBSR 2014.         [ Links ]

4. Liu JX, Chen YS, Chen LF. Accurate and robust extraction of brain regions using a deformable model based on radial basis functions. Journal of Neuroscience Methods, 2009; 183(2):255-266.         [ Links ]

5. Smith SM. Fast robust automated brain extraction. Hum Brain Mapp, 2002; 17(3):143-155.         [ Links ]

6. Hahn H, Peitgen H-O. The skull stripping problem in MRI solved by a single 3D watershed transform. In Proc MICCAI, LNCS, 2000; 1935:134-143.         [ Links ]

7. Beare R, Chen J, Christopher AL, Silk T, Thompson DK, Yang JYM, et al. Brain extraction using the watershed transform from markers. Frontiers in Neuroinformatics, 2013; 7(32):1-15.         [ Links ]

8. Dogdas B, Shattuck DW, Leahy RM. Segmentation of skull and scalp in 3-D human MRI using mathematical morphology. Hum Brain Mapp, 2005; 26(4):273-285.         [ Links ]

9. Sandor S, Leahy R. Surface-based labeling of cortical anatomy using a deformable database. IEEE Trans Med Imaging, 1997; 16(1):41-54.         [ Links ]

10. Ségonne F, Dale AM, Busa E, Glessner M, Salat D, Hahn HK, Fischl B. A hybrid approach to the skull stripping problem in MRI. NeuroImage, 2004, 22(3):1060-1075.         [ Links ]

11. Iglesias JE, Liu CY, Thompson PM, Tu Z. Robust Brain Extraction Across Datasets and Comparison With Publicly Available Methods. IEEE Transactions on Medical Imaging, 2011; 30(9): 1617-1634.         [ Links ]

12. Mendiola-Santibañez JD, Gallegos-Duarte M, Arias-Estrada MO, Santillán-Méndez I, Rodríguez-Reséndiz J, Terol-Villalobos IT. Sequential application of viscous opening and lower leveling for three-dimensional brain extraction on magnetic resonance imaging t1. J. Electron. Imaging. 0001; 23(3):033010. doi: 10.1117/1.JEI. 23.3. 033010, 2014.         [ Links ]

13. Jaccard P. The distribution of the flora in the alpine zone. New Phytologist, 1912;11(2):37-50.         [ Links ]

14. Dice LR. Measures of the amount of ecologic association between species. Journal of Ecology, 1945; 26(3): 297-302.         [ Links ]

15. Rorden C, Brett M. Stereotaxic display of brain lesions. In Behavioural Neurology, 2000; 12(4): 191-200.         [ Links ]

16. Heijmans H. Morphological Image Operators. Academic Press (USA), 1994.         [ Links ]

17. Vincent L. Morphological grayscale reconstruction in image analysis: applications and efficient algorithms. IEEE Transactions on Medical Imaging, 1993; 2(2):176-201.         [ Links ]

18. Serra J. Image Analysis and Mathematical Morphology. Theoretical Advances, vol. 2. Academic Press (San Diego), 1988.         [ Links ]

19. Matheron G. Random sets and integral geometry. En: John Wiley and Sons (New York), 1975.         [ Links ]

20. Mendiola-Santibáñez JD, Terol-Villalobos IR, Herrera-Ruiz G, Fernández-Bouzas A. Morphological contrast measure and contrast enhancement: One application to the segmentation of brain MRI. Signal Processing, 2007, 87(9): 2125-2150.         [ Links ]