Hyperconnected Openings Codified in a Max Tree Structure: An Application for Skull-Stripping in Brain MRI T1

This article presents two procedures involving a maximal hyperconnected function and a hyperconnected lower leveling to segment the brain in a magnetic resonance imaging T1 weighted using new openings on a max-tree structure. The openings are hyperconnected and are viscous transformations. The first procedure considers finding the higher hyperconnected maximum by using an increasing criterion that plays a central role during segmentation. The second procedure utilizes hyperconnected lower leveling, which acts as a marker, controlling the reconstruction process into the mask. As a result, the proposal allows an efficient segmentation of the brain to be obtained. In total, 38 magnetic resonance T1-weighted images obtained from the Internet Brain Segmentation Repository are segmented. The Jaccard and Dice indices are computed, compared, and validated with the efficiency of the Brain Extraction Tool software and other algorithms provided in the literature.


Introduction
Magnetic resonance imaging (MRI) allows the noninvasive assessment of the patient, and it is useful for early diagnosis, medical monitoring, and the detection of many diseases, such as Alzheimer's disease, brain aneurysm, brain tumor, and melanoma of the eye [1]. However, there are many applications related to brain imaging that require accurate brain segmentation to separate the skull, scalp, dura, eyes, etc. This procedure is known as skull stripping [2,3]. Several applications of brain segmentation include brain volume estimation [4], image registration [5], automatic tumor detection [6], the first state in cortical flattening procedures [7], and structural studies [8][9][10].
Algorithms for skull stripping can be placed into six broad categories mentioned in the literature [1,11]. The categories include: mathematical morphology, pixel intensity analysis, deformable surfaces, brain atlas, hybrid approaches, and deep learning focus.
Furthermore, there are several software programs designed to complete the skullstripping task. One of the most widely used software programs is BET (Brain Extraction Tool) [12], which quickly processes complete volumes of images. BET evolves a tessellated mesh of triangles to fit the brain's surface. However, the resulting crude "skull" image contains a relatively large number of false negatives and positives [13]. Considering the disadvantage mentioned above, BET2 is presented as a new software version, an automated tool for extracting mesh surfaces of the brain, the inner and outer skull, and the scalp from an MRI. Ideally, it requires both T1-and T2-weighted anatomical MRIs, each of a <2 mm resolution in each direction. Another BET improvement is presented in [14], where the authors proposed a faster convergence of the algorithm since they enhanced the vertex displacement, added a new search path, and embedded an independent surface reconstruction process. Other popular algorithms used to separate the skull are the following [15]: BSE (Brain Surface Extractor), SPM2 (Statistical Parametric Mapping v2), McStrip (Minneapolis Consensus Strip), ROBEX (robust brain extraction), and TMBE (Threshold Morphologic Brain Extraction). In [16], the authors used the Richardson-Lucy deconvolution, obtaining high-quality results. However, CoLoRS (coherent local intensity rough segmentation) is a new algorithm to segment an MRI that considers intensity inhomogeneity or bias fields presented in MR volumes [17]. The algorithm is based on clustering and a rough sets theory for simultaneous segmentation and bias field correction of brain MR volumes. The clustering technique allows separating or segmenting components in an image, where statistics are used to group and classify considering texture, color, and form factor, among others. In [18], researchers obtained brain tissue using graph theory, supervoxels, and filtering. In [19], the authors explain the evolution of computational methods in human brain connectivity from 1990 to the present, focusing on graph theory. Graph theory has become a powerful approach for brain imaging analysis, mainly because of its potential to study dynamic behavior over time and disease-related brain changes [20]. In using graph theory, the first step in creating brain graphs is to define the nodes and edges connecting them. It is worth mentioning that connectivity is a powerful concept because the processed images under this notion preserve the contours, avoiding the creation of new maxima or minima.
This paper provides two new openings and introduces two procedures for separating the brain in the MRI T1. Regardless of the method used to separate the skull and nonbrain tissue in an MRI, it is necessary to compare the resultant segmentation to evaluate the performance of the proposed method. The Jaccard and Dice [21,22] indices will be used to compare the resulting segmentations.
The proposals use hyperconnectivity [23], the max-tree [24], viscous transformations [25], and lower levelings [26]. To introduce morphological transformations, Section 2 provides some notions of morphological filtering, connectivity, and hyperconnectivity. Such concepts are explained in detail, and several images illustrate how they operate practically. Section 3 gives the mathematical formalism of the openings defined under the hyperconnectivity. The proposals given in this section are also explained widely to follow the ideas related to mathematical morphology. Section 3.1 shows the criterion based on the maxima of the image to detect the maximum hyperconnectivity, and Section 3.2 formalizes the use of lower leveling applied from the higher extreme. Section 4 illustrates the application of hyperconnected viscous transformations, which allow the skull stripping of the 38 MRI T1 images provided by the IBSR [27]. Once the brain is segmented, the indices of Jaccard and Dice [21,22] are computed and compared with the information obtained from the current literature. Section 5 corresponds to the conclusions. Figure 1 shows the elementary structuring element in 3D containing its origin used in this study.B={−x:x ∈ B} denotes the transposed set with respect to its origin, and λ is a homothetic parameter. The morphological opening, closing, erosion, and dilation transformations are expressed as, γ λB , ϕ λB , ε λB , and δ λB , respectively [28]. Furthermore, the openingγ µ ( f ) = R f (ε)(x) and closingφ µ ( f ) = R f (δ)(x) by reconstruction propagate a marker to filter components by size without introducing new contours [29], where R f represents the reconstruction transformation.

Connectivity
One of the most interesting concepts proposed in mathematical morphology is the connected class introduced by Serra [30]. The transformations defined under specific connectivity do not modify the shapes during the processing, and fused regions preserve the contours of the original image. The connected opening γ x (A) is used in practice to separate each one of the components. Figure 2 clarifies this concept.

Viscous Opening
Serra proposed a connection class on the space generated by dilated [31]. The viscous opening belongs to this class and discovers the connected components eroding the image. This transformation is expressed as follows: The number of connected components depends on the viscosity parameter λ.  Figure 3f presents the connected components for viscosity λ = 22. However, by considering disks as the elementary shapes of the image, it is not possible to detect six connected components for any viscosity. The solution to this problem is to select the connected components at different viscosities (scales), and one option is the traditional algorithm known as ultimate erosion. This consists of choosing the connected components at a specific viscosity λ, such that the viscosity λ + 1 will remove them. Another method to select the connected components is to compute the distance function, as shown in Figure 3g, and detect their maxima. The image in Figure 3h contains the ultimate eroded components for viscosities (sizes) 25, 34, 45, 64, 66, and 68. Figure 3i illustrates the connected components in the viscous lattice sense. Viscous connectivity is interesting because it exploits the goal in binary image segmentation, which consists of splitting the connected components into a set of elementary shapes.

Morphologically Connected Filtering in Viscous Lattices
Equation (1) utilizes the marker x to detect the image components using the punctual opening γ x (X). However, the trivial opening γ O (A) uses other criteria, for example, area, the size of the structuring element, or volume. This opening is presented as follows: The operator γ O (A) detects and recovers all image components, fulfilling an increasing criterion. Then, from Equation (1), the connected viscous opening is expressed as follows: For the case of functions, it is denoted as: Equations (3) and (4) permit different viscous openings depending on the increasing criteria.

Hyperconnectivity
Serra introduced the hyperconnectivity concept [23], which permits working with joined or overlapped components. Figure 4 helps to understand this notion for the 2D case. Notice that the eyes link the brain and the skull in Figure 4a, i.e., they are hyperconnected because they form a unique component. Figure 4b presents the regional maxima obtained from Figure 4a. Each of these maxima are located on the brain, eyes, white matter, or skull. From here, we let Max( f ) be the set of the maxima of f , whereas Max k ( f ) denotes the set of maxima of f at the k level.
Serra [23] considered a class of functions f ∈ F admitting a unique maximal connected component |Max( f )| = 1, where the horizontal cross-sections of functions f ∈ F are connected in such away that functions that admit a unique maximum generate a hyperconnection. An example of this situation is illustrated in Figure 5c,d.

Hyperconnected Opening
Similarly to Serra, we consider the class W of those functions admitting a maximal connected component. Particularly, we define the class W as those extracted from a function f , containing only one maximum using the reconstruction transformation R f (g)(x). We let Max( f ) be the set of the maxima of the function f . The marker h M i is expressed as Thus, a particular set of functions recovered from f with a unique maximum is composed of g M i = R( f , h i ) functions. Then, the set W is defined as follows: Figure 5 exemplifies this situation. Figure 5a shows two overlapping functions, g 1 and g 2 . For the MRI case, g 1 represents a maximum on the brain, and g 2 is a maximum on the skull; however, the regions under the intersection of both functions indicate that the brain and the skull overlap, i.e., they are connected or hyperconnected. Figure 5b illustrates the supremum between g 1 and g 2 . Note that it is not possible to recover g 1 or g 2 from g 1 g 2 . This is what we visualize in reality; our eyes would observe how the brain and the skull appear in two places in a certain slice; nevertheless, lower slices connect them.
The reconstructed functions correspond to g M 1 and g M 2 , displayed in Figure 5c,d. These images come from an individual reconstruction using each maximum computed from Figure 5b. Figure 5c,d illustrate how to detect the markers to separate the brain and skull, and Equation (6) represents it formally. The maximum to be treated is selected and subsequently reconstructed using the transformation by reconstruction R. In Figure 5e, the two functions g 1 and g 2 do not overlap; hence, both functions can be retrieved. However, this is not what really happens in image segmentation. Figure 6 illustrates a real example where there are several hyperconnected functions. Figure 6a,b show the original image and its maxima, respectively. In contrast, Figure 6c illustrates the obtention of a hyperconnected function computed from the extreme region marked with a circle of green color. Note that each maximum can be used to produce a hyperconnected function. ) each regional maximum detected corresponds to different regions within the brain, and each maximum can be used to recover other internal structures. However, our work currently investigates this case. The maximum in a green circle is the marker to produce the image in (c); (c) hyperconnected function associated with regional maximum marked with the green circle. Now, we introduce the trivial criterion to build openings as follows: where Vol represents the volume, and µ v denotes the increasing criterion given by volume.
The opening γ M i is expressed as: Then, we define Equation (7) as: withγ Equation (8) specifies thatγ µ v is obtained from those reconstructed maxima that fulfill the increasing criterion; for this, the supremum operator is necessary.
Formally, the next expression considers the highest maximum, called extreme hyperconnectivity: In the present work, the following connected viscous opening is used: Figure 7 shows an example in 2D using the input image in Figure 7a. The erosion ε λ ( f ) is in Figure 7b. The maxima computed using the max-tree can be found in Figure 7c. Figure 7d shows the maximum fulfilling the increasing criterion to compute the most important hyperconnected component, whereas Figure 7e exemplifies the openingγ µ M v . Figure 7f contains the threshold calculated by the Otsu algorithm. The output image corresponds to the image in Figure 7g, and a mask with the input image permits brain recovery.

Hyperconnected Functions and Lower Leveling
Another useful transformation used in this work to segment the brain is the lower [26]. This transformation works similar to the opening by reconstruction; nonetheless, the α ∈ [0, 255] parameter permits controlling the reconstruction of the marker into the original image.
The proposed marker g considers the regional maxima similar to those defined in Section 3.1, and subsequently, it is iterated using the lower leveling transformation. Therefore, following similar steps to deduce Equation (10), Equation (11) is obtained: with In practice, Equation (11) indicates that the erosion is computed on the input image to separate connected components. Posteriorly, the marker is established on a specific maximum selected by an increasing criterion. The marker g is iterated following the lower leveling operator. During the reconstruction process, the α parameter avoids a complete marker reconstruction into the reference image. In the end, the output image is dilated.

Brain Extraction Based on the Maximum Hyperconnected Function
The difficulties found in the skull stripping procedure and its importance have led to the introduction of a wide range of proposals to address the problems with the procedure. A fundamental challenge to segment the brain in an MRI is the connection between the brain and the skull. Based on classical connectivity, the brain and skull make a connected component. Therefore, to overcome this problem, another type of connection must be used: for example, the viscous connectivity to separate them. The MRI dataset of 38 normal subjects processed in this paper comes from the Internet Brain Segmentation Repository (IBSR), developed by the Center for Morphometric Analysis (CMA) at the Massachusetts General Hospital [27]. However, considering that the images of this database contain intense contrast changes between the different three-dimensional image sections, converting our analysis is a real challenge. This problem, linked to the low quality of images, cannot always be avoided. Thus, robust methods must be implemented to process the images to resolve these difficulties. The proposed morphological transformations presented in this work give good results even when working with poor-quality images. The basic idea behind the segmentation problem of MRI images is the use of viscous connected transformations. Instead of computing the classically connected components on the original image, the morphological erosion ε λ provides them (step one in Figure 8), and the following consists of determining the maximum hyperconnected function (step two in the diagram of Figure 8). Morphological dilation δ λ permits the generation of a viscous component in step three, and a threshold based on the Otsu method (step four) enables brain detection. Figure 7 shows the sequence of operators proposed to separate the brain. The Jaccard (JC = |X Y| |X Y| ) and Dice (DC = 2 |X Y| |X|+|Y| ) [21,22] indices are computed to compare our results with those presented in the current literature. The Jaccard and Dice indices were designed to measure the overlap between two given objects and yield a value between zero (no overlap) and one (complete overlap). The metrics are straightforward to compute and interpret, but they should be applied carefully. These indices are not appropriate for use in the following cases [32]: (1) in small segmented structures such as brain lesions, cell images at low magnification, or distant cars; (2) in the presence of noise; (3) they do not have the capability of distinguish differences in shapes; (4) they are not appropriate for detection and localization tasks. Our interest is evaluating the intersection with the ground truth images provided by databases utilized here. Table 1 displays such indices for the BET algorithm in columns "BET Jaccard" and "BET Dice", and our results are in columns "MHF Jaccard" and "MHF Dice". The data reported in Table 1 use 38 volumes of the IBSR repository, and they are plotted in Figure 15. Table 1. Jaccard and Dice indices computed for 38 volumes of an MRI obtained from the IBSR repository. The brain segmentation in this work utilizes the steps provided in Figure 8, and the computed indices are in columns MHF Jaccard and MHF Dice.

Brain Extraction Based on Hyperconnected Functions and Lower Leveling
The procedure to separate the brain using hyperconnectivity and the viscous opening γ * µ v is similar to that presented in Section 4.1. The diagram in Figure 9 presents the new opening γ * µ v in step 2. Because the lower leveling is iterated with some α value, the marker propagation stops in the image minima more quickly. In this case, the dura represents a relevant minimum in the picture. Therefore, the skull will never reach high-intensity levels, and as a result, the lower leveling transformation produces an excellent outcome. The last step consists of separating the complete reconstructed brain using the Otsu threshold. Figure 10 presents a set of processed images following the steps given in Figure 9. The Jaccard and Dice indices are computed and presented in Table 2. Thus, Figure 11 shows several brain slices in different planes to appreciate the performance of the hyperconnected opening on a specific volume. Moreover, Figure 12 illustrates how the brain is reconstructed by iterating the hyperconnected leveling transformation with different slopes α.
In addition, Figure 13 shows the performance of the algorithm presented in Figure 8 using a file obtained from the Neurofeedback Skull-stripped (NFBS) repository [33]. Ten brain IMR volumes of the NFBS are processed using the MHF algorithm displayed in Figure 8, and their respective Jaccard and Dice indices can be observed in Table 3. The  mean values of Table 3 are displayed in Table 4. Such indices are useful to be compared with other results presented in the current literature.

Figure 9.
Steps to segment the brain using the Hyperconnected Lower Leveling (HLL) procedure. Table 2. Jaccard and Dice indices computed for 38 volumes of MRI obtained from the IBSR repository. The brain segmentation in this work utilized the steps provided in Figure 9, and the computed indices are in columns HLL Jaccard and HHL Dice.  10. Set of images processed using the HLL procedure presented in Figure 9: (a) Input images; (b) viscous components obtained from step (3) using α = 3, λ = 3 and the greater volume computed from a max-tree branch; (c) brain segmentation after applying the Otsu threshold; (d) ground truth images obtained from the IBSR repository. Figure 11. Brain slices in axial, sagittal, and coronal planes obtained by applying the hyperconnected viscous opening γ * µ v , with ε λ=3 ( f ). Table 4 displays the mean values of the Dice and Jaccard indices computed from the number of volumes specified in the last column utilizing different methodologies, and most of them utilize the IBSR dataset. Most articles compare with BET, for which Table 4 includes the indices associated with this methodology. (c) α = 10; and (d) α = 5; (e) contours illustrating the reconstruction process of the brain using the leveling, with, α = 22 in green, α = 10 in blue, α = 10 in fuchsia, and α = 5 in purple. Figure 13. File A00028185 belongs to the NFBS repository. This volume was processed by using the algorithm presented in Figure 8. (a) Segmented brain slices in the sagittal plane with their corresponding ground truth images; (b) segmented brain slices in the coronal plane with their corresponding ground truth images; (c) segmented brain slices in the axial plane with their corresponding ground truth images.

Volume
For the SPM8 Seg, SPM8 VBM, SPM8-NewSeg, FSL, and Brainsuite methods [34], the authors did not report the results for BET. However, by comparing such indices related to the MHF and HLL procedures, a better performance is observed for the algorithms proposed in this paper. The same situation is found when compared to the methodology introduced by Somasundaram et al. [35]. The results from the 10 brain MRIs obtained from the NFBS repository are high because the MRI volumes do not present the eyes, facilitating skull segmentation. Nevertheless, similar performances are presented when MHF and HLL procedures are compared to Zhang et al. [36], Jiang et al. [37], Mendiola et. al. [38], and Galdames et al. [39]. In Jiang et al. [37], the authors use a nonlinear speed function in the hybrid level set model to eliminate boundary leakage. When using the method, an active contour neighborhood model is applied iteratively slice by slice until the neighborhood of the brain boundary is obtained. The results show high values for the computed indices because of the edges enhancement. However, there are two problems: (1) the brain extraction requires a semi-global understanding of the image and (2) the weak boundaries between the brain tissues and surrounding tissues.
Additionally, in Mendiola et al. [38], indices of Jaccard and Dice are high, and so is the algorithm's execution time. To overcome this disadvantage, in this paper, the max-tree structure, and the hyperconnectivity notion are utilized. The results displayed in Table 5 show two faster procedures with indices comparable to those obtained in [37,38]. The results indicate minimal differences between the methods [37,38], and the paper presented here. The proposal introduced here has the disadvantage that before loading the flat areas to the max-tree, smoothing filtering is performed to reduce the amount of regions that are loaded to the structure, producing the fusion of elongated and narrow regions. Furthermore, separately comparing the indices of the averages of the 18 (IBSR1) and 20 (IBSR2) volumes, our procedure overcame the indices reported in [37] for the IBSR1 data set and are similar to the IBSR2. Furthermore, the execution time of our proposal can be found in Table 5, together with the times reported in [38,40].
Nonetheless, as noted in Table 4, the BET method presents low indices when compared to other methodologies. To better appreciate the data provided in Table 4, Figure 14 presents such information. BET is used widely around the world for the following reasons: (i) Usually, researchers do not need complete brain segmentation, and (ii) it is quicker. The computation time, which is not commonly reported in papers, should be considered as an important parameter. With this in mind, Table 5 indicates that the computation time of the proposed HLL method is close to that of the BET algorithm, but with a high degree of efficiency in brain segmentation. This situation is illustrated clearly in Figure 15, where MHF and HLL algorithms outperform BET. The BET algorithm produces approximated segmentations. However, sometimes, when images contain important variations in illumination, brain segmentation is similar to a sphere [38].   Table 1 to the IBSR1 data set; (b) Indices of Jaccard and Dice for the methods displayed in Table 1 related to the IBSR2 data set. Figure 15. Indices of Dice and Jaccard to compare procedures BET, MHF, and HLL. Figures 8 and 9 show diagrams corresponding to MHF and HLL algorithms.

Conclusions
The two methods proposed in this paper utilize hyperconnectivity and viscous lattices, which permit separating the brain, even with poor-quality images. A criterion of maximum hyperconnectivity to extract the main component (the brain) and a threshold based on the Otsu method guarantee automation in the process. The efficiency related to segmentation is better than that of the BET algorithm and similar to the results shown in [36][37][38][39]. Moreover, the computation time is comparable to the BET and quicker than that of the method presented in [38]. However, segmenting hyperconnected components entails separating overlapping parts. If the pixel intensity level between the elements is notable, the possibility of obtaining the right segmentation is high; otherwise, it is necessary to use other criteria to separate components mainly during the reconstruction process. Moreover, the computed indices for volumes of teh NFBS repository are high because the IMR volumes lack eyes, and the segmentation is simple. Furthermore, the better execution time to separate brain corresponds to our proposal when compared with those reported in [37,38], whereas the Jaccard and Dice indices are similar.
Three things should be considered in future work: (1) hippocampal segmentation, where the intensity levels among the regions are similar; (2) white and gray matter separation, since they are related to neural damage and memory problems due to aging; and (3) the proposed transformations work adequately for T1 images, and not for other modalities. This inconvenience will be solved.