Glioma Dynamics and Computational Models

15 Pages·2007·0.53 MB·French

262 Current Medical Imaging Reviews, 2007, 3, 262-276
Glioma Dynamics and Computational Models: A Review of Segmentation,
Registration, and In Silico Growth Algorithms and their Clinical
Applications
Elsa D. Angelini, Olivier Clatz, Emmanuel Mandonnet, Ender Konukoglu, Laurent Capelle and
Hugues Duffau*
1GET, Ecole Nationale Supérieure des Télécommunications, CNRS UMR 5141, Paris, France, 2Asclepios Research
Project - INRIA Sophia Antipolis, France, 3Unité 678, Inserm/UPMC, 91 boulevard de l’Hôpital, 75634 Paris Cedex 13,
France, 4Department of Neurosurgery, Hôpital de la Salpêtrière, 47-83 Bd de l’hôpital, 75634 Cedex Paris Cedex 13,
France, 5Department of Neurosurgery, Hôpital Gui de Chauliac, CHU de Montpellier, 80 Avenue Augustin Fliche,
34295 Montpellier Cedex 5, France
Abstract: Tracking gliomas dynamics on MRI has became more and more important for therapeutic management.
Powerful computational tools have been recently developed in this context enabling in silico growth on a virtual brain that
can be matched with real 3D segmented evolution through registration between atlases and patient brain MRI data. In this
paper, we provide an extensive review of existing algorithms for the three computational tasks involved in patient-specific
tumor modeling: image segmentation, image registration, and in silico growth modelling (with special emphasis on the
proliferation-diffusion model). Accuracy and limits of the reviewed algorithms are systematically discussed. Finally
applications of these methods for both clinical practice and fundamental research are also discussed.
INTRODUCTION evolution with virtual in silico dynamics. Such a comparison
would require three different steps (as illustrated in Fig. 1):
The advent of MRI scanning protocols has allowed
segmentation of actual growth, registration on a virtual brain
accurate follow-up of tumor growth through volumetric
atlas, and identification of model parameters corresponding
measurements. Interpretation of the radiological evolution of
to optimal matching between actual and simulated evolution.
the tumor appears of utmost importance for therapeutic
management, especially for low-grade glioma. Indeed, We thus propose in this paper to review existing algo-
patients are most of the time asymptomatic (except in the rithms for segmentation, registration and in silico growth,
case of epilepsy) during the “low-grade” phase, and the with special attention to the accuracy and reliability of the
tumoral evolution is only monitored by MRI follow-up, both methods.
prior and after treatment. However, such information about
the tumor dynamics is usually not fully integrated with the SEGMENTATION
therapeutic strategy, and the assessment of the tumor
An accurate determination of the actual tumor evolution
evolution is still limited to qualitative descriptions (recur-
requires full 3D segmentation on digital images. Manual
rence, progression, regression, stability). Thus, it is expected
segmentation by an expert is still considered as the reference
that bio-mathematical models will help to quantify tumor
method, but is a time consuming task with high inter and
dynamics, to simulate treatment effects, and finally to
intra-observer variability.
optimize therapeutic strategies.
Many automated or semi-automated approaches were
Computational models of gliomas dynamics have been
developed over the past ten years showing great variability in
initiated more than ten years ago [1, 2]. The first studies
results and performance in terms of reproducibility.
aimed at modeling the effect of chemotherapy and surgical
Challenges in the segmentation of gliomas from MRI data
resection on the evolution of high-grade gliomas. If the
are related to i) the infiltration of cells into the tissue,
mathematical framework introduced at that time is still in
inducing unsharp borders with irregularities and disconti-
use, there have been considerable advances in its numerical
nuities (a tumor is not necessary a single connected object),
resolution. In particular, digital brain templates, provided by
ii) the great variability in their contrast uptake (depending on
MRI, enable to implement the biophysics equations onto
their vascularisation) and iii) their appearance on standard
accurate virtual anatomy. This in turn allows to refine the
MRI protocols.
model, by introducing for example different cell motility in
white and gray matter [3], and inside white matter, along and MRI protocols used for brain imaging typically include
orthogonally to axonal fasciculus [4]. However, published Proton density (PD), T1-weighted (e.g. SPGR), T1-weighted
studies have never seriously matched observed radiological enhanced (T1E) with contrast agent (usually Gadolinium),
and T2-weighted (e.g. FLAIR) data. T1 data provides detai-
led anatomical views of the brain along with high signals on
*Address correspondence to this author at the Department of
haemorrhages. T1E data shows strong signal on all vascula-
Neurosurgery, CHU Montpellier, Hôpital Gui de Chauliac, 80
rized structures (including tumors and haemorrhages),
Avenue Augustin Fliche, 34295 Montpellier Cedex 5, France; Tel:
33 4 67 33 02 90; Fax: [phone]; E-mail: h-duffau@chu- whereas usual FLAIR images (with a slice thickness around
montpellier.fr
1573-4056/07 $50.00+.00 ©2007 Bentham Science Publishers Ltd.
Glioma Dynamics and Computational Models Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 263

Fig. (1). Illustration of glioma growth evolution observed on temporal MRI studies.
5 mm) show less anatomical details, but high signal on In 2004 Letteboer et al. [7], proposed an interactive
tumors, infiltrations and edema. segmentation method for three types of tumors: full-
enhancing, ring-enhancing and non-enhancing. After manual
Methods tracing of an initial slice, a series of morphological filtering
operations (based on the watershed algorithm) was applied to
Given an image, the segmentation task can be seen as the
partition the MRI volume data into homogeneous areas. A
partition of the image into homogeneous objects, which
multiscale framework (i.e. analysis of the image data at
correspond to a region-based segmentation approach, or as
different spatial resolutions) was employed to correlate
the detection of object contours within the image,
segmented regions across different scales. The overall
corresponding to an edge-based segmentation approach. The
segmentation process was guided via an interactive user
majority of the MRI-based glioma segmentation methods
interface.
that have been proposed in the literature are region-based.
More recent methods, based on deformable models, also In 2005, Droske et al. [8], proposed to use a deformable
included edge-based information. In the case of MRI model, implemented with a level set formulation, to partition
segmentation, several factors introduce a large amount of the MRI data into regions with similar image properties,
uncertainty in the segmentation process, including partial based on prior intensity-based pixel likelihoods for tumoral
volume effects, integration of multi-protocol image data and tissues. The deformable model optimization was performed
observer variability. In this context, a large set of segmen- on a spatially-adaptive grid, only refined in inhomogeneous
tation methods was designed in a statistical framework, regions. Homogeneity measures included gray value inter-
providing a classification of the image data into different vals, defined from a user input, and image gradient values.
tissue types, while only few were designed with a deter- Some manual supervision of the deformable model was
ministic approach. required, so that incremental segmented areas were proposed
to the user who controlled the final segmentation results.
Deterministic Approaches
More specifically, heterogeneous tumors, involving necrosis
In 1996, Gibbs et al. [5] introduced a morphological edge for example, required successive segmentations by addition
detection technique combined with simple region growing to or removal of intermediate results.
segment enhancing tumors on T1 MRI data. Based on an
Statistical Approaches
initial sample of the enhanced tumor signal and the surroun-
ding tissues, provided manually, an initial segmentation was In 1995, Vaiddynathan et al. [9], compared two super-
performed combining pixel thresholding, fitting to an edge vised multispectral classification methods: k nearest
map of the image data and morphological opening and neighbour (kNN) and spectral fuzzy C-means (FCM). For
closing, inspired by the work proposed by Kennedy et al. [6]. these two classification approaches, nine tissue classes were
The tumor area was defined based on pixel values in the considered (background, CSF, WM, GM, fat, muscle, tumor,
range of 4 standard deviations around the mean value, edema, necrosis). The authors also tested an interactive seed-
constrained by the edge map. growing segmentation approach on T1E MRI data. The seed-
264 Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 Duffau et al.
growing algorithm only segmented tumor tissue based on a In 2001 Moonis et al. [14] proposed a segmentation
sample pixel population manually selected by the user. framework based on fuzzy connectedness (FC) which
optimally clustered voxels into classes of high connectivity
In 1998, Clark et al. [10] introduced a knowledge-based
(analogous to a similarity measure). The method was applied
(KB) automated segmentation method for glioblastomas on
to T1, T1E and T2 data, and initialized with an MRI data
multispectral data combining T1E, PD and T2 weighted data.
standardisation of the gray levels based on non-linear
A training phase was performed on 17 slices from seven
transformation of the histograms [15].
patients, extracting tumor size and enhancement level
characteristics. Slices were first characterized as normal or In 2005, Liu et al. [16], from the same group, used a
abnormal via a fuzzy C-means (FCM) classification and the similar approach based on a volume of interest on co-
analysis of the clustering result through an expert system. registered T1 and T2 data, to process only slices containing
Two examples of knowledge used in the predecessor system the tumor. A set of points inside the tumor were selected to
were: 1) in a normal slice, CSF belongs to the cluster center initialize the statistics used in the FC. The threshold level
with the highest value in the intracranial region; 2) in image applied to the FC maps to define the final segmentation
space, all normal tissues are roughly symmetrical along the result was determined empirically on five datasets and then
vertical axis. fixed once for all. Segmentation was performed separately
on the T2, T1E and subtracted (T1-T1E) data sets in 3D.
After a brain mask was computed, initial tumor segmen-
Manual corrections of the segmentation results were perfor-
tation, generated from vectorial histogram thresholding in
med by experts.
the T1, PD and T2 images, was post-processed with a KB
approach to eliminate non-tumor pixels. In 2001, Fletcher-Heath et al. [17], proposed a combi-
nation of unsupervised classification with FCM and
Tumor heuristics used in the KB system were the
knowledge-based (KB) image processing for segmentation
following: “1) Gadolinium-enhanced tumor pixels occupy
of non-enhancing tumors. The FCM was run on spectral data
the higher-end of the T1 spectrum; 2) Gadolinium-enhanced
(T1, T2, PD). As the authors pointed out, FCM tended to
tumor pixels occupy the higher-end of the PD spectrum,
define clusters with similar sizes, which required an initial
though not with the degree of separation found in T1 space;
classification in ten classes. A KB system was then designed
3) Gadolinium-enhanced tumor pixels are generally found in
to re-cluster the segmentation results into seven classes
the “middle” of the T2 spectrum, making segmentation based
based on a training phase. Difficulties principally arose in
on T2 values difficult; 4) Slices with greater enhancement
the separation of CSF and tumor signals.
have better separation between tumor and non-tumor pixels,
while less enhancement results in more overlap between In 2004, Mazzara et al. [18], compared the kNN appro-
tissue types”. It is important to note that their notion of ach from [9] and the KG-based approach from [10] for
tumor pixel included edema and necrosis. A final processing growth tumor volume (GTV) measurements on eleven
stage was performed, based on histogram analysis of the patients with high and low-grade gliomas. As used in
tumor pixels and heuristics on the “density” of intensity oncology radiation therapy, GTV corresponded to the area
features of non-tumor tissues. Indeed, based on the obser- enclosing several contiguous clusters of enhancing pixels
vation that tumors can show different levels of enhancement (i.e. including non-enhancing pixels within the area). The
and very complex shapes, the final KB approach was study showed severe limitations of the KG-system (which
focused on characterizing non-tumoral tissues. was not trained with the dataset to segment) in handling
particular cases such as non-enhancing tumor margins or the
In 2001, Kaus et al. [11] presented a complete validation
presence of non-enhancing cystic necrotic tissues at the
of an automated segmentation method on T1E data from
center of the tumor. On the other hand, the kNN segmen-
twenty patients with meningiomas and low-grade gliomas.
tation method, trained with sample data from MRI slices to
The segmentation method, called an adaptive template-
segment, lead to robust segmentation results on all patients.
moderated classification, and described in [12, 13] was based
In 2006, Beyer et al. [19], from the same group, presented a
on an iterative process. It alternated between a kNN
similar and more recent comparative study, extracting GTV
classification of voxels into five hierarchical tissue types
with the same two segmentation methods and evaluating the
(background, skin-fat-bone, brain, ventricles, tumor) and a
results in terms of predictive dose measurement for therapy
nonlinear registration of the data with an anatomical atlas
planning.
(manually segmented MRI data of a single subject) to align
the data with the template. The kNN classification used In 2004, Zou et al. [20], proposed a continuous proba-
features from data intensity values and anatomical priors on bilistic segmentation framework, based on mixture modeling
the tissue location from the atlas. This method performed for two classes: tumor and non-tumor tissues. After
extraction of the five tissues in a pre-determined hierarchical initialization of the segmentation with the semi-automated
order. Tissue mean values were learned on the patient’s data method from Kaus et al. [11], the segmentation process
via manual selection of three or four points for each tissue. involved estimation of the distribution parameters and
To handle the presence of the tumor in the registration probability values thresholding. Three metrics were proposed
process, voxels assigned to the tumor class were masked and evaluated to optimize the threshold selection: Receiver
with brain labels prior to registration with the atlas. This operating curve (ROC), which weights the sensitivity versus
method obviously relied on a strong homogeneity assump- the specificity of the segmentation result, a Dice similarity
tion of the tumor’s appearance on MRI data, which was coefficient, which is also a function of sensitivity and
reinforced by the use of anisotropic diffusion filtering. specificity and mutual information that directly compares the
segmentation result to a ground truth.
Glioma Dynamics and Computational Models Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 265
In 2004, Prastawa et al. [21] proposed a segmentation - Manual tracing bears some variability, which has been
framework based on outlier detection on T2 data. The evaluated in several studies as detailed below in this section.
abnormal tumor region was detected via registration on a
In the majority of papers, a standard manual segmen-
normal brain atlas. Statistical clustering of the abnormal
tation of the tumor (and eventually of the brain) was defined
voxels, followed by a deformable model, were then used to
from the segmentations of one or more independent human
isolate the tumor and the edema.
observers.
A summary of the reviewed papers and the specificities
For example, in Kaus et al. [11], pixels were assigned to
of the clinical aspects of the evaluation setup are provided in
the ground truth, if at least three of four observers agreed
Table 1.
with this assignment.
V
Gold Standard Definition t V TP FP
t
Evaluation of tumor segmentation accuracy typically FN V TP FN
requires the definition of a ground truth segmentation. In this
section we call (cid:1) the image data domain (i.e. the entire set
FP TP
of pixels to segment), V the tumor volume (i.e. set of pixels

identified as belonging to the tumor) provided by the
V
evaluated segmentation method, by V the tumor volume
t

identified by the ground truth method and by the
In Fletcher-Heath et al. [17], in Vaiddynathan et al. [9] in
cardinality (i.e. number of voxels included) of a volume. The Droske et al. [8] and in Prastawa et al. [21], ground truth was
complementary set of V (i.e. set of pixels identified as not defined from manual tracing from one expert.
belonging to the tumor) is written as V . To evaluate manual tracing variability, in Kaus et al.
[11], one observer manually segmented four times the same
Static as well as statistical ground truth volumes have
2D MRI slice, over one week.
been defined and used in the reviewed literature to compare
tumor segmentations. We review their definitions and In Liu et al. [16], two sets of manual tracings were used,
constructions in this section. based on 2 operators (a neuroradiologist and a trained expert
in medical imaging assisted by a semi-automated method),
Deterministic Ground Truth
and to evaluate segmentation precision, the two operators
The most widely used approach to define a ground-truth repeated the segmentation three times on ten T2, ten T1E and
tumor object, is via manual tracing of the contours by one or five T1E subtracted studies. To evaluate segmentation
several experts in neuroradiology. There are obvious accuracy, five T2 cases were manually traced three times.
limitations to such approach:
In Stadlbauer et al. [22], the authors compared manual
- Enhancement appearance of tumor on T1E data is very delineation of hyperintense signal on T2 data to high ratio of
variable and depends on the degree of vascularization as well Cho/NAA on MR spectroscopy imaging (MRSI). The
as permeability of the vessels to Gadolinium. ground truth for this study consisted of biopsy samples of
Table 1. Summary of Reviewed Papers and Clinical Setup
Data LGG Glioma HG Year
Gibbs T1E 10 1996
Letteboer 20 2004
Droske T1E ? 2005
Liu FLAIR, T1, T1E 10 2005
Vaidyanathan T1, PD, T2 4 1995
Fletcher-Heath T1, PD, T2 6 2001
Clark T1, PD, T2 (all with Gd) 7 1998
Kaus SPGR-Enh 14 2001
Moonis FLAIR 19 2001
Mazzara T1E, FLAIR (CT) 3 8 2004
Zou T1E (SPGR) 3 3 2004
Prastawa T2 1 1 2004
266 Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 Duffau et al.
brain tissues in the tumor region detected on MRSI but not Error Measurements for a Given Segmentation
on T2 data.
In Kaus et al. [11], validation was performed with
In their early work, Gibbs et al. [5] evaluated their measurements of segmentation accuracy (SA) defined as:
segmentation method on a simple phantom made of a water-
filled glass head containing a balloon of known volume (TP+TN)
(varying from 10 to 16 cm3 with 2cm3 increment) of CuSO4, SA= (2)
which was screened with the same T1E protocol as used for (cid:1)

patients. They also compared their method to a thresholding-
based segmentation method provided by a commercial soft- using the true positive (TP) volume overlap measured as:
ware (ISG Technologies Inc), requiring 2D manual initiali-
zations. Efficiency of the commercial method and the pro- TP= V (cid:1)V (3)
t
posed method was evaluated by repeating the initialization
three times, from a single user. and the true negative (TN) volume overlap measured as :
Statistical Ground Truth
TN = V (cid:1)V (4)
t
In Mazzara et al. [18], a statistical ground truth was
defined from a set of manual tracings as the probability that a It is important to note that this measure is very favorable
pixel is properly considered as part of the tumor. These to small objects such as circumscribed tumors.
probabilities were used as weights in the accuracy measure-
In Clark et al. [10] and in Fletcher-Heath et al. [17], the
ment, so that the true-positive accuracy measurement corres-
authors used two measurements:
ponded to the ratio of the total sum of weights contained in
the segmented area versus the total sum of weights generated percent match (PM) ratio defined as:
from nine manual tracings. A similar approach was used to
TP
measure false-positive values.
PM = (5)
In Zou et al. [20], the authors used the “simultaneous V
t
truth and performance level estimation” (STAPLE) method
introduced by the same group in [23], to compute a proba- This measure corresponds to a TP volume fraction
bilistic estimate of the true segmentation, given multiple (TPVF).
manual tracings, and provided a measure of the performance
correspondence ratio (CR) defined as:
level represented by each segmentation.
FP
In Letteboer et al. [7], the authors considered that manual TP(cid:1)
2
tracing did not provide ground truth per se and compared CR = (6)
V
their segmentation results to such tracings via a Bland and
t
Altman statistical analysis. Considering that differences in
An additional comparison was performed, evaluating the
volume measurements depend on the tumor size, which is
segmentation method for accuracy and precision of longi-
correlated with the fact that the majority of the segmentation
tudinal tumor evolution. In this study they evaluated four
errors occur on the surface of the tumor, they proposed a
normalization of the volume difference values D evaluated patients with at least two scans and failed in measuring a
growing tumor in one case, where a significant amount of
in the test:
fluid at the initial time lead to an overestimation of the initial
(V (cid:1) V ) tumor size.
t
D= , (1) Analogous measures were used by Letteboer et al. [7] as:
((V (cid:1) V )+(V (cid:1) V )) 2
D t t E TP
PM = , (7)
where V is the segmented volume dilated by one pixel and ((V +V ) 2)
D t
V is the segmented volume eroded by one pixel. This
E and
normalized difference value directly correlates the D value
with the number of different voxels. V (cid:1)V
t
CR= (8)
(( ) )
Clinical Validation Methods V +V 2
t
Given a ground truth representation of the tumor, several
In Liu et al. [16] ,ten segmentations of patients with
error measurements can be used to evaluate the accuracy of
glioblastomas were evaluated with the following. Measur-
the segmentation method, including the notions of speci-
ments:
ficity, sensitivity, repeatability and efficiency, which are
discussed in this section. 1. Precision: reproducibility of the segmentation varying all
parameters for the scanning protocol (i.e. using repetitive

scans of individual patients) and for the segmentation
Glioma Dynamics and Computational Models Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 267
algorithm. This measure is analogous to the concept of N
reproducibility or repeatability. (cid:1)Vi (cid:2) (cid:3)Vi
t t
2. Accuracy: comparison of segmented regions with i=1 (12)
respect to manual tracing. Two volume fraction (VF) N
measures were used: (cid:1)Vi / N
t
`FPVF: i=1
Based on this disagreement ratio, inter-observer
V (cid:1)V variability was evaluated among three experts, comparing
t
(9) one expert to the two other ones.
V
t In Letteboer et al. [7], three operators segmented the
tumors twice. Observer variability was assessed through a
and FNVF:
Bland and Altman statistical analysis, measuring deviations
V (cid:1)V from an average volume value. From this analysis,
t (10) variability can be measured as the Bias ± 1.96 SD, where SD
V is the standard deviation of the volume measurements.
t
In Prastawa et al. [21], area overlap and surface distances
Efficiency: computation and manual intervention time. were used to evaluate inter-observer variability.
Similar measurements were also evaluated in Letteboer
et al. [7] on twenty patients. Results
In Gibbs et al. [5], efficiency of the commercial method In Gibbs et al. [5], correlation of volume measurements
and the proposed method was evaluated by repeating the based on phantom data was perfect, with less then 5% errors.
initialization three times, from a single user. Patient data segmentation required 2-3 minutes of user
intervention and around 10 minutes of computational time.
In Zou et al. [20], the three metrics evaluated to optimize
Segmentation with the commercial software required about
the segmentation method were also used to evaluate the
30 minutes of user time. Both methods provided between 9%
segmentation results: area under the receiver operating curve
and 13% precision in volume measurements (for tumor
(AUC) which weights the sensitivity versus the specificity of
volumes in the range of 2 cm3-80 cm3), with a mean diffe-
the segmentation result, a Dice similarity coefficient (DSC), rence in observations of 0.1± 4.5 cm3. Mean difference bet-
which is also a function of sensitivity and specificity and
ween two observers is 0.8± 1.8 cm3. In conclusion, despite
mutual information (MI) that directly compares the segmen-
its simplicity and limitation to well enhanced tumor, not
tation result to a ground truth.
close to the skull, this early method provided fast and
In Prastawa et al. [21], the authors used three error accurate global tumor volume measures.
metrics, including surface comparison of the tumor’s outline.
In Letteboer et al. [7], average intra and inter-observer
Theses metrics were the PM overlap measure, the Hausdorff
volume differences were 3.2% and 9.7%, while average
distance and average surface distance.
manual PM values were 93.5% and 90%. Variability with
Inter and Intra-Observer Variability Bland and Altman analysis was 0.04 ± 1.79 for manual
tracing and -0.01 ± 0.76 for the segmentation method. They
Inter and intra-observer variability is typically measured
concluded that the watershed method was more reproducible
with the coefficient of variation (CV) of volume measure-
than manual tracing. A negative bias showed under-
ment defined as:
segmentation from the watershed and they obtained best PM
similarity measures for enhancing tumors and worst for non-
(cid:2)
CV%= Vt (cid:1)100 (11) enhancing ones.
μ In Droske et al. [8], computation time, including manual
V
t initialization required 3 minutes. They concluded that accu-
racy was high for homogeneous enhanced tumors and lower
where μ is the mean value and (cid:1) is the standard
V V for heterogeneous tumors, without quantitative numbers.
t t
deviation of volume measurement made on V . In Vaidyanathan et al. [9], precision of manual tracing
t
and three segmentation methods was evaluated in terms of
In Moonis et al. [14] the authors used the CV measure to
user input and showed very large variability. Manual tracing
evaluate variability of the user tracings and inputs to the
average variability (inter-intra) was 6%-17%. Segmentation
segmentation method. In Vaiddynathan et al. [9], the authors
methods average variability was around 5%-8% for trained
also used the CV measures to evaluate reproduci-bility of
classification methods and 17%-6% for a seed growing
tumor segmentation with respect to parameter setting of the
method. These results demonstrated the better precision of
different segmentation methods.
trained-based segmentation approaches but also the
In Mazzara et al. [18], intra-observer variability of weakness of these methods in terms of sensitivity to user
manual tracing was assessed with three tracings from one inputs. They concluded that reproducibility was a major fac-
expert, and measured via the ratio of average disagreement: tor affecting tumor volume measurements (with differences
268 Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 Duffau et al.
above 20 ml in their study), with a critical need for less- In Mazzara et al. [18], average intra-observer variability
supervised methods. was very large (20%±16%). Average inter-observer
variability was also high (28%±12%). They observed that
In Fletcher-Heath et al. [17], they obtained a PM bet-
they obtained better average reproducibility in preoperative
ween 53% and 90% and a CR between 0.37 and 0.87 over
cases (15% - 24%) than in postoperative cases (27%-32%).
six cases. For the two cases used for training the KB system,
In terms of accuracy measurements based on a statistical
average PM was 72% and average CR was 0.67. For the four
ground truth, manual tracing provided 85%±7%, KNN
tested cases, average PM was 75% and average CR was 0.56.
provided 56%±6% and KG method provided 52%±7%,
In Clark et al. [10], overall results showed an overesti- while FP values were 8%±11% for the KNN method,
mation of the tumor volume with automated segmentation. 8%±8% for the KG method and 17%±11%, showing a
Regarding longitudinal studies, the methods also showed one tendency of the segmentation method to underestimate tumor
case where tumor shrinkage was falsely measured, volume. The authors concluded that the segmentation
potentially due to high haemorrhage signal in the baseline accuracy was within the manual tracing variability range. It
MRI scan. They also reported 5% inter-observer variability is important to note that the accuracy of the automated
in tumor volume manual measurement. For the proposed KB methods could not have exceeded the one of the manual
segmentation method, repeatability was fully guaranteed. tracing with the proposed accuracy measurements with a
PM was typically above 0.9 (range=0.69:1, average= 0.93), statistical ground truth. Average computational time for the
and CR was more variable (range=0.43:0.85, average=0.66). KNN approach was 30 minutes, including training time. For
They also showed that less than 15% of the FP were really the KG approach, preparing MRI scans for segmentation
FP, not spatially connected to any ground-truth tumor pixel. required 90 minutes of operator time, and segmentation time
Regarding a kNN- based method, they reported poorer PM required 30 minutes. In Beyer et al. [19], the same group
measures (range=0.22:0.99, average= 0.71) and CR showed that expert physician reference volumes were
measures (range=-2.21:0.64, average =0.19). irradiated within the same level of conformity when using
radiation plans generating from automated segmented
In Kaus et al. [11], evaluation of the proposed segmen-
tation method reported extremely high SA (over 99%), with contours.
very low segmentation variability. Intra and inter -observer In Zou et al. [20], high AUC and DSC measures were
variability (CV) for manual and automated segmentation was obtained on 2D slices, from nine cases, but they observed
below 5%-15% and 4%-7%. This again suggested the gain in that the recommended optimal threshold for the probability
precision when using an automated (classification-based) values seemed to be case- and task metric-dependent, largely
segmentation method. Computation time for a three- depending on the clinical goal of the segmentation task. For
dimensional data set was 75 minutes. example, AUC was suited for overall accuracy while MI was
better suited for longitudinal studies of the tumor’s evo-
In Moonis et al. [14], on a set of ten patients, average CV
lution.
for two different users was 0.27-0.21 on FLAIR data and
0.37-0.27 on T1E data. Inter-observer CV values were 0.38 In Prastawa et al. [21], inter-observer variability was very
and 0.29 for FLAIR and T1E data, with non-significant variable, between 59% and 89% of area overlap. Average
differences. The study reported that manual editing of the surface distances remained below 2mm while the Hausdorff
segmentation results lead to smaller tumor volumes, with a distance reached 13mm for one case. Surface overlaps with
median change of volume of -17%. No significant difference the automated segmentation method were comparable,
was found in CV values with or without manual editing. between 68% and 80%, the Hausdorff distance reached
18mm and the average distance 4mm. Computa-tional time
In Liu et al. [16] ,ten patients with glioblastomas were
was around 90 minutes per volume.
evaluated. Manual correction by experts was performed on
all segmentation outputs evaluated. Results showed:
Summary
- Precision: for repeated segmentations by two experts, the
average CV measure for each expert and inter-observer Based on the literature review of tumor segmentation
CV were: (0.33-0.26, 0.43) on FLAIR data. (0.33-0.3, from MRI data it is extremely difficult to conclude on a best
available method. Obviously, automation and minimal user
0.36) on T1E data and (0.37-0.33, 0.41) on subtracted
supervision as in [9], are desirable. Obviously, atlas-based
T1E. Regarding these three experiments, intra- and inter-
methods are limited to tumors without any mass effect. To
observer percentage overlap was above 98%. Mean
volume variation was 1.2% over two repeated FLAIR cover the range of sophistication in segmentation methods,
scans on five patients. simple seed growing showed very poor reproducibility,
while integration of multispectral MRI data, from several
- Accuracy: On five FLAIR cases, average FPVF was protocols, seems critical to mimic visual interpretation of the
4.45% and average FNVF 4.28%. Volume estimation tumor borders from neurologists.
accuracy was above 95%.
Recent methods have included a great amount of
- Efficiency: computational time was 16 minutes (7
interactive manual supervision of the segmentation process,
minutes for registration, and 8 minutes for operator
reflecting multiple observations of high variability from
supervision). manual or automated tracing of tumors. Moreover, MRI data
Manual tracing average inter and intra variability was provides images with varying tumor appearance due to the
measured below 2% on the three protocols. heterogeneity of the tumor physiology as well as important
variations in MRI scanners in terms of image quality.
Glioma Dynamics and Computational Models Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 269
Related to this factor, a recent paper from Angelini et al. [24] tration algorithm. Registration is a primordial step to com-
proposed a longitudinal method for quantification of low- pare the actual segmented evolution with virtual simulated
grade glioma evolution based on non-linear image normali- growth.
zation and direct difference comparison, avoiding the need
Image registration can be defined as the process of
for specific tumor segmentation. As an alternative, notions of
aligning a target image to a source image. It consists in
fuzzy segmentation decisions are used in recent works from
determining the geometrical transformation that maps
Khotanlou et al. and Dou et al. [25, 26].
corresponding structures in the target and the source image.
Explicitly the registration algorithm aims at deforming the
Future Improvements
source images so that similar structures appear at the same
- Performing multispectral image segmentation: we can location in the images.
believe that classification or deformable models incor-
Historically, the primary interest of clinicians was the
porating multispectral data will show superior perfor-
fusion of multimodal images of the same subject. This
mance within the next few years.
registration -usually based on rigid transformations- allows
- Role of MRI spectroscopy imaging (MRSI): in a recent the comparison of different images from a single subject and
study, Stadlbauer et al. [22] showed on ten patients with the visualization of multiple modalities at the same spatial
low-grade (I&II) gliomas that abnormal tumor areas on location (see Fig. 2 for an example with MR T1 and MR T2).
1H-MRSI exceeded by an average of 24% that area It could also be used to finely analyze the progression over
delineated via supervised seed-growing segmentation of time of an evolving process (such as lesion growing) [28].
T2 hyperintense signal. The study also showed that
Once mono-subject rigid registration algorithms were
computation of metabolic ratios on MRSI enabled
considered mature, computer scientists started developing
automated segmentation of tumoral areas. On the other
non rigid registration methods for multi-subjects registration.
hand, spatial resolution of spectroscopic data is still
Those algorithms aim at finding the deformation that maps
limited, with a voxel volume close to 1 mL (compared to
different subjects in the same space. Such methods are
0.005 mL in MRI).
usually used in combination with atlases. An atlas is an
- Different evaluation setups: It is essential to evaluate anatomical image of a subject coupled with a second image.
brain tumor segmentation methods in a clinical setting. In This associated image represents segmented structures of
a study from Pallud et al. [27] biopsy samples isolated interest, local diffusion properties or label probabilities. The
tumoral cells beyond imaging abnormalities observed on registration algorithm is used to compute the deformation
T2 and FLAIR MRI data in 17 patients with low-grade field from the anatomical image of the atlas to the image of
gliomas. In the context of radiation therapy, Beyer et al. the subject. This deformation field is then applied to the
[19] have recently shown the superiority of knowledge- property image, so that this information can be mapped to
based segmentation methods in automatically deter- the patient’s anatomy.
mining brain gross tumor volume (GTV), with higher
reproducibility than manual tracing. Methods
Registration Algorithm for Healthy Subject
REGISTRATION
The majority of registration algorithms share the same
Introduction
three components. First the type of transformation, which
Since numerical simulations are usually computed in a controls the geometric flexibility of the algorithm, needs to
different space from the one encompassing the patient’s be defined. The transformation T() defines the mathematical
anatomy, modeling results have to be warped back to match formulation that relates a point X in the first (floating) image
the specific patient geometry. This deformation of the refe- and a point X’ in the second (fixed) image:X'=T(X).
rence space to the patient space is achieved through a Regi-

Fig. (2). Rigid registration of MR T1 (left) and T2 (middle) images of the same patient. The registered T2 image (right) now have
homologous structures displayed at the same location in the image.
270 Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 Duffau et al.
Usually, the transformation is assumed to be rigid between similarity criteria [31, 44]. The displacement of the tumor
images of the same subject. This transformation is region is then guided on its boundary by the surrounding
parameterized by 6 degrees of freedom:X'= RX+D, 3 healthy structures, and the regularization criterion ensures
for translation (vector D) and 3 for rotation (matrix R). the continuity of the displacement inside the tumoral region.
The registration algorithm could be considered as passive in
Different subjects are registered using a non-rigid
the tumor region, in the sense that it follows the motion of its
transformation. Various transformation models have been
environment.
proposed in the literature. These models were chosen for
their mathematical simplicity, and none of them relies on a An alternative approach consists in defining a model of
priori brain variability. However, they allow for multi- the tumor-induced displacement. Initially, a statistical
resolution registration, meaning that the user could define the method was proposed, trained on a dataset of possible
level of complexity of the transformation through a specified displacements [45]. These methods make the assumption that
number of degrees of freedom. Most common transformation similar tumors have similar deformation patterns; so that the
models are splines [29], cosine basis [30], tetrahedral mesh deformation induced by a new tumor can be deduced with a
[31, 32], multi-affine [33-35] and free (one vector per voxel) linear combination of other displacement induced by tumors
at the same location. However, these statistical methods need
[36-38].
a large number of patient images and the corresponding true
The second component of a registration algorithm is a displacements. To overcome this problem, it was proposed to
similarity metric between the two images I and I’. This train the statistical model on tumor growth simulations [46,
metric will depend on the assumptions made on image 47]. Numerical simulation could then be used to train the
intensities distribution for similar structures in the two model for any tumor and at any location in the brain. In this
images [39, 40]. For instance, if similar structures are case, the registration algorithm could be considered as
considered to have homologous intensities (as for mono- active, in the sense that it tries to find an appropriate model
modal registration), the appropriate metric will be the sum of to fit to observed growth in the image. Nevertheless, using
squared intensities: (cid:1)(I(X)(cid:2)I(X'))2 . If the relationship such registration methods for tumor growth simulation in an
atlas image is highly controversial.
between the two intensities is supposed linear (as for mono Clinical Validation Methods
modal registration, but from different MR scanners), the
Validation of non-rigid registration algorithms is a
correlation coefficient should be used:
research topic on its own. Because it is a relatively recent
(cid:1)(I(X)(cid:2)I)(I'(X')(cid:2)I')
research domain, emphasize has mostly been put on the
, where I defines the average
development of new algorithms. In addition, validation
(cid:1)I(X)I'(X)(cid:2)I'I
methods highly depend on the application of the registration
intensity of image I. For more complex and probabilistic algorithm.
relationship, mutual information is the natural choice [41, When the displacement observed between the two images to
42]. Because of its adaptability to multi-modal images, be registered is due to a mechanical deformation (for
mutual information is now the default similarity measure example: a tumor induced mass effect), comparison of
used in most registration algorithms. manually identified landmarks seems to be the standard
Usually, the transformation model does not explicitly method of evaluation. Points on similar structures are
impose the smoothness of the deformation. This smoothness manually identified in the moving image (Pm) and the
is then ensured through a regularization component in the reference image (Pr). The computed displacement of the
global registration formulation. Popular smoothing energies moving point T(P ) is then compared to the manually
m
rely on second order derivative of the displacement field like evaluated displacement D=(P -P) and the error is defined as
m r
thin plates splines [29] or Laplacian [37] and continuum the sum of the squared differences: (cid:1)(T(P )(cid:2)D)2 . This
mechanics based energy [31]. It has recently been proposed m
to take into account inter-subject variability to constrain the method is usually considered as the standard method for
regularization energy [35, 43]. mono-subject registration validation. It is however subject to
inter-expert variability and human error in identifying
These components are integrated into the registration
corresponding points.
algorithm through an optimization process based on an
energy formulation. The algorithm tries to estimate the In the context of multi-subject registration, the objective
parameters of the transformation that minimizes an energy is different: the deformation now represents the variability
composed of the similarity criteria and the regularization between the two individuals and is more difficult to evaluate.
energy. In addition, distinguishing between variability-induced
deformation and tumor-induced deformation is challenging
Registration Algorithm in the Presence of Tumor
for pathological cases. However, the registration algorithm is
The presence of a tumor in images represents a challenge usually used in combination with an atlas to import the
for registration algorithms: the assumption of an intensity associated property image in the patient space. It then makes
relationship between homologous structures does not hold, more sense to validate the quality of this final matching. For
as the presence of a tumor in the image changes the MR example, the registration of a segmented atlas will be
signal in invaded areas. The initial approach to tackle this evaluated on the quality of the segmentation on the new
problem consisted in discarding the tumor region from the image [48].
Glioma Dynamics and Computational Models Current Medical Imaging Reviews, 2007, Vol. 3, No. 4 271
Recently, new methods for validation based on atlas these algorithms. To the best of our knowledge, LesionMask
building have been proposed [49, 50]. In [50], multiple [52] is the only non-rigid registration algorithm available
registration algorithms are used to build an average diffusion that could be used in the presence of tumors. More recent
tensor image. The quality of the registration is then evaluated tools based on ITK have been proposed for the registration
on the variability (variance) around the mean image: the of healthy subjects [54, 55]. It is probable that those
better the registration, the closer every deformed diffusion algorithms will soon be adapted to account for the tumor
tensor images should be to the mean. In [49], an atlas of deformation in the images.
resection probability is build. The performance of the
registration algorithm is then evaluated using the predictive IN SILICO GLIOMA GROWTH
power of the atlas on new cases.
The mathematical description of tumor growth can be
formulated at different spatial scales: either one tries to
Results
simulate the growth at the cellular level (cellular automata),
In [31], the registration error is evaluated using manually or one defines on a macroscopic scale how the tumor density
identified landmarks. Results are obtained on mono-subject will evolve (with partial differential equations). In both
registration of patients in the course of tumor resection. cases, the models attempt to predict the mathematical law of
Images are acquired with an open MR system. The average the tumor growth.
registration error measured on 6 cases and 54 landmarks is
The cellular automata approach has been proposed for
bellow 1mm, and the maximum error bellow 3mm. It is
different tumors, including high grade gliomas.[56-58]. The
mentioned that an accuracy decrease is observed in the
rules of division and invasion are the key elements of this
regions very close to the tumor. The algorithm is not
approach. Since these models describe the tumor growth at
distributed.
the microscopic scale, their prediction could be also
In [46], the impact of the active tumor model is evaluated validated in clinical practice by microscopic histological
on atlas-to-subject image registration of real and simulated analysis of tumor samples (spatial correlation in graph cells
tumors. The benefit of the active model over the usual [59]) or in vitro, with dynamic study of glial cells migrations
passive registration algorithm is evaluated using landmark [60]. On the contrary, the partial differential equations
errors: the average accuracy improvement is 58% for real approach does not tell anything about the spatial ordering of
tumors, and 39% for simulated tumors. The average error is the tumor cells: it simulates a coarse-grained cell density.
bellow 4mm in both cases. Part of the algorithm is
There is of course a link between these two spatial scales
distributed [51], but it does not include the active tumor
of description: partial differential equations can be solved by
model.
stochastic methods, mimicking the cellular phenomena
A first step towards the validation of registration (Gaussian random walk for example). Fractals models [61-
algorithms for healthy subject has been proposed in [50]. 63] also propose to make a link between these two scales. By
The quantitative study of Sanchez et al. shows that best determining the fractal dimension of the tumor boundaries, a
registration algorithms available for DTI mapping (where scaling analysis gives critical exponents from which
splines and demons are compared) do not have statistical dynamics law can be inferred. This approach seems to be in
difference in their respective measured errors. good agreement with both in vivo and in vitro data of bulky
tumors, but its relevance for infiltrative tumors like glioma is
In [49], the software LesionMask [52] (available for
not established.
download) is used to build a tumor resection probability
atlas. Retrospective evaluation of the predictive power of the We will now focus on the macroscopic model of tumor
atlas for the preoperative classification of subtotal versus growth, which seems to be the most appropriate for compa-
partial tumor resection shows a correct prediction in 82% of rison with clinical MRI data, and is the most widely used. Of
cases. Such method could be used in the very near future to interest is a phenomenological approach very recently
assess the relative performance of registration algorithms in proposed [64]: in this case, machine learning algorithms are
the presence of tumors. used to determine, based on observed growth patterns, a 3D
probabilistic classification of diffusion patterns. This
Summary promising approach, yet preliminary, will help in the future
to refine the deterministic model of proliferation-diffusion
To summarize, different registration methods have been
that we will now present.
proposed in the literature to tackle the problem raised by the
presence of tumors in the images. There is however today no
Methods
consensus on a satisfactory registration method able to map a
patient image with a tumor on an (healthy) atlas. Indeed, Proliferation-Diffusion Model for Tumor Cell Density
validation of registration algorithm seems to be very
This modeling framework is based on an equation
application dependant. Most promising validation methods
governing (coarse-grained) tumor cell density (denoted c,
seem to be based on atlas building [49, 50].
expressed in cells/mm3) increase by time unit [4]. The
The mathematical formulation of the problem seems to generic form of the equation is the following:
be now better understood. It is then expected that recent (cid:1)c
=(cid:2)c+(cid:3).(D(cid:3)c), (13)
efforts in the open source community [53], implementing (cid:1)t
most popular algorithms and making them available, will stating that, at a given spatial location, new tumor cells
improve the research efforts in comparing and validating appear either by division (proliferation), or by moving from