1. Introduction

The term brain tumor refers to an abnormal growth of cells in the brain. It could be cancerous or not and could originate directly from the brain tissues or have a metastatic nature. Benign brain tumors have a uniform structure and do not contain active cells; on the contrary, malignant brain tumors have a non-uniform structure and contain active (cancerous) cells. All the types of brain tumors may produce symptoms that vary according to the parts of the brain involved. Moreover, brain tumors may vary in size and location, besides presenting a lot of abnormalities which make it difficult to characterize and identify them.

Magnetic Resonance Imaging (MRI) is a very common tool used in neurology for visualizing the brain anatomy along three different planes: axial, coronal, and sagittal. When protons are placed in a magnetic field, they align. This alignment (magnetization) is then perturbed by introducing an external radio frequency energy. The protons start oscillating and, at this state, they can absorb energy and hence the nuclei can release it or re-radiate it. The nuclei then return to their original equilibrium state: the time that it takes for the relaxation of the component parallel to the magnetic field to go back to the equilibrium is called longitudinal relaxation time or T1; on the other hand, the transversal relaxation time is abbreviated with T2. The transmitted signals are measured; then, the Fourier transform is used to convert the frequency information to corresponding gray scale intensity values: this makes it possible to construct an image. By varying the sequence of the radio frequency pulses applied and collected, different types of images can be created. In particular, calling repetition time (RT) the time interval between successive pulse sequences applied to the same slice, and calling time to echo (RE) the time interval between the delivery of the pulse and the receiving of the echo signal, the so called T1-weighted and T2-weighted images can be created. Images T1-weighted are produced using short TE and TR, while images T2-weighted are produced using longer TE and TR.

The objective of this paper is to describe an automatic methodology to identify any possible brain tumor given either a T1-weighted or a T2-weighted MRI. By screening the given image, the proposed approach locates the tumor area as well as identifies those scans that do not present any abnormal tissues, and hence, can be categorized as healthy brains. The method described here relies on the use of a finite-elements approximation of the input image and on the use of the Zienkiewicz-Zhu (Z2) error estimator [1,2,3] to adaptively and locally refine the resolution of the image at specific sites. In particular, any area of strong change in the gray intensity values will be easily identified by the Z2 error estimator, allowing for highlighting the abnormal tissues.While the use of finite-elements and adaptive strategies is not new in the image segmentation context, see [4,5,6,7,8], these methods are rather complex as they are framed in a variational formulation, and their objective is different from the one proposed here: the given image is split into several parts or components that can be easily identified via the location of edges. Hence, these methodologies aim at locally reconstructing the profiles of these regions, but they do not address the problem of the classification of such regions.In general, regarding the brain tumor detection task, many methods are developed by using statistical and machine learning techniques; see [9,10]. In particular, advanced architectures like convolutional neural networks (NN) [11], probabilistic NN [12,13], deep NN [14], and U-Nets [15] have been successfully applied to the task of image segmentation and brain tumor classification; see also [16] for a comprehensive review. On the other hand, many unsupervised methods exist as well; some are based on clustering techniques [17,18,19,20,21,22]; others are based on automatic thresholding methods and morphological transformations [23,24,25,26,27]. Some techniques apply a suitable projection in lower dimensional spaces in order to get rid of unwanted redundancies, and so classical matrix decompositions like the SVD or the non-negative matrix factorization are applied; see [28]. A histogram-based gravitational optimization algorithm [29], Gabor-wavelets [30], and other big data analytics techniques are also employed; see the review [31] and references therein.Many other approaches are based on supervised anomaly detection: the main goal is to identify abnormalities, given only normal data in the training set, see [32,33,34,35,36,37]. The method introduced here can be interpreted as an unsupervised anomaly detection technique. Besides using the Z2 error estimator, it also exploits specific γ corrections which help to increase or to reduce brightness of the normal tissues surrounding the tumor, according to which modality of MRI is analyzed, i.e., either T1-weighted or T2-weighted. As a final step, the use of morphological operators allows to get rid of redundant noise and to locate the affected, i.e., anomalous, areas inside the brain tissues.Recent approaches have also exploited particle swarm optimization [38], genetic and evolutionary algorithms; see [39] and references therein. The main contributions proposed in this paper can be summarized as follows:

A suitable γ correction is applied to pre-process the data;

The use of finite-elements and Z2 error estimator with isotropic mesh refinement allows for automatically detecting the anomalous areas, if present;

A post-processing phase, based on the use of morphological transformations, allows for getting rid of the cortex and better locating the anomalous tissues, if present.

The paper is organized as follows: In Section 2, the finite-elements mathematical framework is described. In Section 3, the adopted pre-processing and post-processing are discussed. In Section 4, the method is tested and compared with other available techniques. Section 5 provides some conclusive remarks. 2. Mathematical FrameworkProviding as input an MR image I having size m×n pixels, the following methodology produces a finite-elements approximation on a locally refined regular triangulation. In particular, a denser number of triangles will be automatically located where the anomalous tissues are present, allowing for an easy detection of such regions. The first step consists of constructing a regular quadrangular grid associated with I, i.e., every pixel is a node in such a grid. Then, every regular square is diagonally divided in two triangles. In this way, a non-overlapping regular triangulation Th can be obtained, where the subscript h denotes the maximal diameter of the triangular elements. At this early stage, hT=diam(T)=h constant for every triangle T∈Th. More in general, diam(T):=maxx,y∈T∥x−y∥ and hence h:=maxT∈ThhT. On the constructed triangulation, the space of linear finite-elements can be defined:

Vh1:=vh∈C0(Ω-):vh|T∈P1,∀T∈Th,

(1)

where Ω:=(0,m)×(0,n); P1 is the space of polynomials of a degree less than or equal to 1, and hence Vh1 denotes continuous piecewise linear polynomial functions on every triangle T. Every function vh∈Vh1 can be expressed as a linear combination of the basis functions  as

vh(x)=∑k=1Nhvkφk(x),

(2)

with vk values of the function vh at the grid nodes and Nh total number of grid nodes in Th. In our case, a continuous approximation image I is constructed by interpolating the original image I at the provided pixel values, i.e., via using expression (2),

I(x)=∑k=1NhIkφk(x),

(3)

where Nh denotes the total number of nodes in Th, and a linear global index k is used to identify every pixel (i,j) in I, according to the mapping shown in Figure 1. In particular, k can be obtained as As a second step, a downsampling of the constructed approximation I in Equation (3) is performed. In particular, I is approximated itself via using a coarser triangulation T˜h˜ having N˜h˜ vertices, where the set of these nodes is not necessarily a subset of the original triangulation nodes. This new approximation shall be denoted as I˜h˜. The third step consists of computing the Zienkiewicz-Zhu error estimator ηT˜ on every triangle T˜∈Th˜:

ηT˜:=∥G(I˜h˜)−∇I˜h˜∥L2(T˜).

(4)

In our case, the function G(I˜h˜) in Equation (4) is constructed by evaluating the gradient ∇I˜h˜ on the barycenter of each triangle and then by averaging on those triangles which share a common vertex. In particular, the reconstruction Gk˜ of the gradient at every node k˜ of the triangulation T˜h˜ is done by

Gk˜:=∑T˜∋k˜|T˜||Tk˜|∇I˜h˜|T˜,

(5)

where the symbol |·| in formula (5) denotes the area of the considered element, and Tk˜ denotes the union of the triangles which share the vertex k˜. The elements T˜ having a bigger ηT˜ are the ones where the intensity value of the pixels varies the most. Then, the estimators for every triangle are sorted from the smallest to the biggest one and the triangles that are marked, and hence refined, are the ones such that with t:=pNt, Nt the total number of triangles and p percentage of triangles that will be considered, while 0<cref<1. The marking strategy is fundamental to guarantee a specific order of convergence for the chosen approximant; see [40,41]. However, this is not the main concern to be addressed for the current application. The parameters setting for Equation (6) in the tests has been experimentally carried out.Remark 1. Although more accurate formulas exist for the gradient recovery estimation (see [42,43,44,45]), in this work, the choice of Equation (5) is intentional: a more rough approximation is indeed preferable in this context rather than an extremely accurate one, as the main goal of the proposed methodology is to estimate the region extension of any possible anomalous tissue, rather than just delineating the profile or the edges of such regions.The refinement strategy is repeated for several iterations, until the maximum ηT˜ does not reach the chosen threshold value, producing every time an isotropic triangulation. At the end, smaller triangles will be located in the regions of interest. Therefore, a final binary map b is produced by constructing b∈Vh˜0 as a piece-wise constant polynomial function, which is defined as

b|T˜:=1if|T˜|is minimal,0otherwise.

(7)

It should be noticed that having an isotropic triangulation is a fundamental step in order to produce feasible results with the described strategy. Indeed, having an anisotropic triangulation does not necessarily guarantee that the regions of interests are covered with triangles of a small area. The shape of the triangles can be rather elongated and flattened, preventing the use of the rule in Equation (7) to produce a suitable binary map b. A Metastatic Brain Tumor ExampleTo make the described procedure clearer to the reader, the following example shows the various employed steps. An image 283×295 pixels of a metastatic brain tumor is selected, Figure 2a. The first interpolant, I, obtained with linear finite-elements and defined with 27,828 degrees of freedom (dof), is shown in Figure 2b, while the starting approximation I˜h˜ for the conducted analysis, on a down-sampled triangulation consisting of Nt=1200 triangles and 217 dof, is shown in Figure 2c. The triangulation is locally refined by using the error estimator in Equation (4) and setting the parameters p=0.25 and cref=0.60 for the Equation (6). The output of this process in shown in Figure 2d,e, where the produced approximant I˜h˜ is defined with 936 dof. In Figure 2f, only for visualization purposes, the original MRI bright regions are superimposed onto the final triangulation in order to show the coherence of the obtained results. Finally, in Figure 2g, the final triangulation is shown again, but the triangles are colored according to their area value. In Figure 2h, the produced binary map b is constructed from the final triangulation by selecting only the triangles with a minimal area and setting the corresponding pixels to 1, while setting any other pixel value equal to 0. It is evident how there is still some noise affecting the quality of the produced results. Therefore, this final image will be post-processed as described in Section 3.2. 5. ConclusionsIn this paper, a method based on finite-elements and on the Z2 error estimator is presented to address the task of tumor detection in brain MR images. The proposed approach is called Z2-γ as it also relies on a suitable preliminary gamma correction which helps to better highlight the anomalous regions, whenever present. This initial study shows promising results in the application of the Z2 error estimator, together with isotropic triangulations, and specific tools of image processing in the automatic detection of brain tumors. Future work is currently devoted to provide fully automatic settings for the parameter γ, according to which input image is given, and to the study of brain–tumor segmentation in its various parts, i.e., active cells, necrotic core, and edema, see [58,59,60] and references therein, a task which requires a better discriminative power of the error estimator, as well as suitable input images; see [61]. Moreover, in the current approach, the final validation should always be carried out by the radiologist, as the method achieves 84.6% accuracy, and hence some MR images are still not correctly classified. In order to make the participation of the medical staff more active, for future upgrades, the segmentation task could be performed by combining both classical MRIs and functional-MRIs, and the initial screening of possible affected areas could be completed by the radiologist by drawing specific curves, as it is proposed in [62]. Although the resulting method would still not be fully automatic, the author believes that it is always important to rely on the supervision of competent medical staff as well as it is fundamental to develop suitable technologies to improve the accuracy and precision of any initial formulated diagnosis.