Standardized assessment of automatic segmentation of white matter hyperintensities; results of the WMH segmentation challenge

1 Standardized Assessment of Automatic Segmentation of White Matter Hyperintensities; Results of the WMH Segmentation Challenge Hugo J. Kuijf, J. Matthijs Biesbroek, Jeroen de Bresser, Rutger Heinen, Simon Andermatt, Mariana Bento, Matt Berseth, Mikhail Belyaev, M. Jorge Cardoso, Adri` a Casamitjana, D. Louis Collins, Mahsa Dadar, Achilleas Georgiou, Mohsen Ghafoorian, Dakai Jin, April Khademi, Jesse Knight, Hongwei Li, Xavier Llad´ o, Miguel Luna, Qaiser Mahmood, Richard McKinley, Alireza Mehrtash, S´ ebastien Ourselin, Bo-yong Park, Hyunjin Park, Sang Hyun Park, Simon Pezold, Elodie Puybareau, Leticia Rittner, Carole H. Sudre, Sergi Valverde, Ver´ onica Vilaplana, Roland Wiest, Yongchao Xu, Ziyue Xu, Guodong Zeng, Student Member, IEEE, Jianguo Zhang, Guoyan Zheng, Member, IEEE, Christopher Chen, Wiesje van der Flier, Frederik Barkhof, Max A. Viergever, Fellow, IEEE, and Geert Jan Biessels H.J. Kuijf and M.A. Viergever are with the Image Sciences Institute, UMC Utrecht, Utrecht University, the Netherlands. J.M. Biesbroek, R. Heinen, and G.J. Biessels are with the Brain Center Rudolf Magnus, UMC Utrecht, Utrecht University, the Netherlands. J. de Bresser is with the Department of Radiology, UMC Utrecht, the Netherlands and with the Department of Radiology, LUMC, Leiden, the Netherlands. S. Andermatt and S. Pezold are with the Department of Biomedical Engineering, University of Basel, Allschwil, Switzerland. M. Bento is with the Radiology and Clinical Neuroscience, Hotchkiss Brain Institute, University of Calgary, AB, Canada. M. Berseth is with NLP Logix. M. Belyaev is with the Skolkovo Institute of Science and Technology. M.J. Cardoso is with the School of Biomedical Engineering and Imaging Sciences, King’s College London and the Centre for Medical Image Computing, University College London. A. Casamitjana and V. Vilaplana are with the Signal Theory and Communications Department, Universitat Polit` ecnica de Catalunya, BarcelonaTech, Barcelona, Spain. D.L. Collins and M. Dadar are with the McGill University, Canada. A. Georgiou is with the Computational Statistics and Machine Learning MSc, University College London. M. Ghafoorian is with TomTom, Amsterdam, the Netherlands. D. Jin and Z. Xu are with the Department of Radiology and Imaging Science, National Institutes of Health, USA. A. Khademi is with the Ryerson University, Canada. J. Knight is with the University of Guelph, Canada. H. Li is with the Sun Yat-sen University, University of Dundee, and Technical University of Munich. X. Llad´ o and S. Valverde are with the Research institute of Computer Vision and Robotics, University of Girona, Spain. M. Luna and S.H. Park are with the Department of Robotics Engineering, Daegu Gyeongbuk Institute of Science and Technology, Daegu, South Korea. Q. Mahmood is with the Pakistan Institute of Nuclear Science and Technology. R. McKinley and R. Wiest are with the Support Center for Advanced Neuroimaging, Institute for Diagnostic and Interventional Neuroradiology, Inselspital, University of Bern, Switzerland. A. Mehrtash is with the Electrical and Computer Engineering Department, University of British Columbia, Vancouver, BC and the Department of Radiology, Brigham and Women’s Hospital, Harvard Medical School, Boston, MA. S. Ourselin is with the is with the School of Biomedical Engineering and Imaging Sciences, King’s College London. B. Park is with the Department of Electronic, Electrical and Computer Engineering, Sungkyunkwan University and the Center for Neuroscience Imaging Research, Institute for Basic Science (IBS), Suwon, Korea. H. Park is with the Center for Neuroscience Imaging Research, Institute for Basic Science (IBS) and the School of Electronic and Electrical Engineering, Sungkyunkwan University, Suwon, Korea. E. Puybareau is with the EPITA Research and Development Laboratory (LRDE), France. L. Rittner is with the School of Electrical and Computer Engineering, University of Campinas, SP, Brazil. C.H. Sudre is with the School of Biomedical Engineering and Imaging Sciences, King’s College London; and the Centre for Medical Image Computing and the Dementia Research Centre, Institute of Neurology, University College London.. Yongchao Xu is with the EPITA Research and Development Laboratory (LRDE) and the LTCI, T´ el´ ecom ParisTech, Universit´ e Paris-Saclay, France, and the Huazhong University of Science and Technology, China. G. Zeng and G. Zheng are with the Institute for Surgical Technology & Biomechanics, University of Bern, Bern, Switzerland. J. Zhang is with the University of Dundee. C. Chen is with the Memory Aging and Cognition Center, NUHS, Singapore. W. van der Flier is with the Alzheimer Center, VU Amsterdam, the Netherlands. F. Barkhof is with the Department of Radiology & Nuclear Medicine, VU University Medical Center, Amsterdam, the Netherlands and the UCL institutes of Neurology and Healthcare Engineering, London, UK. Abstract—Quantification of cerebral white matter hyperintensities (WMH) of presumed vascular origin is of key importance in many neurological research studies. Currently, measurements are often still obtained from manual segmentations on brain MR images, which is a laborious procedure. Automatic WMH segmentation methods exist, but a standardized comparison of the performance of such methods is lacking. We organized a scientific challenge, in which developers could evaluate their method on a standardized multi-center/-scanner image dataset, giving an objective comparison: the WMH Segmentation Challenge (https://wmh.isi.uu.nl/). Sixty T1+FLAIR images from three MR scanners were released with manual WMH segmentations for training. A test set of 110 images from five MR scanners was used for evaluation. Segmentation methods had to be containerized and submitted to the challenge organizers. Five evaluation metrics were used to rank the methods: (1) Dice similarity coefficient, (2) modified Hausdorff distance (95th percentile), (3) absolute log-transformed volume difference, (4) sensitivity for detecting individual lesions, and (5) F1-score for individual lesions. Additionally, methods were ranked on their inter-scanner robustness. Twenty participants submitted their method for evaluation. This paper provides a detailed analysis of the results. In brief, there is a cluster of four methods that rank significantly better than the other methods, with one clear winner. The inter-scanner robustness ranking shows that not all methods generalize to unseen scanners. The challenge remains open for future submissions and provides a public platform for method evaluation. Index Terms—Magnetic resonance imaging (MRI), Brain, Evaluation and performance, Segmentation I. INTRODUCTION White matter hyperintensities (WMH) of presumed vascular origin are one of the main manifestations of cerebral small vessel disease and play a key role in stroke, dementia, and ageing [1], [2]. On T2-weighted and fluid-attenuated inversion recovery (FLAIR) brain MR images, WMH are clearly visible Copyright (c) 2019 IEEE. Personal use of this material is permitted. However, permission to use this material for any other purposes must be obtained from the IEEE by sending a request to [email protected]. 2 as hyperintense regions within the white matter [3]. An example image is shown in Figure 1, with the manual segmentation shown in Figure 1(c). Quantification of WMH is of importance in clinical research studies, where measures of WMH volume, shape, and location are obtained from detailed segmentations. These measures are associated with the presence and severity of clinical symptoms, such as cognitive impairment and gait disturbances, and are likely to find their way into daily clinical practice, supporting diagnosis, prognosis, and treatment monitoring [2], [4]. However, manual delineation of WMH is a time-consuming and observer-dependent procedure. Automatic WMH segmentation methods have been developed, but a review by Caligiuri et al. [5] revealed a key issue: it is hard to compare the various methods that are described in the literature. Each proposed segmentation method has been evaluated on a different ground truth (different number of subjects, different experts, different protocols), using different evaluation criteria. A further challenge of automatic WMH segmentation methods is the deployment of such a method within a new institute that might have different scanners or imaging protocols. Many (deep) machine learning methods require some form of transfer learning or fine-tuning on the target images [6], which in practice is not always feasible. These issues are not unique to the task of automatic WMH segmentation, but occur in many medical image analysis tasks. Organizing a scientific challenge is a way to address this, having a number of competing methods perform the same task on the same data. This has been successfully applied to various tasks, such as liver segmentation [7], image registration [8], coronary calcium scoring [9], or gland segmentation in histology images [10]. In the past, a number of challenges have been organized that included abnormalities on brain MR images, such as the multiple sclerosis (MS) lesions [11], [12], tumour [13], or tissue [14] segmentation challenges1. However, none of these challenges focuses on WMH of presumed vascular origin (although MS lesions share some characteristics with such WMH; and the brain tissue segmentation challenge included WMH lesions, but not as a separate task). The WMH Segmentation Challenge described in this paper provides a standardized assessment of automatic methods for the segmentation of WMH. The task for the challenge was defined as: “the segmentation of white matter hyperintensities of presumed vascular origin on brain MR images”2[3]. Key features of the challenge include: Participants have to submit their method to the organizers for independent evaluation on a test set. The test set includes data from two additional scanners not in the training data, to evaluate generalizability of segmentation methods across scanners. The dataset was derived from patients with various degrees of ageing related degenerative and vascular pathologies, which is important for the generalizability since segmentation methods should be able to deal with this variation. Evaluation is performed using five 1For a more complete overview, visit: https://grand-challenge.org/ challenges/ 2https://wmh.isi.uu.nl/details/ TABLE I OVERVIEW OF THE NUMBER OF IMAGES AVAILABLE FOR TRAINING (TR.) AND TEST (TE.). Institute Scanner Tr. Te. UMC Utrecht 3 T Philips Achieva 20 30 NUHS Singapore 3 T Siemens TrioTim 20 30 VU Amsterdam 3 T GE Signa HDxt 20 30 1.5 T GE Signa HDxt 0 10 3 T Philips Ingenuity (PET/MR) 0 10 different metrics and participants are ranked relative to each other. In this paper, the organization of the challenge, its results, and a detailed evaluation are presented. II. METHODS A. Training and test data A total of 60 training and 110 test images were used in this challenge. Imaging data was acquired from five different scanners, from three different vendors, in three different institutes: the University Medical Center (UMC) Utrecht, VU University Medical Centre (VU) Amsterdam, both in the Netherlands, and the National University Health System (NUHS) in Singapore. For each subject, a 3D T1-weighted and a 2D multi-slice FLAIR image were provided. The training data consisted of sixty images: twenty 3 T images of a single scanner of each institute. The test set included ninety images (three times thirty) of those same scanners and additionally twenty images (two times ten) of scanners that were not in the training data set. An overview of the data set is given in Table I. Subjects included from UMC Utrecht and VU Amsterdam were selected from the memory clinic patients of both institutes [15]. Subjects included from the NUHS Singapore were selected from the Memory Ageing and Cognition Centre Cohort recruited from the memory clinics of the National University Hospital and St. Luke’s Hospital in Singapore [16]. For each scanner, subjects were randomly picked from all subjects and randomly placed into the training or test sets. 1) MRI parameters: All 3D sequences were acquired in the sagittal direction and all 2D multi-slice sequences in the transversal direction. UMC Utrecht, 3 T Philips Achieva: 3D T1-weighted sequence (192 slices, voxel size: 1.00 ×1.00 ×1.00 mm3, repetition time (TR)/echo time (TE): 7.9/4.5ms), 2D FLAIR sequence (48 slices, voxel size: 0.96 ×0.95 ×3.00 mm3, TR/TE/inversion time (TI): 11,000/125/2,800 ms) NUHS Singapore, 3 T Siemens TrioTim: 3D T1-weighted sequence (voxel size: 1.00 ×1.00 ×1.00 mm3, TR/TE/TI: 2,300/1.9/900 ms), 2D FLAIR sequence (voxel size: 1.00 × 1.00 ×3.00 mm3, TR/TE/TI: 9,000/82/2,500 ms) VU Amsterdam, 3 T GE Signa HDxt: 3D T1-weighted sequence (176 slices, voxel size: 0.94×0.94×1.00 mm3, TR/TE: 7.8/3.0ms), 3D FLAIR sequence (132 slices, voxel size: 0.98 ×0.98 ×1.20 mm3, TR/TE/TI: 8,000/126/2,340 ms) 3 (a) T1-weighted image (b) FLAIR image (c) Manual WMH segmentation Fig. 1. Example brain MR images of a subject with white matter hyperintensities (WMH) of presumed vascular origin. On the T1-weighted image (a), WMH show as hypointense regions within the white matter. On the FLAIR image (b), WMH are clearly visible as hyperintense regions within the white matter. The corresponding manual WMH segmentation is shown in (c). VU Amsterdam, 1.5 T GE Signa HDxt: 3D T1- weighted sequence (172 slices, voxel size: 0.98 ×0.98 × 1.50 mm3, TR/TE: 12.3/5.2ms), 3D FLAIR sequence (128 slices, voxel size: 1.21 ×1.21 ×1.30 mm3, TR/TE/TI: 6,500/117/1,987 ms) VU Amsterdam, 3 T Philips Ingenuity (PET/MR): 3D T1-weighted sequence (180 slices, voxel size: 0.87 × 0.87 ×1.00 mm3, TR/TE: 9.9/4.6ms), 3D FLAIR sequence (321 slices, voxel size: 1.04 ×1.04 ×0.56 mm3, TR/TE/TI: 4,800/279/1,650 ms) All 3D FLAIR sequences were resampled into the transversal direction with slices of 3 mm thickness for two reasons: (1) to save time on the manual annotation of WMH and (2) to become more similar to the 2D multi-slice sequences. An example FLAIR image of each scanner is shown in Appendix A Figure 73. 2) Data pre-processing: All images were bias-corrected using SPM12 [17]. Using the elastix toolbox for image registration [18], the 3D T1-weighted images were aligned with the (resampled) FLAIR images. The transformation parameters were provided with the data. The faces of the subjects were manually removed from all sequences and the masks used for that were provided as well. Data before and after preprocessing is provided on the challenge website for registered participants: https://wmh.isi. uu.nl/data/. 3) Manual reference standard: WMH and other pathologies (i.e. lacunes and non-lacunar infarcts, (micro) hemorrhages) were manually segmented in accordance with the STandards for ReportIng Vascular changes on nEuroimaging (STRIVE) criteria [3]. The outline of WMH and other pathology was delineated using a contour drawing technique by an expert observer (O1). This observer had extensive prior experience with the manual segmentation of WMH and had segmented 1000+ cases before this dataset. Manual delineations were peer-reviewed by a second expert observer (O2) with eleven years of experience in quantitative neuroimaging and clinical neuroradiology. In case of mistakes, 3Available in the supplementary files / multimedia tab. errors, or delineations that were not according to the STRIVE criteria, O1 corrected the manual segmentation in a consensus meeting with O2. Hence, the provided reference standard is the corrected segmentation of O1, after peer review by O2. The contours were converted to binary masks, whereby all voxels whose volume was within the manual delineation for >50 %, were considered WMH. Background received label 0 and WMH label 1. Other pathology was converted to binary masks as well, receiving label 2. These masks were dilated by 1 pixel in-plane (with a 3×3×1voxel kernel). In case of overlap between labels 1 and 2 (after dilation), label 1 was assigned. Two additional observers segmented the sixty training images to obtain inter-observer agreement measures. Observer O3 was trained for WMH segmentation, but had no extensive prior experience. Observer O4 was trained for WMH segmentation and had prior experience. B. Set-up of the challenge Participants could register on the challenge website and download the training data. Methods had to be containerized with Docker4[19] and submitted for evaluation. Containerization eases deployment of methods and guarantees that the method will produce identical output when run on a different platform. To ensure this, the output of the containerized method on the first training subject was sent back to the participants for verification. During testing, the containerized method was run on each test subject one by one. No identifiers were present that would indicate from which of the five scanners the current subject originated. After processing a subject, the container was completely destroyed and reloaded. Full details on how the containers would be run, including a Python and MATLAB example container, were provided on the challenge website5. An NVIDIA Titan Xp GPU was available for methods that needed one. 4https://www.docker.com/ 5https://wmh.isi.uu.nl/methods/ 4 C. Participants Twenty teams submitted their method before the deadline and participated in the challenge. A brief summary of each method is given below, in alphabetical order. achilles a neural network similar to HighResNet [20] and DeepLab v3 [21], utilizing atrous (dilated) convolutions, atrous spatial pyramid pooling, and residual connections. The network is trained only on the FLAIR images, taking random 713sized patches, and applying scaling and rotation augmentations [22]. cian a network based on multi-dimensional gated recurrent units (MD-GRU) was trained on 3D patches using data augmentation techniques including random deformation, rotation and scaling [23]–[25]. hadi a random forest classifier trained on multi-modal image features. These include intensities, gradient, and Hessian features of the original images, after smoothing, and of generated super-voxels [26]. ipmi-bern a two-stage approach that uses fully convolutional neural networks to first extract the brain from the images and second identifies WMH within the brain. Both stages implement long and short skip connections. The second stage produces output at three different scales. Data augmentation was applied, including rotations and mirroring [27]. k2 a 2D fully convolutional neural network with an architecture similar to U-Net [28]. A number of models were trained for the whole dataset, as well as for each individual scanner. During application, first the type of scanner was predicted and next that specific model was applied together with the model trained on all data [29]. knight a voxel-wise logistic regression model that is fitted independently for each voxel in the FLAIR image. Images were transformed to the MNI-152 standard space [30] for training and at test time the parameter maps were warped to the subject space [31], [32]. lrde a modification of the pre-trained 16-layer VGG network [33], where the FLAIR, T1, and a high-pass filtered FLAIR are used as multi-channel input. The VGG network had its fully connected layers replaced by a number of convolutional layers [34]–[37]. misp a 3D convolutional neural network with 18 layers using patches of 27 ×27 ×9voxels. The first eight layers were trained separately for the FLAIR and T1 images and had skip-connections [38] [39]. neuro.ml a neural network using the DeepMedic [40] architecture, having two parallel branches that process the images at two different scales. The network used 3D patches, which were sampled such that 60 % of the patches contained a WMH [41]. nic-vicorob a 10-layer 3D convolutional neural network architecture previously used to segment multiple sclerosis lesions [42]. A cascaded training procedure was employed, training two separate networks to first identify candidate lesion voxels and next to reduce false positive detections. A third network re-trains the last fully connected layer to perform WMH segmentation [43]. nih cidi a fully convolutional neural network modified from the U-Net architecture [28] was used to segment WMH on the FLAIR images. Next, another network was trained to segment the white matter from T1 images, and the segmented white matter mask is applied to remove false positives from the WMH segmentation results. The original U-Net architecture was trimmed to keep only three pooling layers [44]. nist a random descision forest classifier trained on location and intensity features [45]–[47]. nlp logix a multiscale deep neural network similar to [48], with some minor modifications and no spatial features. The network was trained in ten folds and the three best performing checkpoints on the training data were selected. These were applied on the test set and the results averaged [49]. scan a densely connected convolutional network using dilated convolutions [50], [51]. In each dense block, the output is concatenated to the input before passing it to the next layer. Two classifiers were trained: one to apply brain extraction and the second to find lesions within the extracted brain [52]. skkumedneuro an intensity-based thresholding method with region growing approach to segment periventricular and deep WMH separately, and two random forest classifiers for false positive reduction. Per imaging modality, 19 texture and 100 “multi-layer” features were computed. The “multi-layer” features were computed using a feed-forward convolutional network with fixed filters (e.g. averaging, Gaussian, Laplacian); consisting of two convolutional, two max-pooling, and one fully connected layer [53]. sysu media a fully convolutional neural network similar to U-Net [28]. An ensemble of three networks was trained with different initializations. Data normalization and augmentation was applied. To remove false positive detections, WMH in the first and last 1 8th slices was removed [54], [55]. text class a random forest classifier trained primarily on texture features. Features include local binary pattern, structural and morphological gradients, and image intensities [56], [57]. tig a three-level Gaussian mixture model, slightly adapted from [58]. The model is iteratively modified and evaluated, until it converges. After that, candidate WMH is selected and possible false positives are pruned based on their location [59]. tignet a neural network with the HighResNet architecture [20]. The network was trained on 2,660 images segmented using the previous method of team tig [58], [60]. upc dlmi a neural network modified from the V-Net architecture [61]. An additional network with convolutional layers is trained on upsampled images and then concatenated with the output of the V-Net [62]. Detailed information on each method can be found online at https://wmh.isi.uu.nl/results/results-miccai-2017/. D. Evaluation and Ranking Methods were evaluated according to five criteria: (1) the Dice Similarity Coefficient (DSC), (2) a modified Hausdorff distance (95th percentile; H95), (3) the absolute percentage volume difference (AVD), (4) the sensitivity for detecting 5 individual lesions (recall), and (5) F1-score for individual lesions (F1). For recall and F1, individual lesions are defined as 3D connected components within an image. The exact implementation of each metric was put online6beforehand and could be used by participants for self-evaluation during development. During evaluation of the results, it was discovered that the AVD metric had a slight flaw. A method could undersegment WMH by at most 100%, but could oversegment WMH almost infinitely. Therefore, in this manuscript, the AVD metric was replaced by the absolute log-transformed volume difference (lAVD, (1)). lAVD =|log segmented volume true volume |(1) The final ranking was based on the five metrics and each method received a rank relative to the performance of all methods. This was computed in a number of steps. First, each metric was averaged over all test scans per method. For each metric, the methods were sorted from best to worst. Next, the best method received a rank of 0and the worst method a rank of 1; all other methods were ranked relatively in the range (0,1). Finally, the five ranks were averaged into the overall rank. 95% confidence intervals on each individual metric and the final ranking were computed using bootstrapping. The bootstrap distribution included 2,000 samples taken randomly from the test set with replacement. Non-overlapping confidence intervals indicate a significant difference between methods, with α= 0.05. It is expected that methods might have more difficulties detecting and segmenting small lesions compared to large lesions. For each subject, the recall will be computed separately for individual lesions smaller than or equal to the median lesion size and for lesions larger than the median lesion size. Additionally, a ranking was computed based solely on the inter-scanner differences. This ranking highlights which methods have the most robust performance across various scanners. For each method and scanner, the median performance of each metric was computed. Next, the standard deviation of those medians per scanner was averaged; giving a single value per metric: the standard deviation of the median per scanner. Methods were then ranked based on this value: first per metric and then averaged over all five metrics. A lower standard deviation across the median performance on all scanners (for all metrics) indicates a better inter-scanner robustness. Finally, the Simultaneous Truth And Performance Level Estimation (STAPLE) algorithm [63] was applied to all methods and to the top-ranking methods. STAPLE takes multiple segmentations as input and produces a combined segmentation, which was evaluated and ranked separately. It has been shown for other applications, e.g. brain tumour segmentation [13], that fusing the output of multiple methods can outperform all individual methods. III. RESULTS The subjects included in the challenge were (mean ±sd) 70.1±9.3years old and 50 % were male. The WMH volume 6https://github.com/hjkuijf/wmhchallenge/blob/master/evaluation.py Fig. 2. Histogram showing the WMH volume distribution throughout the dataset. The ticks on the x-axis represent each individual subject. Fig. 3. Histogram showing the WMH count distribution throughout the dataset. The ticks on the x-axis represent each individual subject. An individual lesion is defined as a 3D connected component within an image. in the dataset was (mean ±sd): 16.9±21.6ml (min: 0.78 ml, Q1: 3.24 ml, median: 11.18 ml, Q3: 23.00 ml, max: 195.15 ml; see Figure 2). The WMH count in the dataset was (mean ± sd): 62 ±35 lesions (min: 12 lesions, Q1: 36 lesions, median: 57 lesions, Q3: 81 lesions, max: 194 lesions; see Figure 3). The distribution of lesions throughout the dataset is shown in the top row of Figure 4. There were no significant differences between the training and test sets for age (p= 0.45), gender (p= 0.87), WMH volume (p= 0.74), WMH count (p= 0.75), the presence of lacunes (p= 0.86), nor for the volume of other pathology (p= 0.62). Tests for age and volumes were performed using Welch’s unequal variances t-test [64]. Tests for gender and presence of lacunes were performed using Fisher’s exact test. The bottom row of Table II shows the inter-observer agreement of observers O3 and O4 compared with the manual reference standard of the sixty training images. Additionally, the associated positions of O3 and O4 with respect to all methods in the ranking is provided. This is the position these observers would have achieved if they had participated as method in the challenge. The mean performance of each participating method on each 6 Fig. 4. The MNI-152 standard brain template [30], showing different overlays. Top row: WMH distribution throughout the dataset, where the colour indicates the percentage of subjects that have a lesion in that specific voxel. Middle row: false negative rate, showing the percentage of lesions that were missed in a specific voxel. Bottom row: false positive rate, showing the percentage of false positives in a specific voxel. All voxels where only one subject has a lesion are shown half translucent. individual metric is shown in Table II, together with the 95% confidence intervals. The method of sysu media performed best on the DSC, H95, and recall metrics. The method of cian performed best on the lAVD metric. The method of nlp logix performed best on the F1 metric. Figure 5 shows boxplots of all results of each method on each metric. The final ranking is shown in Table III, together with the 95% confidence intervals. The middle and bottom rows of Figure 4 show spatial maps of the false negative rate and false positive rate, respectively, of all methods combined. Appendix C Figures 8–277show these spatial maps per method, ordered by their final ranking. Figure 6 shows the (relative) difference in recall between small and large lesions. All methods perform worse in recalling small lesions compared to large lesions. For example, the method of sysu media has a recall for large lesions of 94% and a recall for small lesions of 76%, resulting in a relative difference of −20%. Overall, the drop in recall ranges from −20% (sysu media) to −87% (text class), as indicated by the solid lines in the figure. Table IV highlights various properties of all methods, sorted by their final ranking. The top 11 methods all employ some form of deep learning, with a U-Net-like architecture [28] being the overall most common. Amongst the non-deep learning methods, the use of a random forest classifier is most common. Almost all methods apply various kinds of pre- 7Available in the supplementary files / multimedia tab. processing, where normalizing intensities of an image to a standardized range is applied by most methods. Some methods apply post-processing techniques, mainly aimed at reducing the number of false positive detections. The H95 and F1 metrics are most sensitive to false positive detections, but the methods that apply post-processing do not have a notable better score on these metrics than nearby ranking methods. When considering only deep learning methods, most use data augmentation to generate more training samples. Scaling, rotating, and mirroring an image are quite common, but the top 2 methods also apply shearing or non-linear deformations. The last columns of Table IV highlight some properties of deep learning methods, in which a few clusters can be distinguished. Top ranking methods have applied dropout during training, some form of hard negative mining, and use an ensemble of networks. Three methods use dilated convolutions, but these cluster in the middle of the ranking. Most methods that use 3D convolutions appear to rank at the bottom. Using batch normalization, multi scale approaches, or learning rate schedules does not seem to influence the ranking; and neither does the choice of loss function. The inter-scanner robustness was determined as follows: Appendix C Figures 8–27 show the median performance of each method per metric per scanner (the line in the individual boxplots). Per metric, the standard deviation of the median values per scanner is computed. Next, methods are ranked based on those values, where a lower standard deviation 7 TABLE II MEAN PERFORMANCE AND 95 % CONFIDENCE INTERVALS OF EACH PARTICIPATING METHOD ON EACH INDIVIDUAL METRIC. METRICS INCLUDE: (1) DICE SIMILARITY COEFFICIENT (DSC), (2) MODIFIED HAUSDORFF DISTANCE (95TH PERCENTILE; H95), (3) ABSOLUTE OF THE PERCENTAGE VOLUME DIFFERENCE (AVD), (4) SENSITIVITY FOR DETECTING INDIVIDUAL LESIONS (RECALL), AND (5) F1-SCORE FOR INDIVIDUAL LESIONS (F1). BOLD INDICATES THAT A METHOD HAS THE BEST SCORE ON THAT METRIC. METHODS ARE SORTED BASED ON THE FINAL RANKING AS SHOWN IN TABLE III. THE BOTTOM ROWS INCLUDE THE RESULTS OF THE SIMULTANEOUS TRUTH AND PERFORMANCE LEVEL ESTIMATION (STAPLE) ALGORITHM APPLIED ON ALL METHODS AND ON THE TOP 4RANKING METHODS,AND THE RESULTS OF OBSERVERS O3 AND O4, TOGETHER WITH THE ASSOCIATED POSITIONS IN THE RANKING IF STAPLE, O3, AND O4 WOULD HAVE PARTICIPATED IN THE CHALLENGE. NOTE THAT O3 AND O4 SEGMENTED THE SIXTY TRAINING IMAGES. # Team DSC H95 (mm) lAVD Recall F1 1 sysu media 0.80 (0.78 - 0.82) 6.30 (4.75 - 7.93) 0.193 (0.165 - 0.224) 0.84 (0.82 - 0.86)∗0.76 (0.73 - 0.78)† 2 cian 0.78 (0.76 - 0.80) 6.82 (4.92 - 9.22) 0.193 (0.162 - 0.228) 0.83 (0.81 - 0.84)∗0.70 (0.67 - 0.73)‡ 3 nlp logix 0.77 (0.75 - 0.80) 7.16 (5.61 - 8.82) 0.219 (0.174 - 0.271) 0.73 (0.71 - 0.76)∗∗ 0.78 (0.76 - 0.80)† 4 nic-vicorob 0.77 (0.74 - 0.79) 8.28 (6.60 - 10.06) 0.248 (0.201 - 0.303) 0.75 (0.73 - 0.77)∗∗ 0.71 (0.68 - 0.73)‡ 5 k2 0.77 (0.74 - 0.79) 9.79 (7.72 - 12.28) 0.246 (0.187 - 0.310) 0.59 (0.56 - 0.61) 0.70 (0.68 - 0.72)‡ 6 misp 0.72 (0.69 - 0.75) 14.88 (10.52 - 19.41) 0.258 (0.167 - 0.388) 0.63 (0.60 - 0.65) 0.68 (0.65 - 0.70)‡ 7 lrde 0.73 (0.70 - 0.76) 14.54 (10.32 - 19.31) 0.309 (0.218 - 0.442) 0.63 (0.60 - 0.66) 0.67 (0.65 - 0.69)‡ 8 nih cidi 0.68 (0.65 - 0.70) 12.82 (10.54 - 15.16) 0.281 (0.200 - 0.394) 0.59 (0.56 - 0.62) 0.54 (0.51 - 0.57) 9 ipmi-bern 0.69 (0.67 - 0.72) 9.72 (7.98 - 11.56) 0.225 (0.178 - 0.275) 0.44 (0.42 - 0.46) 0.57 (0.55 - 0.58) 10 scan 0.63 (0.59 - 0.66) 14.34 (12.25 - 16.50) 0.277 (0.223 - 0.336) 0.55 (0.52 - 0.58) 0.51 (0.48 - 0.53) 11 achilles 0.63 (0.60 - 0.66) 11.82 (9.80 - 13.94) 0.276 (0.226 - 0.331) 0.45 (0.42 - 0.47) 0.52 (0.50 - 0.53) 12 skkumedneuro 0.58 (0.54 - 0.61) 19.02 (16.64 - 21.58) 0.384 (0.292 - 0.503) 0.47 (0.44 - 0.49) 0.51 (0.48 - 0.54) 13 tignet 0.59 (0.56 - 0.63) 21.58 (18.15 - 25.33) 0.533 (0.450 - 0.623) 0.46 (0.41 - 0.51) 0.45 (0.42 - 0.49) 14 tig 0.60 (0.56 - 0.63) 17.86 (15.57 - 20.20) 0.400 (0.333 - 0.474) 0.38 (0.36 - 0.41) 0.42 (0.40 - 0.44) 15 knight 0.70 (0.67 - 0.72) 17.03 (14.48 - 19.88) 0.352 (0.290 - 0.427) 0.25 (0.22 - 0.27) 0.35 (0.32 - 0.38) 16 upc dlmi 0.53 (0.48 - 0.58) 27.01 (22.25 - 31.99) 0.612 (0.481 - 0.762) 0.57 (0.53 - 0.60) 0.42 (0.38 - 0.46) 17 nist 0.53 (0.49 - 0.57) 15.91 (14.44 - 17.42) 0.581 (0.469 - 0.695) 0.37 (0.34 - 0.40) 0.25 (0.22 - 0.27) 18 neuro.ml 0.51 (0.45 - 0.56) 37.36 (33.70 - 40.89) 1.033 (0.836 - 1.241) 0.71 (0.68 - 0.75)∗∗ 0.21 (0.19 - 0.24) 19 text class 0.50 (0.45 - 0.54) 28.23 (24.15 - 32.68) 0.605 (0.492 - 0.724) 0.27 (0.25 - 0.29) 0.29 (0.26 - 0.31) 20 hadi 0.23 (0.19 - 0.27) 52.02 (49.25 - 54.82) 1.685 (1.448 - 1.939) 0.58 (0.52 - 0.63) 0.11 (0.09 - 0.12) 4 STAPLE (all) 0.77 (0.74 - 0.80) 5.74 (4.26 - 7.43) 0.315 (0.249 - 0.393) 0.77 (0.75 - 0.79) 0.74 (0.71 - 0.76) 2 STAPLE (top 4) 0.80 (0.78 - 0.82) 6.43 (4.48 - 8.81) 0.171 (0.144 - 0.201) 0.80 (0.78 - 0.82) 0.76 (0.74 - 0.78) 5 O3 0.77 (0.74 - 0.80) 6.79 (5.32 - 8.54) 0.176 (0.135 - 0.222) 0.65 (0.62 - 0.69) 0.74 (0.71 - 0.76) 4 O4 0.79 (0.76 - 0.81) 7.22 (5.36 - 9.36) 0.195 (0.148 - 0.245) 0.66 (0.63 - 0.70) 0.76 (0.73 - 0.78) ∗sysu media and cian perform significantly better on the recall metric than all other teams. ∗∗ nic-vicorob,nlp logix, and neuro.ml perform significantly better on the recall metric than all remaining teams. †nlp logix and sysu media perform significantly better on the F1 metric than all other teams. ‡nic-vicorob,k2,cian,misp, and lrde perform significantly better on the F1 metric than all remaining teams. indicates better inter-scanner performance. The result of this inter-scanner ranking is shown in the last column of Table III, together with the new position of that method in the ranking. The method of ipmi-bern achieves the highest inter-scanner rank and sysu media is just behind on the second rank. The methods of achilles and knight enter the top 4 of the ranking. STAPLE was applied on all methods and on the top 4 ranking methods, since these rank significantly higher than all other methods. The results are shown on the bottom rows of Table II and in Appendix C Figures 28 and 29. STAPLE on all methods would rank fourth in the challenge and achieves the best H95. STAPLE on the top 4 ranking methods would rank second in the challenge and achieves the best DSC and lAVD. When re-computing the inter-scanner robustness, both STAPLE methods outperform all other methods. STAPLE compared with the top 3 methods in the inter-scanner ranking is shown separately in Table V, because the relative ranking values change when including STAPLE. Finally, it could be hypothesized that low ranking methods suffer from training set overfitting [65] or poor generalization. This was evaluated by applying all submitted methods to the training data and comparing the performance on the training data to the performance on the test data. This analysis shows excellent correlation (R-squared: 0.94, with p < 0.001), suggesting that there is no indication for overfitting of methods on the training data. IV. DISCUSSION We have presented a standardized assessment of automatic methods for the segmentation of white matter hyperintensities of presumed vascular origin. This assessment was performed in the context of the WMH Segmentation Challenge, hosted at the 20th International Conference on Medical Image Computing and Computer Assisted Intervention (MICCAI) in 2017, Qu´ ebec City, Quebec, Canada. The manual reference standard was created in consensus by two skilled observers with extensive prior experience in WMH segmentation, which resulted in high quality WMH segmentations. Two additional observers individually segmented the sixty training images, without a consensus reading, to 8 Fig. 5. Boxplots showing all five metrics per method. The box indicates the interquartile range (IQR) with a line at the median. The whiskers extend up to 1.5 times the IQR and the fliers indicate the remaining data points. Note for the Hausdorff distance that hadi did not produce any output for 10 subjects and hence their boxplot is based on only 100 subjects (see Appendix C Figure 27 for full details). Note for the log-transformed volume difference that for visibility purposes, this figure is clipped at 3.0. Teams hadi,lrde, misp,neuro.ml,nih cidi,nist,skkumedneuro,text class, and upc dlmi have lAVD values above 3.0. For full details, see Appendix C Figures 27, 14, 13, 25, 15, 24, 19, 26, and 23, respectively. determine inter-observer agreement. The top-ranking methods achieve similar or superior performance as these two individual observers, which suggests that automatic methods might be able to replace individual observers in WMH segmentation. The moderate recall of the individual observers is mainly caused by not segmenting or missing small WMH. The F1 TABLE III FINAL RANKING OF THE METHODS THAT PARTICIPATED IN THE CHALLENGE. THE COLUMN RANK SHOWS THE RELATIVE PERFORMANCE OF EACH METHOD,BASED ON ALL FIVE METRICS LISTED IN TABLE II, TOGETHER WITH THE 95 % CONFIDENCE INTERVALS. THE COLUMN INTER-SCANNER RANK SHOWS THE RANKING WHEN IT IS COMPUTED SOLELY BASED ON INTER-SCANNER ROBUSTNESS. THE SYMBOLS BETWEEN BRACKETS INDICATE WHETHER A TEAM IS RANKED ON THE SAME POSITION (−), LOWER (∨), OR HIGHER (∧)COMPARED WITH THE ORIGINAL RANKING;WITH THE NEW POSITION INDICATED AS WELL. DOTTED LINES INDICATE CLUSTERS OF METHODS THAT RANK SIGNIFICANTLY DIFFERENT FROM METHODS RANKED ABOVE/BELOW, BECAUSE OF NON-OVERLAPPING CONFIDENCE INTERVALS. # Team Rank (95 % CI) Inter-scanner rank 1 sysu media 0.0068 (0.0019 - 0.0161)†0.0375 (∨2) 2 cian 0.0357 (0.0248 - 0.0539)‡0.0831 (∨5) 3 nlp logix 0.0520 (0.0365 - 0.0744)‡0.1111 (∨7) 4 nic-vicorob 0.0785 (0.0577 - 0.1045)‡0.1629 (∨11) 5 k2 0.1437 (0.1188 - 0.1711) 0.1174 (∨8) 6 misp 0.1740 (0.1356 - 0.2273) 0.1915 (∨12) 7 lrde 0.1782 (0.1395 - 0.2290) 0.3510 (∨17) 8 nih cidi 0.2376 (0.2131 - 0.2680) 0.1570 (∨10) 9 ipmi-bern 0.2537 (0.2391 - 0.2727) 0.0345 (∧1) 10 scan 0.2836 (0.2631 - 0.3099) 0.2252 (∨14) 11 achilles 0.3058 (0.2896 - 0.3276) 0.0714 (∧3) 12 skkumedneuro 0.3649 (0.3325 - 0.4044) 0.1105 (∧6) 13 tignet 0.4090 (0.3765 - 0.4481) 0.2969 (∨15) 14 tig 0.4097 (0.3795 - 0.4454) 0.1289 (∧9) 15 knight 0.4320 (0.4082 - 0.4598) 0.0785 (∧4) 16 upc dlmi 0.4429 (0.3903 - 0.5016) 0.7415 (∨20) 17 nist 0.5040 (0.4724 - 0.5404) 0.3052 (∧16) 18 neuro.ml 0.5615 (0.5193 - 0.6084) 0.6110 (∨19) 19 text class 0.5961 (0.5539 - 0.6430) 0.2117 (∧13) 20 hadi 0.8886 (0.8687 - 0.9103) 0.4974 (∧18) †sysu media ranks significantly higher than all other participants. ‡cian, nlp logix, and nic-vicorob rank significantly higher than all remaining participants. is higher than the recall, which is opposite for most automatic methods, and indicates that both O3 and O4 have hardly segmented any false positive WMH. The organizers have chosen not to disclose the test set, contrary to what is common in medical image analysis challenges. By keeping the test set secret, a high reliability of the results can be ensured because it obviates the possibility of (visual) self-evaluation by participants. The rapidly increasing popularity of (deep) neural networks as methodology of choice for analysing medical images [66] is noticeable in this challenge as well. Fourteen of the twenty submitted methods employ some form of (deep) neural networks, including all methods in the top ten. Nevertheless, the use of deep learning methodology is not a guaranteed recipe for success, since a number of low-ranking methods use it as well. Ensemble methods appear to do very well in this challenge. The methods of sysu media (# 1), nlp logix (# 3), and k2 (# 5) use an ensemble of separately trained neural networks to achieve top-ranking results. Furthermore, the STAPLE algo- 9 TABLE IV OVERVIEW OF VARIOUS PROPERTIES OF ALL METHODS. METHODS ARE SORTED BASED ON THE FINAL RANKING AS SHOWN IN TABLE III. Neural network features5 # Team Pre1Method Post2Data DL3Aug4Loss function Dim Dil BN Drop MS LR HN Ens 1 sysu media I,RU-Net SL T,FXH,R,SDSC 2D X 2 cian F,IMDGRU T,FXD,R,Smultinom. log. 2D†X 3 nlp logix I,SCNN T,FXcross-entropy 2D X X X X X 4 nic-vicorob CNN SM T,FXM,Rcross-entropy 3D X X X 5 k2 I,R,SU-Net T,FXMDSC 2D X X X 6 misp I,RCNN T,FXmean sq. error 3D X X X 7 lrde F,IVGG-16 T,FXR,Smultinom. log. 2D X 8 nih cidi SU-Net G T,FXM,Rcross-entropy 2D X X 9 ipmi-bern I,SU-Net T,FXM,Rcross-entropy 2D X X X 10 scan SDenseNet T,FXcross-entropy 2D X 11 achilles I,RHighResNet FXR,SDSC 3D X X X 12 skkumedneuro I,SRF T,F 13 tignet B,I,THighResNet T,F∗XDSC 3D X X X 14 tig B,S,TGMM FP T,F 15 knight B,I,S,TVLR SM F∗M,T,YDSC 16 upc dlmi IU-Net T,FXMDSC 3D X X X 17 nist B,I,TRF T,F 18 neuro.ml DeepMedic T,FXcross-entropy 3D X X 19 text class I,RRF SM T,F 20 hadi RF T,F 1Pre-processing: B= bias field correction, F= morphological filter to enhance small lesions, I= intensity normalization, R= resizing or resampling to a predefined grid, S= skull stripping, and T= transformation to a standard space. 2Post-processing: FP location based false positive reduction, Ggraph-based segmentation refinement, SLremove slices prone to false positives, and SM remove small segmentation results. 3DL indicates whether this method uses deep learning. 4Augmentation of training data: D= non-linear deforming, H= shearing, M= mirroring, R= rotating, S= scaling, T= translating/moving, and Y= generating synthetic lesions. 5Features used in the neural networks. Dim: 2D or 3D convolutions. Dil: dilated convolutions. BN: batch normalisation. Drop: dropout. MS: multi scale approaches (e.g. separate paths at different resolutions). LR: use of a learing rate schedule (e.g. reducing the learning rate during training). HN: hard negative mining. Ens: an ensemble of multiple networks. ∗additional data from other sources was used to train this method. †the convolutions are 2D, but the third dimension is processed within an RNN that incorporates all dimensions. TABLE V THE RE-COMPUTED RESULTS OF THE INTER-SCANNER ROBUSTNESS RANKING WHEN INCLUDING THE SIMULTANEOUS TRUTH AND PERFORMANCE LEVEL ESTIMATION (STAPLE) ALGORITHM APPLIED ON ALL METHODS OR ON THE TOP 4RANKING METHODS. STAPLE OUTPERFORMS ALL METHODS AND THEREFORE THE RELATIVE RANKING VALUES CHANGE. HERE, STAPLE IS COMPARED TO THE TOP 3METHODS IN THE ORIGINAL INTER-SCANNER RANKING IN TABLE III. THE SYMBOLS BETWEEN BRACKETS INDICATE WHETHER A TEAM IS RANKED ON THE SAME POSITION (−), LOWER (∨), OR HIGHER (∧)COMPARED TO THE ORIGINAL RANKING;WITH THE NEW POSITION INDICATED AS WELL. #Team Inter-scanner rank STAPLE (top 4) 0.0152 (1) STAPLE (all) 0.0390 (1) 1 ipmi-bern 0.0400 (∨2) 0.0402 (∨2) 2 sysu media 0.0433 (∨3) 0.0434 (∨3) 3 achilles 0.0769 (∨4) 0.0768 (∨4) rithm that combines all methods or the top 4 ranking methods achieves good results as well. On the inter-scanner robustness ranking, both results of the STAPLE algorithm outperform all other participating methods. Combining the results of various methods has been performed in other challenges as well, for example in the brain tumour segmentation challenge (BRATS) [13]. However, in that challenge the combination of methods always outperformed all individual methods, whereas in this challenge the method of sysu media remains the winner. This seems to be mainly caused by the good performance of sysu media on the recall metric compared to the STAPLE results. Both STAPLE methods perform less well in recalling small lesions below the median size, as can be seen in Figure 6. The use of dropout during training is another characteristic of top-ranking methods. Random dropouts prevent units in neural networks from co-adapting too much [67] and introduces some redundancy in the network. A larger network trained with dropout might behave like an ensemble of smaller networks; and ensemble methods also rank at the top. However, the deep learning methods trained with dropout have a considerably lower inter-scanner rank (Table III). They drop more in the inter-scanner ranking than methods trained without dropout, suggesting that these methods might not generalize well to unseen data from unseen scanners; but instead only to APPENDIX C SUMMARIES OF RESULTS Detailed summaries of all results per participant are given in the following appendices. All figures show the performance of each method on the five scanners described in Section II-A1 for the following five criteria: (1) the Dice Similarity Coefficient (DSC), (2) a modified Hausdorff distance (95th percentile; H95), (3) the absolute log-transformed volume difference (lAVD), (4) the sensitivity for detecting individual lesions (recall), and (5) F1-score for individual lesions (F1). Next to that are two columns that show spatial maps of the false negative rate (left) and false positive rate (right). The figures are presented in the order of the final ranking, as shown in Table III. Fig. 8. Detailed results of sysu media. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 9. Detailed results of cian. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 10. Detailed results of nlp logix. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 11. Detailed results of nic-vicorob. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 12. Detailed results of k2. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the H95 boxplot contains values out of range (max = 108.18 mm). Fig. 13. Detailed results of misp. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the H95 and lAVD boxplots contain values out of range (max H95 = 111.09 mm; max lAVD = 5.27). Fig. 14. Detailed results of lrde. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the H95 and lAVD boxplots contain values out of range (max H95 = 118.57 mm; max lAVD = 5.60). Fig. 15. Detailed results of nih cidi. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the lAVD boxplot contains values out of range (max = 5.22). Fig. 16. Detailed results of ipmi-bern. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 17. Detailed results of scan. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 24. Detailed results of nist. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the lAVD boxplot contains values out of range (max = 3.02). Fig. 25. Detailed results of neuro.ml. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the lAVD boxplot contains values out of range (max = 4.80). Fig. 26. Detailed results of text class. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: the H95 and lAVD boxplots contain values out of range (max H95 = 152.66 mm; max lAVD = 3.48). Fig. 27. Detailed results of hadi. The two columns on the right show the false negative rate (left) and false positive rate (right). Note: some output was empty and the H95 and lAVD were not evaluated, and the lAVD boxplot contains values out of range (max = 4.83). Fig. 28. Detailed results of STAPLE applied on all methods. The two columns on the right show the false negative rate (left) and false positive rate (right). Fig. 29. Detailed results of STAPLE applied on the top 4 ranking methods. The two columns on the right show the false negative rate (left) and false positive rate (right).