# A Probabilistic Atlas of Diffuse WHO Grade II Glioma Locations in the Brain

Sarah Parisot, Amélie Darlix, Cédric Baumann, Sonia Zouaoui, Yordanka Yordanova, Marie Blonski, Valérie Rigau, Stéphane Chemouny, Luc Taillandier, Luc Bauchet, et al.

<figure><img src="image_1.png" /></figure>

► To cite this version:

Sarah Parisot, Amélie Darlix, Cédric Baumann, Sonia Zouaoui, Yordanka Yordanova, et al. A Probabilistic Atlas of Diffuse WHO Grade II Glioma Locations in the Brain. PLoS ONE, 2016, 11 (1), pp.e0144200. 10.1371/journal.pone.0144200. hal-01882093

HAL Id: hal-01882093

https://hal.science/hal-01882093

Submitted on 26 May 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
RESEARCH ARTICLE

# A Probabilistic Atlas of Diffuse WHO Grade II Glioma Locations in the Brain

Sarah Parisot¹,²,³✉*, Amélie Darlix⁴,⁵✉, Cédric Baumann⁶, Sonia Zouaoui⁷, Yordanka Yordanova⁸, Marie Blonski⁹, Valérie Rigau¹⁰, Stéphane Chemouny³, Luc Taillandier⁹, Luc Bauchet⁵,¹¹, Hugues Duffau⁵,¹¹, Nikos Paragios¹,²

1 Center for Visual Computing, Ecole Centrale Paris, Chatenay Malabry, France, 2 INRIA, Galen Team, Saclay-Ile-de-France Center, Chatenay Malabry, France, 3 Intrasense SAS, Montpellier, France, 4 INSERM U1051, Montpellier Neurosciences Institute, University Hospital, Montpellier, France, 5 Department of Medical Oncology, Institut Régional du Cancer de Montpellier - Val d'Aurelle, Montpellier, France, 6 Department of Clinical Epidemiology and Evaluation, University Hospital, Nancy, France, 7 Department of Epidemiology, Groupe de Neuro-Oncologie du Languedoc-Roussillon, Registre des Tumeurs de l'Hérault, Institut Régional du Cancer de Montpellier - Val d'Aurelle, Montpellier, France, 8 Department of Neurosurgery, HIA Val de Grâce, Paris, France, 9 Neuro-oncology Unit, Department of Neurology, University Hospital, Hôpital Central, Nancy, France, 10 Department of Pathology, University Hospital, Hôpital Gui de Chauliac, Montpellier, France, 11 Department of Neurosurgery, University Hospital, Hôpital Gui de Chauliac, Montpellier, France

✉ These authors contributed equally to this work.

✉ Current address: Biomedical Image Analysis Group, Department of Computing, Imperial College London, 180 Queen's Gate, London SW7 2AZ, United Kingdom

* s.parisot@imperial.ac.uk

<figure><img src="image_3.png" /></figure>

CrossMark
click for updates

<figure><img src="image_4.png" /></figure>

OPEN ACCESS

**Citation:** Parisot S, Darlix A, Baumann C, Zouaoui S, Yordanova Y, Blonski M, et al. (2016) A Probabilistic Atlas of Diffuse WHO Grade II Glioma Locations in the Brain. PLoS ONE 11(1): e0144200. doi:10.1371/journal.pone.0144200

**Editor:** Joseph Najbauer, University of Pécs Medical School, HUNGARY

**Received:** June 15, 2015

**Accepted:** November 13, 2015

**Published:** January 11, 2016

**Copyright:** © 2016 Parisot et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** Data have been deposited to Figshare: http://dx.doi.org/10.6084/m9.fishare.1550871.

**Funding:** This work was supported by ANRT (Grant 147/2010), Intrasense SAS, the European Research Council Starting Grant Diocles (ERC-STG- 259112). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. The commercial funder Intrasense SAS (www.intrasense.fr) provided support in the form of salaries for authors SP and SC, but did not have any additional role in the study design, data collection and analysis, decision to publish, or

## Abstract

Diffuse WHO grade II gliomas are diffusively infiltrative brain tumors characterized by an unavoidable anaplastic transformation. Their management is strongly dependent on their location in the brain due to interactions with functional regions and potential differences in molecular biology. In this paper, we present the construction of a probabilistic atlas mapping the preferential locations of diffuse WHO grade II gliomas in the brain. This is carried out through a sparse graph whose nodes correspond to clusters of tumors clustered together based on their spatial proximity. The interest of such an atlas is illustrated via two applications. The first one correlates tumor location with the patient's age via a statistical analysis, highlighting the interest of the atlas for studying the origins and behavior of the tumors. The second exploits the fact that the tumors have preferential locations for automatic segmentation. Through a coupled decomposed Markov Random Field model, the atlas guides the segmentation process, and characterizes which preferential location the tumor belongs to and consequently which behavior it could be associated to. Leave-one-out cross validation experiments on a large database highlight the robustness of the graph, and yield promising segmentation results.
preparation of the manuscript. The specific roles of
these authors are articulated in the “author
contributions” section.

**Competing Interests:** Authors SP and SC are affiliated to a commercial company: Intrasense SAS (www.intrasense.fr). This does not alter the authors' adherence to PLOS ONE policies on sharing data and materials.

**Introduction**

Diffuse WHO grade II gliomas (DLGG) are pre-malignant brain tumors characterized by a continuous growth, a migration along the white matter tracts and an unavoidable anaplastic transformation [1]. They affect people in their thirties or forties [2–4]. Two types of tumors are observed depending on their growth speed and the age at the time of diagnosis [5]. An improved knowledge on their natural history has led to shift from a “wait and see” policy to more aggressive therapeutic strategies, with early surgery as the first option, in order to delay malignant transformation while preserving the quality of life of the patients [6, 7].

The management of these tumors is highly dependent on their location. First, several studies have suggested a variability in terms of molecular biology among DLGG with regards to the tumor location [8–11], with a higher rate of 1p deletion in the anterior part of the brain (in particular in the frontal lobe) [8] and a lower rate in the insula [9], or the absence of IDH1 mutation within the insula [10] and its presence for tumors located within the frontal lobe [11]. These differences among DLGGs regarding molecular abnormalities could be one of the explanations to the variability in terms of date of occurrence, growth speed and thus, age at the time of the firsts symptoms or at diagnosis. Second, the dynamical interactions between DLGG and the brain may also vary depending on the eloquence of the areas in which the tumor is located. Indeed, it was shown that slow-growing DLGG might induce cerebral plasticity, explaining the absence of neurological deficit for most patients despite a voluminous tumor, even in the so-called functional regions [12, 13]. Nonetheless, a recent atlas of resectability of DLGG demonstrated that some cerebral areas had low compensatory abilities, constituting a “minimal common brain” among patients [14]. Consequently, the extent of the surgical resection (and thus the median survival) is correlated with the glioma's location, with a better tumor removal in non-eloquent areas rather than in eloquent ones and in compensable regions rather than non compensable ones [15, 16]. Finally, prognosis can differ according to tumor location. In a recently published series on more than 1000 DLGG patients, frontal locations were associated with a better prognosis compared to other locations, and this could not be entirely explained by the possibility of better resection [4].

Interestingly, it has been suggested that DLGG have preferential locations within the brain, with 82.6% of them located within functional areas in a study analysing 132 DLGG [17]. The supplementary motor area and the insular lobe seem to be their most frequent locations [4, 17, 18]. Nonetheless, to the best of our knowledge, there is no probabilistic map of DLGG locations available in the literature.

The aim of this paper is the introduction of a probabilistic atlas computed on a homogeneous series of patients with a DLGG. The idea is to produce a statistical and compact representation of the tumors' preferential locations in the brain to facilitate subsequent location-dependent analyses. We will achieve this using binary maps indicating the positions of the tumors obtained by manual segmentation of a series of MRI images. Those maps are registered in the same reference coordinates. A statistical measure evaluates the tumors' relative position and proximity, and enables to construct a complete graph where each node is a tumor and the arcs' strength corresponds to the proximity measure. We then aim at reducing the graph's size to a handful representative nodes situated in the densest areas (i.e. areas where tumors preferentially appear) using unsupervised clustering with unknown number of populations. The quality of the clustering is evaluated using conventional cluster validation methods and cross validation. The interest with regards to a better understanding of the origins of DLGG and their interaction with the central nervous system will be discussed on the basis of this original atlas, as well as possible clinical applications allowing the optimization of the therapeutic management.
The potential and clinical relevance of this atlas are illustrated through two different applica-
tions. The first one highlights the interest of location specific studies and how they can be facili-
tated by our atlas through a comparison between tumour location and the patient's age. We
show how the cluster configuration of our atlas allows for simplified location specific statistical
analyses and provide preliminary results suggesting a possible link between age and location.

Our second application highlights how this atlas can be integrated in DLGG segmentation
methods and improve their performance by exploiting the knowledge of preferential locations.
Knowing the size and extent of the tumor is of key importance for follow up and therapy plan-
ning. Currently, the size of a DLGG is approximated by using manual segmentation, a time
consuming and tedious task subject to inter expert variability. Automatic tumor segmentation
is a difficult task, due to the tumors' heterogeneous shapes and appearances, overlapping inten-
sities with the surrounding healthy tissue and fuzzy boundaries.

A popular approach for their segmentation is to rely on learning statistical classifiers to sep-
arate tumor voxels from healthy voxels. The Support Vector Machines (SVM) classifiers have
been extensively studied for this purpose [19–21]. However, this kind of methods have encoun-
tered limited success due to the underlying assumption that all voxels are independent and
identically distributed (i.i.d), meaning that all voxels are segmented independently from their
neighbors. Neighborhood dependent features or morphological filters have been considered to
introduce spatial information with limited success. More recently, methods have modeled the
spatial dependencies via the use of random fields [22–25]. In this setting, a learned classifier is
coupled with a regularization term that penalizes segmentation discontinuities on a defined
local neighborhood. The main drawback of those methods is that they encode the spatial
dependencies in a local manner, lacking global information on the brain structures’ boundaries.
Gering et al. [26] adopted a multilayer MRF approach where an intensity based voxel wise seg-
mentation is progressively refined by incorporating high-level neighboring information such
as the distance between the different brain structures. The tumor is detected as a deviation to
healthy voxels intensity. User interaction is required for initialization. Corso et al. [27] intro-
duce more global information by coupling Gaussian Mixture Model classification with a multi
level graph structure, where the edges of the graph have an affinity that characterizes the simi-
larity between the neighboring nodes.

Atlas-based segmentation methods are endowed with global properties, based on a healthy
brain's expected structures. Kaus et al. [28] alternate statistical classification based on intensity
difference with registration of the data with an anatomical atlas, where the tumor voxels are
reclassified as healthy. This assumes strong homogeneity in the tumor appearance. Prastawa
et al., Moon et al. [29, 30] use a registered probabilistic atlas in which probabilities for tumor
and edema are encoded as prior information, based on contrast enhancement. Spatial and geo-
metric constraints are added to avoid false detections. In [31, 32] the tumors are detected as
outliers with respect to the normal brain tissue characteristics of a registered healthy atlas. Spa-
tial constraints are modeled via Markov Random Fields [31] or level sets [32]. In [33], we pro-
posed a MRF based coupled tumor segmentation and registration framework where atlas based
prior knowledge is introduced through the registration task.

Our method combines local and global informations. Tumor detection is based on a learned
classifier and spatial constraints are enforced using a Markov Random Field model, while the
statistical atlas of uneven repartition in the brain adds global prior information on the most
probable location of tumor voxels. The method is tested on a challenging clinical database with
variable quality and poor resolution, and yields promising results that demonstrate the clinical
potential of the atlas.

The remainder of this paper is organized as follows: section 2 describes the methodology
towards mapping DLGG's preferential locations in the brain and its application for tumor
segmentation. Experimental validation and obtained results are part of section 3 while discussion and perspectives conclude the paper.

# Materials and Methods

## Database Construction and Preprocessing

We conducted a retrospective study at a single institution on a series of 210 DLGG patients. This study was reviewed and approved by the Institutional Review Board (Institut Regional du Cancer de Montpellier - Val d'Aurelle, ID number ICM-URC-2015/35). The patient information was anonymized and de-identified prior to analysis. Patients aged over 18 at diagnosis were included if they fulfilled the following criteria: surgery for a DLGG performed by one of the 3 neurosurgeons of the neuro-oncology unit (HD, LB or YY) between January 1st, 2006 and August 30th, 2012; FLAIR weighted MRI images performed before any oncological treatment (chemotherapy or radiation therapy) available and exploitable. Exclusion criteria were as follows: no histological diagnosis available; low-grade gliomatosis cerebri; tumor volume over 150 cm³ or multicentric; DLGG located in the spinal cord or in the brainstem; delay >3 years between the first MRI performed (at diagnosis) and the first MRI available. The age at diagnosis (i.e. at the time of the first MRI performed), the age at the time of the first symptoms and at the time of the first MRI available were collected for each patient. For incidental tumors (n = 16) the age at the time of the first symptoms was considered to be the age at the time of the first MRI.

The MRI performed at diagnosis was collected for each patient. When that MRI was not available, a later MRI could be used if the patient had not yet received any oncological treatment at the time when it was performed (biopsy possible). The OSIRIX® software was used to perform the image manual segmentation on the FLAIR weighted sequences for each patient. Segmentations were performed by HD and AD. A region of interest (ROI) was defined for each patient, determining for each voxel of each slide of the FLAIR sequence if it belongs to the ROI (i.e. to the tumor) or not. The complete database and expert segmentations are available online at http://db-gliomas-gradeii.net/ and on figshare at http://dx.doi.org/10.6084/m9.fishshare.1550871.

Due to our patient selection constraints for atlas construction, all images are acquired in a clinical setting and have therefore variable quality (often low resolution) and progression of the disease. This makes all image processing tasks (tumor segmentation, deformable registration) particularly challenging.

Preprocessing involves intensity regularization for the segmentation task. One of the main issues with MRI imaging is that the intensities of the same anatomical tissues can be very different from one image to another, and even within the same image. Many algorithms have been developed to correct the intensity inhomogeneity [34]. We adopt a simple regularization process:

$$I_{\text{reg}} = \frac{I - \text{Med}(I)}{\text{IQR}(I)} \qquad (1)$$

where Med(I) is the median intensity of the image and IQR(I) is the interquartile range. Those values are computed without taking into account the background pixels.

## Probabilistic Atlas Construction

A High anatomical variability exists between the brains of different individuals. In order to compare the tumors' positions in the brain, we need to compensate such variability and bring all the MRI images in the same reference coordinate space. To this end, we perform affine registration towards the same reference pose on all the volumes. While deformable registration
would yield better correspondences, the presence of the tumor alters the quality of the registration. The most straightforward approach to deformable registration in such a setting, which we used in a previous work on a much smaller database [35], is to mask the pathology while performing the registration. It has been shown that it can lead to important errors in registration [33] and therefore cannot be used reliably on a large database. Existing methods developed specifically for registration in the presence of tumours [33, 36] currently lack reliability and scalability to be used systematically on our large database, notably due to the variable quality of our data. Let us consider N images $I_i(\mathbf{x})$, $i \in [1, \dots, N]$ featuring a DLGG, and the healthy reference pose $A(\mathbf{x})$. All images are affinely registered to the reference pose. We call $M_i(\mathbf{x})$ the corresponding binary tumor segmentation map obtained by manual segmentation, after affine registration.

Registration brings all the tumors (segmentation maps) in the same coordinate system. The next step is to build a graph expressing the statistical relationships between the tumors' positions. Using an arbitrary distance measure, we can build a complete graph $G$ where the nodes correspond to the tumors and the edges' strength corresponds to the distance between the tumors. The nodes' probability of being at the center of a preferential location increases with the amount of tumors being at proximity (node centrality in the network). As a distance measure, we adopt the euclidean distance between the centers of gravity of the tumors, computed as the mean coordinate among all tumor voxels. This measure has the advantage of being independent from the size of the tumor and of being more robust to registration errors than the use of surface based distances. This strongly reduces the impact of using affine registration over deformable registration. We will now consider $\mathbf{d} = d(M_i, M_j)$, the set of distances between all possible pairs of tumors.

The complete graph obtained with this distance measure is shown in Fig 1, where the color and thickness of the edges illustrate the distances values (from thick and red (low distance) to thin and blue) while the size of the nodes increases with their centrality.

This graph already suggests some areas with a high density of tumors (red edges) and others where we seldom find tumors. To confirm those observations and identify preferential areas, we aim at reducing the dimensionality of the graph by regrouping close tumors in clusters. If there exists indeed preferential locations, we should obtain compact and well separated clusters. To perform such a clustering, we need to identify (i) the number of clusters $K$, (ii) the centers of the clusters $\bar{c}_1, \dots, \bar{c}_K$, (iii) the remaining nodes' assignments $l_i, i \in [1, \dots, N]$ to the different clusters $c_1, \dots, c_K$. Popular clustering algorithms such as K-means have two main drawbacks: they are very sensitive to initialization and can easily get trapped in a local minimum, while the number of clusters needs to be predefined. To tackle these issues, we make use of a recently proposed unsupervised clustering method [37] that is able to overcome the aforementioned limitations of conventional methods like K-means. This method automatically determines the optimal number of clusters and is independent from initialization.

$$ \min_{k} \min_{c_1, \dots, c_K} \min_{\bar{c}_1, \dots, \bar{c}_K} \left( \sum_{i=1}^{N} \sum_{j=1}^{K} \delta(l_i, c_j) d(\bar{c}_j, M_i) + \alpha \sum_{i=1}^{K} \mu(\bar{c}_i) \right) \quad (2) $$

Where $\delta(\cdot, \cdot)$ is the Kronecker delta function. First, the nodes assignments should be determined so as to minimize their distance to the closest center (first term). The second term is a penalty introduced in order to avoid the trivial solution of selecting all the nodes as clusters. It is defined as the mean distance value between the considered node and all the remaining nodes:

$$ \mu(M_i) = \frac{1}{N} \sum_{j=1}^{N} d(M_i, M_j) \quad (3) $$
<figure><img src="image_5.png" /><figcaption>**Fig 1. Network representation of the whole data-set before clustering superimposed to the mean registered image.** The nodes are located in the center of gravity of the tumors and bigger nodes have bigger network centrality (i.e. they are strongly connected to many other nodes). The edges colors and strength represent the distance between nodes (from red and thick (short distance) to blue and thin). For visibility reasons, arcs corresponding to a distance greater than 35 are not displayed.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g001

This term enables the nodes with the highest centrality (i.e with the highest amount of close
tumors) to be selected as centers. We control the relative importance of the 2 terms using the α
coefficient. For each α value, we evaluate the quality of the clustering using the 3 following
indexes:

The Dunn index [38] is very commonly used for cluster evaluation. It compares the biggest
distance intra cluster (cluster “diameter”) to the smallest distance inter cluster. The formulation
proposed in [38] is very sensitive to noise and outliers as the distances between all nodes are
computed. Here, we follow the more robust graph theoretic approach of [39] to define the
diameter and distance inter cluster. Each cluster is represented as a Minimum Spanning Tree
(MST). Given the set of nodes $V_j$ that belong to the cluster $c_j$, a spanning tree is the minimal
subgraph to the complete graph that connects all nodes. The weight of the spanning tree corre-
sponds to the sum of the weights of all edges in the tree. The MST $G_j^{MST} = (V_j, \mathcal{E}_j^{MST})$ is the
spanning tree of minimal weight. The cluster diameter is then defined as the maximum edge
weight among all edges in the MST graph.

$$D = \min_{i \in [1:K]} \left\{ \min_{j \in [1:K], j \neq i} \left\{ \frac{d(\bar{c}_i, \bar{c}_j)}{\max_{k,k' \in \varepsilon_j^{\text{MST}}} d(M_k, M_{k'})} \right\} \right\} \quad (4)$$

A high Dunn index corresponds to compact and well separated clusters (high $d_{inter}$ and low $d_{intra}$), i.e. to a good clustering. As a result, we aim at finding a clustering that maximizes the Dunn index.

The Davies-Bouldin [40] index defines a measure of similarity between clusters:

$$R_{i,j} = \frac{\sigma_i + \sigma_j}{d(\bar{c}_i, \bar{c}_j)} \quad (5)$$

and then computes the maximum similarity for each cluster:

$$DB = \frac{1}{k} \sum_{i=1}^{k} \max_{j \in [1:k], j \neq i} R_{i,j} \quad (6)$$

where $\sigma_i$ is the average distance of all points in cluster $c_i$ to its center. A good clustering corresponds to a low DB index (little similarities between different clusters). Like the Dunn index, DB identifies compact and well separated clusters.

Last, we compute the Silhouette index [41]. It computes for each node a score of confidence with respect to its cluster assignment.

$$s(M_i) = \frac{b(M_i) - a(M_i))}{\max(a(M_i), b(M_i))} \quad (7)$$

$a(M_i)$ is the average distance between $M_i$ and all the remaining elements assigned to the same cluster while $b(M_i)$ is the average distance between $M_i$ and all the elements in the closest cluster. The Silhouette index takes values between -1 and 1. If the value is close to 1, the node has been assigned correctly. A value close to zero suggests that the node is equally far away from 2 clusters, while a silhouette index close to -1 infers that the node has been misclassified. To evaluate the quality of the clustering, we compute the global silhouette index:

$$GS = \frac{1}{K} \sum_{j=1}^{K} \frac{1}{n_j} \sum_{i=1}^{n_j} s(M_i) \quad (8)$$

The highest GS corresponds to the best clustering, where the individual Silhouette indexes are closest to 1.

We select the clustering result (i.e. $\alpha$ value) that corresponds to the best indexes value.

## Atlas Application: Tumor Segmentation and Characterization with Position Priors

The probabilistic atlas constructed through clustering provides interesting insight on the locations where tumors are likely to appear and on the tumors' location dependent properties. Considering a new image I featuring a tumor, assigning it to a cluster could enable to predict its future behavior and spatial repartition based on the cluster properties as well as provide powerful spatial prior information for automatic segmentation.

In this section, we consider that the atlas has been constructed as a K clusters graph. For each cluster $c_i$, we now consider the spatial extension of all tumors within this cluster through their associated binary segmentation maps. We build a distribution map as $P(x|c_i)$, $i \in [1:K]$
<figure><img src="image_6.png" /><figcaption>**Fig 2.** Examples of cluster probability maps describing the spatial repartition of tumors in the cluster. The maps are superimposed to the mean registered image.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g002

describing the conditional probability of tumors appearances depending on their location (with respect to the reference frame coordinate system). This can be done simply by counting the number of tumors in the corresponding cluster appearing at each voxel voxel. Considering there are N tumors in cluster cᵢ, we defined P(x|cᵢ) as:

$$P(\mathbf{x}|c_i) = \frac{1}{N} \sum_{j \in \{1:N\}} M_j(\mathbf{x}) \qquad (9)$$

To compensate for the possible small amount of tumors per cluster, all probability maps are smoothed using a Gaussian filter. Examples of such probability maps are shown in Fig 2.

We formulate the problem of segmentation and characterization as a labeling problem where we seek to assign a characterization label $\theta_x \in C = \{c_1, ..., c_k\}$ and a segmentation label $\omega_x \in S = \{0, 1\}$ to each voxel of the image. The characterization label's role is to determine
which preferential location the tumor belongs to by assigning it to one of the sparse graph's clusters while the segmentation label aims at separating tumor ($\omega_x = 1$) and background voxels ($\omega_x = 0$). We adopt a Markov Random Field (MRF) model on a graph $G = (V, E)$, where the graph nodes $V$ correspond to the image voxels $x$ and the graph's edges $E$ define a first order neighborhood system $N$ by connecting a node to its 6 immediate neighbors. The optimal labeling $l = [\theta, \omega]$ is recovered by minimizing the MRF energy:

$$E(1) = \sum_{x} V_{x}(l_{x}) + \sum_{x} \sum_{y \in N(x)} V_{x,y}(l_{x}, l_{y}) \qquad (10)$$

$V_x(l_x)$ is the unary cost and describes the likelihood of assigning a specific label to voxel $x$, while $V_x, y(l_x, l_y)$ describes the pairwise interactions between two neighbouring voxels. The unary and pairwise costs will be details in the subsequent subsections.

**Detection.** Our first step towards tumor segmentation and characterization is the construction of a classifier that will distinguish tumorous and healthy pixels. We learn such a classifier using the boosting algorithm Gentle Adaboost [42] and a set of features extracted from the image. The training set is made of an ensemble of pairs $\{x_i, y_i\} \in \mathbb{R}^3 \times \{-1, 1\}, i \in [1:N], x_i$ is an observed voxel from a training volume and $y_i$ is its corresponding label. The idea behind boosting is that a strong classifier can be created from an ensemble of weak classifiers. Each pair is associated to a feature vector $f(x_i)$ and a weight $D_i = 1/N$. At each iteration t, a weak classifier $h_t(x_i)$ is constructed as a decision stump in order to minimize the training error $\sum_{i=1}^N D_i(h_t(x_i) - y_i)^2$. The weights are then updated as $D_i = D_i \exp(-y_i h_t(x_i))$, putting the emphasis on misclassified sample. Eventually, the strong classifier is computed as a linear combination of the weak classifiers, yielding a score of confidence $H(I(x))$ for each location of the image to be segmented given the intensity map.

The features used for classification rely on intensity, texture and symmetry. The main feature differentiating healthy and tumorous tissues is the intensity enhancement of tumor pixels on FLAIR images. We exploit this characteristic by using $9 \times 9 \times 5$ intensity patches centered on the pixel of interest. The underlying idea is also to get information about the neighborhood of the pixels since there exists overlapping intensities between tumor and background. Rotation invariant intensity statistics are also used, as we computed the median, standard deviation and entropy of intensity patches of size $k \times k \times 3$, where $k = [3, 5, 7]$, all centered on the pixel of interest. Gabor filters [43] are particularly useful for distinguishing objects of different textures. The filters are used on 2 scales and 10 orientations in order to detect the main structures of the brain (such as the skull). The brain's hemispheres have the interesting property of being roughly symmetrical. This can be taken advantage of in our case, since the tumor introduces dissymmetry between the 2 hemispheres. Thanks to affine registration, the symmetry plane of all images is roughly equivalent to the approximated reference pose's $\Pi$. For each pixel $x$ in the left half of the image, we compute a symmetry measure as follows:

$$S(\mathbf{x}) = \frac{1}{N} \sum_{\mathcal{N}_s(\mathbf{x})} I(\mathbf{x}) - \frac{1}{N} \sum_{\mathcal{N}_s(\mathbf{x})} I(\mathbf{x}_{\Pi}) \qquad (11)$$
$$S(\mathbf{x}_{\Pi}) = -S(\mathbf{x})$$

where $\mathbf{x}_{\Pi}$ is the symmetric of $\mathbf{x}$ with respect to $\Pi$, $\mathcal{N}_s(\mathbf{x})$ is a neighborhood of $\mathbf{x}$ and $N$ the number of pixels in $\mathcal{N}_s(\mathbf{x})$. We use a neighborhood to compensate the fact that the symmetry plane is approximate.

**Pairwise Potentials.** The MRF energy's pairwise potential acts as smoothing priors, imposing spatial consistency of the labellings. We adopt a classic Potts model, penalizing
neighboring nodes that are assigned different labels:

$$V_{x,y}([\theta_x, \omega_x], [\theta_y, \omega_y]) = \begin{cases} 0, & \text{if } \theta_x = \theta_y \text{ and } \omega_x = \omega_y \\ \beta, & \text{if } \theta_x = \theta_y \text{ and } \omega_x \neq \omega_y \\ \infty, & \text{otherwise} \end{cases} \qquad (12)$$

where $\beta$ is a constant parameter imposing the strength of the penalization. Based on the hypothesis that there is only one tumor per patient, the characterization label has to be the same for all grid nodes. Setting the penalization to infinity imposes the same characterization label on all nodes.

**Unary Potentials.** The unary potential represents the nodes likelihoods and seeks the most likely label configuration. It is composed of three terms:

$$V_x([\theta_x, \omega_x]) = V_{S,x}(\omega_x) + V_{C,x}(\theta_x) + V_{C,S,x}(\theta_x, \omega_x) \qquad (13)$$

The segmentation unary potential $V_{S,x}(\cdot)$ seeks the most likely segmentation configuration based on the classification likelihoods recovered from the boosting classification decisions. We follow the approach of [44] and define it as:

$$V_{S,x}(\omega_x) = -w_x \log P_{tm}(I(\mathbf{x})) - (1 - w_x) \log P_{bg}(I(\mathbf{x})) \qquad (14)$$

where the prior probabilities are computed from the boosting classification score as:

$$P_{tm}(I(\mathbf{x})) = (1 + \exp(-2H(I(\mathbf{x}))))^{-1} \qquad (15)$$

$$P_{bg}(I(\mathbf{x})) = 1 - P_{tm}(I(\mathbf{x}))$$

The role of the characterization term $V_{C,x}(\cdot)$ is to identify the closest cluster the tumor to be segmented belongs to. In order to identify the most likely cluster, we measure the degree of non overlap of pixels with a high boosting classification score $H(I(\mathbf{x}))$ with the center of each cluster.

$$V_{C,x}(\theta_x) = \max(0, H_{bin}(I(\mathbf{x})) - M_{\theta_x}(\mathbf{x})) \qquad (16)$$

where $M_{\theta_x}(\cdot)$ is the binary map associated to the center of the corresponding cluster and $H_{bin}(\cdot)$ is the boosting score converted into a thresholded (with threshold T) binary map:

$$H_{bin}(I(\mathbf{x})) = \begin{cases} 1, & \text{if } H(I(\mathbf{x})) > T \\ 0, & \text{otherwise} \end{cases} \qquad (17)$$

Basically, $V_{C,x}(\cdot)$ counts the number of pixels with a boosting score beyond a given value that do not overlap with the center of the cluster. The value of $\theta_x$ that minimizes $V_{C,x}$ is the closest cluster (highest overlap with the cluster center).

The coupling term $V_{C, s, x}(\cdot, \cdot)$, linking the segmentation and characterization, is the key element of our framework. The selected cluster gives information on where the tumor voxels are expected to appear and not to appear, according to the spatial distribution of tumors in the cluster. This term plays the part of a spatial position prior and is defined as:

$$V_{S,C,x}(\theta_x, \omega_x) = -w_x \log P(\mathbf{x}|\theta_x) - (1 - w_x) \log (1 - P(\mathbf{x}|\theta_x)) \qquad (18)$$

This term ensures that the segmentation labels are consistent with the location corresponding to the assigned type of tumor: tumor detections that are not in accordance with the expected spatial position are penalized.
The MRF optimization is done using Fast PD [45], an optimization method based on linear programming, which offers a great compromise between speed and performance.

## Second Application: Correlation with Patient's Age: Statistical Analysis

Patients' characteristics, including age at the time of the first symptoms and age at MRI diagnosis (i.e. age at the time of the first MRI), were described using percentages for categorical variables, and mean, standard deviation, median and range for continuous variables. Ages in each cluster (at the time of the first symptoms and at MRI diagnosis) were described with boxplots. An analysis of variance (ANOVA) and a Bonferroni correction for multiple tests were used to compare the mean ages of each cluster. Statistical analysis was performed using the SAS (Statistical Application System) v9.3 software.

## Results

Two hundred and ten patients were included, among them 97 women (46.2%) and 113 men (53.8%). The median age at the time of the first symptoms, at the time of the first MRI and at the time of the MRI analyzed were 37 years (range 15-70), 37 years (range 17-70) and 37 years (range 17-70), respectively. The mean age at the time of the first symptoms, at the time of the first MRI and at the time of the MRI analyzed were 38 years. A surgical resection of the tumor was performed for 202 patients (96.2%) and a biopsy for 8 of them (3.8%). The median and the mean tumor volumes were 44.7 cm³ (range 0.1-147.9 cm³) and 52.5 cm³, respectively. The median and the mean delays between the first MRI performed (i.e. at diagnosis) and the MRI analyzed were 37 days (range 0-2.8 years) and 115 days, respectively. The initial symptoms were epilepsy in a majority of patients (86.2%), a neurological deficit (speech disturbance, visual deficit or cognitive impairment) in 2.4%, and were unknown in 1.9%. The diagnosis of DLGG was incidental for 9.5% of patients (the brain imaging was performed for another disease or for symptoms that cannot be related to the DLGG).

The image size varied from 144 to 576 in the (x, y) plane and 12 to 213 in the z plane, and the voxel resolution from 0.4 × 0.4 to 1 × 1 mm² in the (x, y) plane and 0.9 to 10 mm in the z plane. The tumor size ranged from 0.3 to 180 cm³. The reference pose used for registration is a 256 × 256 × 24 FLAIR image of a healthy brain, with resolution 0.9 × 0.9 × 5.5.

## Sparse Graph Validation

In order to select the best clustered graph, we performed clustering of the complete graph for $\alpha$ values ranging from 0.3 to 7 producing clustering results involving 5 to 49 clusters. We computed the 3 indexes for all the values $\alpha$. Combining the indexes after normalization as $I_c = (D, (1 - DB).GS)^{1/3}$, we observed a global maximum for $\alpha = 2.1$, corresponding to a 11 clusters graph. Fig 3 shows the different indexes values and their combination, while Fig 4 shows the 11 clusters graph (i.e. the best clustering). We can observe some symmetry in the graph topology. Five out of the 11 clusters are located within the left hemisphere (clusters 1, 4, 5, 6 and 11) while six are located within the right hemisphere (clusters 2, 3, 7, 8, 9 and 10). Clusters 1 and 2 are located in the parietal lobe (left and right parietal lobes, respectively). Clusters 4 and 7 are located within the left and right temporal lobes, respectively, while clusters 8, 9 and 11 are located within the insula. The remaining clusters involve the frontal lobes, with clusters 6 and 10 located in the left and right anterior frontal lobes, respectively, and cluster 3 located within the right prefrontal lobe. The anatomical repartition of clusters is summarized in Table 1 and the positions of the centers of each cluster in the MNI atlas [46, 47] are shown in Fig 5. Positions in the MNI atlas are obtained by non-rigid registration [48] of the FLAIR reference pose to the MNI atlas.
<figure><img src="image_7.png" /><figcaption>**Fig 3.** Cluster validity indices with respect to the value of α (a, b, c, d) and the number of clusters (e). Dunn index (a), Davies Bouldin index (b), Silhouette index (c) and combined indices (d, e).</figcaption></figure>

doi:10.1371/journal.pone.0144200.g003

Next, we performed Leave-One-Out Cross Validation to assess the quality of our clustering. For 210 experiments (one per sample), we removed one sample and learned the optimal graph's topology on the remaining 209 nodes. We then studied the impact of the 210th sample on the structure of the graph and the quality of the clusters. For each of the 210 experiments, we computed the 3 indexes. In most cases (207 samples, 98.5% of the data), the graph was a 11 nodes graph which corresponds to the graph obtained with the full data-set. For 3 samples, the optimal number of clusters was 10, which was due to the fact that they belonged to a loose cluster (center sample or close to center sample) that was merged with its closest neighbor when the sample was removed.

We then studied the robustness of the graph by analyzing the influence of each sample. For each clustering obtained by cross validation, we evaluated the impact of the 210th sample on the structure of the graph and the quality of the clusters. To evaluate how well the clustering represented the data, we tried to assign each removed sample to a cluster. We used three criteria to do so:

1. *The distance from the center of the cluster.* The sample has to be assigned to the closest cluster and should be as close to the corresponding center as the nodes constituting the cluster.
<figure><img src="image_8.png" /><figcaption>**Fig 4.** Visualization of the complete clustered graph superimposed to the mean registered image. The numbers correspond to the number of nodes in each cluster.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g004

<table><caption>**Table 1.** Anatomical location of the different clusters.</caption><thead><tr><th>Anatomical location</th><th>Number of clusters</th><th>Percentage of data</th></tr></thead><tbody><tr><td>Parietal lobe</td><td>2 (symmetric)</td><td>12%</td></tr><tr><td>Supplementary Motor Area (frontal lobe)</td><td>1</td><td>9%</td></tr><tr><td>Temporal lobe</td><td>2 (symmetric)</td><td>18%</td></tr><tr><td>Frontal lobe</td><td>1</td><td>9%</td></tr><tr><td>Anterior frontal (frontal lobe)</td><td>2 (symmetric)</td><td>15%</td></tr><tr><td>Insula</td><td>3 (2 symmetric clusters)</td><td>37%</td></tr></tbody></table>

doi:10.1371/journal.pone.0144200.t001

To determine if the distance to the cluster center follows the cluster's distance distribution or if it is an outlier we use rely on the distribution's quartiles:

$$d(M_i, \bar{c}_k) = \min_{j \in [1:K]} d(M_i, \bar{c}_j) \text{ and } d(M_i, \bar{c}_k) \le Q_3(c_k) + 1.5(IQR(c_k)) \quad (19)$$

where $Q_3(c_k)$ and $IQR(c_k)$ are respectively the third quartile and interquartile range of cluster $c_k$'s distance distribution.
<figure><img src="image_9.png" /><figcaption>**Fig 5.** Positions of the cluster centers on the MNI atlas. The clusters are organized in numerical order (from cluster 1 to 11).</figcaption></figure>

doi:10.1371/journal.pone.0144200.g005
<figure><img src="image_10.png" /><figcaption>**Fig 6. Examples of graph matching results.** (a) Complete match, (b, c) Partial match. Positive matches correspond to blue edges and mismatched samples to red edges. The nodes' locations correspond to the coordinates of the cluster centres of each clustering. The green nodes have been translated along the x axis for visualisation purposes.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g006

2. The individual Silhouette index. The global silhouette index is indicative of the quality of the clustering. The individual index tells how well the sample fits in its cluster. Considering the sample belongs to the closest cluster, we compute the corresponding silhouette index. If the index is close to 1, there is no doubt about the fact that the sample belongs to this cluster. A value close to zero does not necessarily mean that the clustering is bad. Indeed, the main drawback of the silhouette is that it can penalize large clusters (large diameter) when compared to smaller ones. To avoid penalizing large clusters, we also compare the individual index to the distribution of silhouette values:

$$s_k(M_i) = \max_{j \in [1:K]} s_j(M_i) \text{ and } s_k(M_i) \ge Q_1(c_k) - 1.5(IQR(c_k)) \quad (20)$$

where $Q_1(c_k)$ and $IQR(c_k)$ are respectively the third quartile and interquartile range of cluster $c_k$'s silhouette index distribution.

3. The last criterion verifies the consistency between the clusters selected by the first two criteria:

$$\arg\min_k d(M_i, \bar{c}_k) = \arg\max_k s_k(M_i) \quad (21)$$

Using those criteria, we manage to classify 88% of the whole data-set without affecting the quality of the cluster (mean distance intra cluster, mean silhouette index).

We also studied the topology correspondences of the graphs with respect to the original graph $G_0$. This was only possible for the 207 samples that yielded a 11 clusters graph. Using the algorithm presented in [49], we seek a matching between the nodes and arcs' strength of the respective graphs. We observed a complete match for 70% of the graphs, and a one node difference for 26%. The remaining graphs differed from 2 nodes (9 cases) and 4 nodes (1 case). In the case of partial matching, the graphs were still very similar as the maximum distance between the different cluster centers was 12.6 which corresponds to the mean intracluster distance of tight clusters. Fig 6 illustrates partial and complete matches, the worst matching case (4 different samples) is illustrated in Fig 6c. We can see that the graphs are still very similar in the case of partial match.
<table><caption>**Table 2.** Mean and median (in parentheses) values of the different scores and methods for tumor segmentation.</caption><thead><tr><th>Method</th><th>Dice</th><th>True Positive</th><th>False Positive</th><th>MAD</th></tr></thead><tbody><tr><td>Boosting</td><td>50 (55)%</td><td>77 (83)%</td><td>59 (58)%</td><td>22 (22) mm</td></tr><tr><td>Boosting regularized</td><td>66 (74)%</td><td>70 (78)%</td><td>29 (21)%</td><td>9.5 (6) mm</td></tr><tr><td>Our method</td><td>69 (77)%</td><td>66 (74)%</td><td>20 (12)%</td><td>6 (3.5) mm</td></tr></tbody></table>

doi:10.1371/journal.pone.0144200.t002

## Segmentation results

The segmentation was evaluated on the 210 volumes through Leave One Out Cross Validation. For each sample, an atlas is constructed using the data's optimal indexes value and a boosting classifier is learned by randomly selecting 35 volumes among the 209 training volumes. Experiments showed that the potential assigning a cluster to the image ($V_{clus}$) had to be given a much bigger importance than the others in order to properly label the whole training set to the right cluster (especially if the tumor is small). T was set at 1.5, $\lambda$ at 10000, $\gamma$ at 2 and $\beta$ at 1.

We evaluate the quality of our segmentation framework by using the manual segmentations as ground truth. Let us define M the ensemble of pixels labeled as tumor manually, and A the pixels automatically labeled. For each image, we compute the Dice coefficient, False (FP) and True (TP) Positive rates and the Mean Absolute Distance between contours (MAD):

$$ 
\begin{align*}
D &= \frac{2 \| A \cap M \|}{\| A \cup M \| + \| A \cap M \|} \\
TP &= \frac{\| A \cap M \|}{\| M \|} \\
FP &= \frac{\| A \| - \| A \cap M \|}{\| A \|} \\
MAD &= \frac{1}{2} \left( \frac{1}{n_A} \sum_{x \in A} d_{min}(x, M) + \frac{1}{n_M} \sum_{x \in M} d_{min}(A, x) \right)
\end{align*}
\qquad (22) $$

where $n_A$ and $n_M$ are the number of voxels in A and M respectively and $d_{min}(x, M) = \inf(d(x, y)|y \in M)$ is the smallest distance of a point x to any point in M. The distance $d(.,..)$ used here is the euclidean distance.

Results are shown using boosting classification, boosting with pairwise regularization, and boosting with pairwise regularization and spatial prior. Table 2 shows the mean and median scores for the different configurations. Boxplots of the values are shown in Fig 7, while visual examples of segmentation results are illustrated in Fig 8. We observe an increase of the Dice score (+ 3%) as well as a strong decrease of the false positive rate (- 9%) with respect to the regularized method. The MAD score provides the most interesting evaluation since it is not biased by the size of the tumor (contrarily to the other measures). We can observe a strong improvement of the MAD score using the probabilistic atlas. A decrease of the true positive rate is observed, that is mostly associated to small and poorly detected tumors. In that case, the selected cluster may be inadequate and cause a decrease in the quality of the segmentation.

## Correlation between Tumor Topography and Patient's Age

The ages at the time of the first symptoms and at the time of MRI diagnosis are described for each of the 11 clusters in Tables 3 and 4 and Fig 9 (box-plots). The analysis of variance found no significant difference in terms of age at the time of the first symptoms or age at the time of
Dice Score

<figure><img src="image_11.png" /></figure>

(a)

True Positive Rate

<figure><img src="image_12.png" /></figure>

(b)

False Positive Rate

<figure><img src="image_13.png" /></figure>

(c)

Mean Absolute Distance

<figure><img src="image_14.png" /></figure>

(d)

Fig 7. Boxplots of the Dice score (a), True Positive rate (b), False positive rate (c), and MAD score (d) between the automatic and manual tumor segmentation for the three different methods.

doi:10.1371/journal.pone.0144200.g007
<figure><img src="image_15.png" /><figcaption>**Fig 8.** Visual Segmentation results. (a) boosting score, (b) Boosting classification (thresholding), (c) Pairwise MRF, (d) MRF with spatial prior.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g008
<table><caption>Table 3. Mean and median ages at the time of the first symptoms for each cluster.</caption><thead><tr><th rowspan="2">Cluster</th><th rowspan="2">Nsamples (%)</th><th colspan="5">Age at the time of the first symptoms (in years)</th></tr><tr><th>Mean</th><th>Standard Deviation</th><th>Median</th><th>Min</th><th>Max</th></tr></thead><tbody><tr><td>1</td><td>9 (4.3)</td><td>39.1</td><td>12.2</td><td>40.4</td><td>20.5</td><td>61.4</td></tr><tr><td>2</td><td>17 (8.1)</td><td>39.3</td><td>10.7</td><td>39.8</td><td>23.1</td><td>62.8</td></tr><tr><td>3</td><td>18 (8.6)</td><td>34.2</td><td>8.6</td><td>33.0</td><td>22.3</td><td>52.4</td></tr><tr><td>4</td><td>13 (6.2)</td><td>39.4</td><td>7.5</td><td>38.6</td><td>28.6</td><td>52.7</td></tr><tr><td>5</td><td>19 (9.0)</td><td>37.2</td><td>10.3</td><td>37.3</td><td>15.9</td><td>56.5</td></tr><tr><td>6</td><td>12 (5.7)</td><td>34.7</td><td>10.8</td><td>34.5</td><td>21.0</td><td>51.7</td></tr><tr><td>7</td><td>25 (11.9)</td><td>40.7</td><td>10.0</td><td>41.2</td><td>20.3</td><td>63.2</td></tr><tr><td>8</td><td>19 (9.0)</td><td>35.6</td><td>10.5</td><td>36.6</td><td>18.5</td><td>53.9</td></tr><tr><td>9</td><td>30 (14.3)</td><td>39.0</td><td>12.6</td><td>35.8</td><td>22.0</td><td>70.0</td></tr><tr><td>10</td><td>19 (9.0)</td><td>38.6</td><td>11.2</td><td>37.9</td><td>17.7</td><td>59.1</td></tr><tr><td>11</td><td>29 (13.8)</td><td>39.1</td><td>12.1</td><td>37.3</td><td>18.2</td><td>67.1</td></tr></tbody></table>

doi:10.1371/journal.pone.0144200.t003

<table><caption>Table 4. Mean and median ages at the time of MRI diagnosis for each cluster.</caption><thead><tr><th rowspan="2">Cluster</th><th colspan="5">Age at the time of MRI diagnosis (in years)</th></tr><tr><th>Mean</th><th>SD</th><th>Median</th><th>Min</th><th>Max</th></tr></thead><tbody><tr><td>1</td><td>39.5</td><td>12.5</td><td>40.4</td><td>20.7</td><td>61.4</td></tr><tr><td>2</td><td>39.3</td><td>10.7</td><td>39.8</td><td>23.1</td><td>62.8</td></tr><tr><td>3</td><td>34.4</td><td>8.8</td><td>33.1</td><td>22.3</td><td>52.7</td></tr><tr><td>4</td><td>39.7</td><td>7.4</td><td>38.6</td><td>28.9</td><td>52.9</td></tr><tr><td>5</td><td>37.4</td><td>10.2</td><td>37.3</td><td>17.3</td><td>56.6</td></tr><tr><td>6</td><td>36.2</td><td>10.7</td><td>35.1</td><td>21.0</td><td>56.4</td></tr><tr><td>7</td><td>41.0</td><td>10.2</td><td>41.2</td><td>20.3</td><td>63.2</td></tr><tr><td>8</td><td>35.8</td><td>10.6</td><td>36.7</td><td>18.5</td><td>54.0</td></tr><tr><td>9</td><td>39.6</td><td>12.7</td><td>37.2</td><td>22.0</td><td>70.3</td></tr><tr><td>10</td><td>38.6</td><td>11.2</td><td>38.0</td><td>17.8</td><td>59.1</td></tr><tr><td>11</td><td>40.2</td><td>12.2</td><td>37.4</td><td>18.2</td><td>67.1</td></tr></tbody></table>

doi:10.1371/journal.pone.0144200.t004

MRI diagnosis (p = 0.73 and p = 0.72, respectively). Nevertheless, patients whose MRI belongs to cluster 3 seem to be younger than the others, with a median age at 33,0 years at the time of the first symptoms and at 33.1 years at the time of the MRI diagnosis. When comparing this group of patients (n = 18) with all the others (n = 192), we found no significant difference as regards to age, with p = 0.11 for the age at the time of the first symptoms and p = 0.91 for age at the time of MRI diagnosis, respectively. However, it is likely that actual differences exist since the 95% confidence intervals for the difference between both mean ages were very unbalanced around the zero value: [-0.99; 9.48] for the age at the time of the first symptoms and [-0.74; 9.79] for the age at the time of the MRI diagnosis.

Discussion

In this paper, we introduced a probabilistic atlas that identifies the preferential locations of DLGGs in the brain. This atlas is based on a homogeneous series of non-treated DLGG, which is, to our knowledge, the biggest series of “pure” DLGG. We analyzed MRIs only if they had
<figure><img src="image_16.png" /><figcaption>Fig 9. Patient's age at the time of the first symptoms (a) and MRI diagnosis (b) for each cluster.</figcaption></figure>

doi:10.1371/journal.pone.0144200.g009

been performed before the patient had received any oncological treatment. Because the MRIs
in DICOM format and with FLAIR sequences were necessary for our analysis, we could not use
the first MRI performed (i.e. the MRI at diagnosis) for some patients (n = 122). However for
these patients the median delay between the first MRI and the MRI analyzed was short (125
days), and was over a year for 19 patients only. As DLGG are slowly growing tumors, this delay
is unlikely to cause dramatic changes in the tumor volume or to bias the determination of the
tumor center of gravity.

The atlas provides a precise description of the tumors locations and spatial extension in
each cluster (prior probability maps and clusters positions). The insula is undoubtedly a pre-
dominant location for the DLGG as more than a third of the tumors are situated in this area
(37% of the data). This result, as well as the quasi-total absence of tumors in or near the occipi-
tal and prefrontal lobes, is in accordance with what was observed in [17]. The proposed atlas
and complete database are available for download at http://db-gliomas-gradeii.net/.

Several applications for this atlas are imaginable. The proposed atlas constitutes a robust
method for correlations between DLGG topography and their characteristics. The patient's age
at diagnosis and first symptoms was considered in this paper, and a study with respect to the
tumors' molecular biology is currently ongoing.

We were not able to show any variability regarding both the age at the moment of the first
symptoms or the age at MRI diagnosis according to the tumor location. Such variability has
been recently suggested by Gerin et al, who identified two types of tumors among DLGGs, the
first one occurring during teenage years with a very slow growth, and the second one occurring
in early adulthood with a slow growth, [5]. In our series there is a trend towards a younger age
for patients with a tumor belonging to the cluster 3 involving the SMA compared to the other
locations, but it did not reach significance. However, it is likely that there is actually a difference
with regard to the age since the 95% confidence intervals of the difference between mean ages
of the “SMA group” and the “non-SMA” group, even though they contain the value zero, are
unbalanced towards zero ([-0.99; 9.48]) for the age at the time of the first symptoms and
[-0.74; 9.79] for at the time of MRI diagnosis, respectively). The small effective in the “SMA group” (n = 18) could partly explain our lack of statistic power. If such a difference does exist for DLGG involving the SMA, it would reinforce the fact that this area has specific biological, architectonic and functional properties, with particular interactions between glial cells and neurons that could, maybe, explain in parts earlier gliomatogenesis [17].

This study should be pursued with a larger database in order to establish whether those tendencies are confirmed. Furthermore, we have now designed a pipeline for correlating the position of the tumor with any parameter through the use of the probabilistic atlas. It is straightforward to extend the study to investigate the correlation of other parameters (molecular biology for instance) to the tumor location.

Furthermore, it would be interesting to couple this atlas with information obtained from functional imaging and study the relationship between tumors' position and brain functions. This could enable to identify compensable regions and predict the future functional reorganization of the brain depending on the location of the tumor. All of this would be extremely useful for surgery planning. Eventually, the same method could also be applied to other kinds of lesions, perhaps by changing the distance measure for one that is more adapted to the type of lesion.

We presented in this paper an application as a spatial position prior for tumor detection and segmentation. We managed to efficiently detect and characterize (i.e. assign to a cluster) the tumors from a sizable dataset and despite some unclear boundaries and heterogeneous appearances. As stated in the introduction, DLGG have different behaviors depending on their location. We could be able to infer the evolution of the tumor depending on the cluster it is assigned to.

Obtaining the best possible tumor segmentation is not the main focus of this paper. The proposed method was used to highlight the clinical relevance of the probabilistic atlas. We therefore used simple unary costs (boosting) while much more elaborated approaches from the literature could have been considered. It is therefore difficult to compare our results with the current state of the art from the BRATS challenge [50]. On top of that, we are using a very different and challenging database (monomodal data, variable quality, low resolution) due to our clinical setting.

Segmentation can be harder for small tumors, that do not offer an important contrast and are sometimes even difficult to detect visually. The impact of the spatial prior on the segmentation quality is directly dependent on the accuracy of the cluster assignment as well as how well the tumor fits in the cluster. As stated above, 12% of the data could not be represented by any cluster and therefore should not be assigned to any. The segmentation results would most likely be increased by the introduction of a label identifying outlier tumors, possibly by identifying if the tumor is equally distant to two or more clusters.

Despite the fact that we are already working on a large data-set, the quality and precision of the graph as well as the clusters probability distributions (i.e. the prior probability maps), could benefit from an increase of the size of the data set. Currently, we evaluate the segmentation's quality with respect to a manual segmentation by comparing the pixels detected as tumor to the pixels manually labeled as tumors. Manual segmentation is, unfortunately, a subjective approach and dependent on the operator [28]. Results could differ with another operator, especially in the case of DLGG where boundaries of the tumors can be very fuzzy.

The registration step also leaves room for improvement. First, the fact that all registrations are done towards an arbitrary reference pose introduces a bias in the atlas' structure. We should get better results using population registration methods [51] that register all the images together without the use of a reference pose. Second, the use of an affine registration scheme maintains the different brains anatomy but lacks precision in the matching of the anatomical structures. Increased precision would be obtained by aligning all volumes through a
deformable registration scheme that takes into account the missing correspondences induced by the presence of the tumor [33]. Due to the systematic processing of a potentially large scale database, deformable registration methods have to be made fully reliable before being used in our context.

## Author Contributions

Conceived and designed the experiments: SP NP HD AD CB. Performed the experiments: SP AD CB. Analyzed the data: SP AD CB. Contributed reagents/materials/analysis tools: SP HD MB LT YY LB SZ VR. Wrote the paper: SP AD HD NP. Provided the funding: SC.

## References

1. Duffau H (2009) Surgery of low-grade gliomas: towards a 'functional neurooncology'. Curr Opin Oncol 21: 543-549. doi: 10.1097/CCO.0b013e3283305996 PMID: 19606033

2. Rigau V, Zouaoui S, Mathieu-Daudé H, Darlix A, Maran A, et al. (2011) French brain tumor database: 5-year histological results on 25 756 cases. Brain Pathology 21: 633-644. doi: 10.1111/j.1750-3639.2011.00491.x PMID: 21554472

3. Dolecek TA, Propp JM, Stroup NE, Kruchko C (2012) Cbtrus statistical report: primary brain and central nervous system tumors diagnosed in the united states in 2005-2009. Neuro-Oncology 14: v1-v49. doi: 10.1093/neuonc/nos218 PMID: 23095881

4. Capelle L, Fontaine D, Mandonnet E, Taillandier L, Golmard JL, et al. (2013) Spontaneous and therapeutic prognostic factors in adult hemispheric world health organization grade ii gliomas: a series of 1097 cases. Clinical article. Journal of Neurosurgery 118: 1157-1168. doi: 10.3171/2013.1.JNS121 PMID: 23495881

5. Gerin C, Pallud J, Grammaticos B, Mandonnet E, Deroulers C, et al. (2012) Improving the time-machine: estimating date of birth of grade ii gliomas. Cell Proliferation 45: 76-90. doi: 10.1111/j.1365-2184.2011.00790.x PMID: 22168136

6. Soffietti R, Baumert BG, Bello L, von Deimling A, Duffau H, et al. (2010) Guidelines on management of low-grade gliomas: report of an efns-eano task force. Eur J Neurol 17: 1124-33. doi: 10.1111/j.1468-1331.2010.03151.x PMID: 20718851

7. Jakola AS, Myrmel KS, Kloster R, Torp SH, Lindal S, et al. (2012) Comparison of a strategy favoring early surgical resection vs a strategy favoring watchful waiting in low-grade gliomassurgical resection vs waiting in low-grade gliomas. JAMA 308: 1881-1888. doi: 10.1001/jama.2012.12807 PMID: 23099483

8. Laigle-Donadey F, Martin-Duverneuil N, Lejeune J, Crinière E, Capelle L, et al. (2004) Correlations between molecular profile and radiologic pattern in oligodendroglial tumors. Neurology 63: 2360-2362. doi: 10.1212/01.WNL.0000148642.26985.68 PMID: 15623700

9. Gozé C, Rigau V, Gibert L, Maudelonde T, Duffau H (2009) Lack of complete 1p19q deletion in a consecutive series of 12 WHO grade II gliomas involving the insula: a marker of worse prognosis? J Neurooncol 91: 1-5. doi: 10.1007/s11060-008-9680-8 PMID: 18726074

10. Metellus P, Coulibaly B, Colin C, de Paula AM, Vasiljevic A, et al. (2010) Absence of idh mutation identifies a novel radiologic and molecular subtype of who grade ii gliomas with dismal prognosis. Acta Neuropathologica 120: 719-729. doi: 10.1007/s00401-010-0777-8 PMID: 21080178

11. Stockhammer F, Misch M, Helms H, Lengler U, Prall F, et al. (2012) Idh1/2 mutations in who grade ii astrocytomas associated with localization and seizure as the initial symptom. Seizure: the journal of the British Epilepsy Association 21: 194. doi: 10.1016/j.seizure.2011.12.007

12. Duffau H (2005) Lessons from brain mapping in surgery for low-grade glioma: insights into associations between tumour and brain plasticity. Lancet Neurol 4: 476-486. doi: 10.1016/S1474-4422(05)70140-X PMID: 16033690

13. Desmurget M, Bonnetblanc F, Duffau H (2007) Contrasting acute and slow-growing lesions: a new door to brain plasticity. Brain 130: 898-914. doi: 10.1093/brain/awl300 PMID: 17121742

14. Ius T, Angelini E, de Schotten M, Mandonnet E, Duffau H (2011) Evidence for potentials and limitations of brain plasticity using an atlas of functional resectability of who grade ii gliomas: Towards a "minimal common brain". NeurolImage 56: 992-1000. doi: 10.1016/j.neuroimage.2011.03.022 PMID: 21414413

15. Chang EF, Gabriel RA, Potts MB, Berger MS, Lawton MT (2011) Supratentorial cavernous malformations in eloquent and deep locations: surgical approaches and outcomes. Clinical article. J Neurosurg 114: 814-827. doi: 10.3171/2010.5.JNS091159 PMID: 20597603
16. Yordanova YN, Moritz-Gasser S, Duffau H (2011) Awake surgery for WHO Grade II gliomas within "noneloquent" areas in the left dominant hemisphere: toward a "supratotal" resection. J Neurosurg 115: 232-239. doi: 10.3171/2011.3.JNS101333 PMID: 21548750

17. Duffau H, Capelle L (2004) Preferential brain locations of low-grade gliomas. Cancer 100: 2622-2626. doi: 10.1002/cncr.20297 PMID: 15197805

18. Ren X, Cui X, Lin S, Wang J, Jiang Z, et al. (2012) Co-deletion of chromosome 1p/19q and idh1/2 mutation in glioma subsets of brain tumors in chinese patients. PloS ONE 7: e32764. doi: 10.1371/journal.pone.0032764 PMID: 22427879

19. García C, Moreno JA (2004) Kernel based method for segmentation and modeling of magnetic resonance images. In: Advances in Artificial Intelligence-IBERAMIA 2004, Springer. pp. 636-645.

20. Zhang J, Ma KK, Er MH, Chong V, et al. (2004) Tumor segmentation from magnetic resonance imaging by learning via one-class support vector machine. In: International Workshop on Advanced Image Technology (IWAIT'04). pp. 207-211.

21. Verma R, Zacharaki EI, Ou Y, Cai H, Chawla S, et al. (2008) Multi-parametric tissue characterization of brain neoplasms and their recurrence using pattern classification of mr images. Academic Radiology 15: 966. doi: 10.1016/j.acra.2008.01.029 PMID: 18620117

22. Wels M, Carneiro G, Aplas A, Huber M, Hornegger J, et al. (2008) A discriminative model-constrained graph cuts approach to fully automated pediatric brain tumor segmentation in 3-d mri. In: Medical Image Computing and Computer-Assisted Intervention-MICCAI, Springer. pp. 67-75.

23. Bauer S, Nolte LP, Reyes M (2011) Fully automatic segmentation of brain tumor images using support vector machine classification in combination with hierarchical conditional random field regularization. In: Medical Image Computing and Computer-Assisted Intervention-MICCAI, Springer. pp. 354-361.

24. Lee CH, Schmidt M, Murtha A, Bistritz A, Sander J, et al. (2005) Segmenting brain tumors with conditional random fields and support vector machines. In: Computer Vision for Biomedical Image Applications, Springer. pp. 469-478.

25. Lee CH, Wang S, Murtha A, Brown MR, Greiner R (2008) Segmenting brain tumors using pseudo-conditional random fields. In: Medical Image Computing and Computer-Assisted Intervention-MICCAI, Springer. pp. 359-366.

26. Gering DT, Grimson WEL, Kikinis R (2002) Recognizing deviations from normalcy for brain tumor segmentation. In: Dohi T, Kikinis R, editors, Medical Image Computing and Computer-Assisted Intervention-MICCAI. Springer, volume 2488 of Lecture Notes in Computer Science, pp. 388-395.

27. Corso JJ, Sharon E, Dube S, El-Saden S, Sinha U, et al. (2008) Efficient multilevel brain tumor segmentation with integrated bayesian model classification. IEEE Transactions on Medical Imaging 27: 629-640. doi: 10.1109/TMI.2007.912817 PMID: 18450536

28. Kaus M, Warfield S, Nabavi A, Black P, Jolesz F, et al. (2001) Automated segmentation of MR images of brain tumors. Radiology 218: 586-591. doi: 10.1148/radiology.218.2.r01fe44586 PMID: 11161183

29. Prastawa M, Bullitt E, Moon N, Leemput K, Gerig G (2003) Automatic Brain Tumor Segmentation by Subject Specific Modification of Atlas Priors. Academic Radiology 10: 1341-1348. doi: 10.1016/S1076-6332(03)00506-3 PMID: 14697002

30. Moon N, Bullitt E, Leemput KV, Gerig G (2002) Model-based brain and tumor segmentation. In: ICPR (1). pp. 528-531.

31. Menze BH, Van Leemput K, Lashkari D, Weber MA, Ayache N, et al. (2010) A generative model for brain tumor segmentation in multi-modal images. In: Medical Image Computing and Computer-Assisted Intervention-MICCAI, Springer. pp. 151-159.

32. Prastawa M, Bullitt E, Ho S, Gerig G (2004) A brain tumor segmentation framework based on outlier detection. Medical Image Analysis 8: 275-283. doi: 10.1016/j.media.2004.06.007 PMID: 15450222

33. Parisot S, Wells W. III, Chemouny S, Duffau H, Paragios N (2014) Concurrent tumor segmentation and registration with uncertainty-based sparse non-uniform graphs. Medical Image Analysis 4: 647-659. doi: 10.1016/j.media.2014.02.006

34. Vovk U, Pernus F, Likar B (2007) A review of methods for correction of intensity inhomogeneity in mri. IEEE Trans Med Imaging 26: 405-421. doi: 10.1109/TMI.2006.891486 PMID: 17354645

35. Parisot S., Duffau H., Chemouny S., Paragios N. Graph Based Spatial Position Mapping of Low-grade Gliomas. In: MICCAI, pp.508-515, 2011.

36. Gooya A., Pohl K.M., Bilello M., Biros G., Davatzikos C., 2011. Joint segmentation and deformable registration of brain scans guided by a tumor growth model.: In: MICCAI, Part II, LNCS, vol. 6892, pp. 532-540

37. Komodakis N, Paragios N, Tziritas G (2008) Clustering via Ip-based stabilities. In: NIPS. pp. 865-872.
38. Dunn JC (1974) Well separated clusters and optimal fuzzy-partitions. Journal of Cybernetics 4: 95-104. doi: 10.1080/01969727408546059

39. Pal NR, Biswas J (1997) Cluster validation using graph theoretic concepts. Pattern Recognition 30: 847-857. doi: 10.1016/S0031-3203(96)00127-6

40. Davies DL, Bouldin DW (1979) A cluster separation measure. IEEE Transactions on Pattern Analysis and Machine Intelligence PAMI-1: 224-227. doi: 10.1109/TPAMI.1979.4766909

41. Rousseeuw PJ (1987) Silhouettes: A graphical aid to the interpretation and validation of cluster analysis. Journal of Computational and Applied Mathematics 20: 53-65. doi: 10.1016/0377-0427(87)90125-7

42. Friedman J, Hastie T, Tibshirani R (2000) Additive logistic regression: a statistical view of boosting. The Annals of Statistics 38: 337-407. doi: 10.1214/aos/1016218223

43. Manjunath BS, Ma WY (1996) Texture features for browsing and retrieval of image data. Pattern Analysis and Machine Intelligence, IEEE Transactions on 18: 837-842. doi: 10.1109/34.531803

44. Boykov Y, Funka-Lea G (2006) Graph cuts and efficient n-d image segmentation. International Journal of Computer Vision 70: 109-131. doi: 10.1007/s11263-006-7934-5

45. Komodakis N, Tziritas G, Paragios N (2008) Performance vs computational efficiency for optimizing single and dynamic mrfs: Setting the state of the art with primal-dual strategies. Computer Vision and Image Understanding 112: 14-29. doi: 10.1016/j.cviu.2008.06.007

46. Fonov VS, Evans AC, McKinstry RC, Almli CR, Collins DL (2009) Unbiased nonlinear average age-appropriate brain templates from birth to adulthood. NeuroImage 47: S102. doi: 10.1016/S1053-8119(09)70884-5

47. Fonov V, Evans AC, Botteron K, Almli CR, McKinstry RC, et al. (2011) Unbiased average age-appropriate atlases for pediatric studies. NeuroImage 54: 313-327. doi: 10.1016/j.neuroimage.2010.07.033 PMID: 20656036

48. Glocker B, Komodakis N, Tziritas G, Navab N, Paragios N (2008) Dense image registration through mrfs and efficient linear programming. Medical Image Analysis 12: 731-741. doi: 10.1016/j.media.2008.03.006 PMID: 18482858

49. Torresani L., Kolmogorov V., Rother C. (2008) Feature correspondence via graph matching: Models and global optimization. In: ECCV (2). pp. 596-609.

50. Menze B., Jakab A., Bauer S., Kalpathy-Cramer J, Keyvan Farahani, et al. The Multimodal Brain Tumor Image Segmentation Benchmark (BRATS). IEEE Transactions on Medical Imaging, Institute of Electrical and Electronics Engineers (IEEE), 2014, pp. 33.

51. Sotiras A, Komodakis N, Glocker B, Deux JF, Paragios N (2009) Graphical models and deformable diffeomorphic population registration using global and local metrics. In: MICCAI (1). pp. 672-679.