Int J CARS DOI 10.1007/s11548-013-0961-0
ORIGINAL ARTICLE
PNist: interactive volumetric measurements of plexiform neurofibromas in MRI scans Lior Weizman · Dina Helfer · Dafna Ben Bashat · Li-tal Pratt · Leo Joskowicz · Shlomi Constantini · Ben Shofty · Liat Ben Sira
Received: 20 June 2013 / Accepted: 4 November 2013 © CARS 2013
Abstract Purpose Volumetric measurements of plexiform neurofibromas (PNs) are time consuming and error prone, as they require the delineation of the PN boundaries, which is mostly impractical in the daily clinical setup. Accurate volumetric measurements are seldom performed for these tumors mainly due to their great dispersion, size and multiple locations. This paper presents a semiautomatic method for segmentation of PN from STIR MRI scans. Methods Plexiform neurofibroma interactive segmentation tool (PNist) is a new tool to segment PNs in STIR MRI scans. Electronic supplementary material The online version of this article (doi:10.1007/s11548-013-0961-0) contains supplementary material, which is available to authorized users. L. Weizman (B) · D. Helfer · L. Joskowicz School of Engineering and Computer Science and The Edmond and Lily Safra Center for Brain Sciences (ELSC), The Hebrew University of Jerusalem, Jerusalem, Israel e-mail: [email protected] D. Helfer e-mail: [email protected] L. Joskowicz e-mail: [email protected] D. Ben Bashat · L. Pratt · S. Constantini · B. Shofty · L. Ben Sira Tel Aviv Medical Center, Tel Aviv, Israel e-mail: [email protected] L. Pratt e-mail: [email protected] S. Constantini e-mail: [email protected] B. Shofty e-mail: [email protected] L. Ben Sira e-mail: [email protected]
The method is based on histogram tumor models computed from a training set. Results Experimental results from 28 datasets show an average absolute volume difference of 6.8 % with an average user time of approximately 7 min versus more than 13 min with manual delineation. In complex cases, the PNist user time is less than half in compared to state-of-the-art tools. Conclusions PNist is a new method for the semiautomatic segmentation of PN lesions. Its simplicity and reliability make it unique among other state-of-the-art methods. It has the potential to become a clinical tool that allows the reliable evaluation of PN burden and progression.
Keywords Segmentation · Neurofibromatosis-1 · Plexiform neurofibromas · PNS tumors · Shwannoma
Introduction Tumor volume estimation in MRI is a key indicator for disease progression and response to treatment [9,14]. The most common measurement of radiological follow-up and response to treatment are based on linear measurements of the tumor defined by WHO, RECIST and MacDonald [2,7,14]. The disadvantages of these criteria are that they disregard the three-dimensional information available in modern imaging techniques, as well as suffering from limited reproducibility [11]. Furthermore, some tumors involve multiple parts of the body and have asymmetric growth, making them unsuitable for linear measurements. In these cases, the major factor for decision-making is the tumor burden, which is estimated based on volumetric quantification of the tumor. Manual delineation is too time consuming and labor intensive to be practical. As a result, the tumor burden is
123
Int J CARS
“eyeballed” for comparative purposes in most cases. This may compromise the assessment quality of clinical course and treatment decisions. Therefore, automatic and semiautomatic methods for tumor segmentation using medical image processing are essential to make volumetric measurements clinically feasible. In this paper, we focus on plexiform neurofibromas (PNs). These tumors are one of the hallmarks of neurofibromatosis1 (NF1) [5], a common genetic disorder that involves the central and peripheral nervous systems. The reported incidence of PN ranges between 30 and 40 % of NF1 patients. This is probably an underestimation since internal PNs are often undetected without appropriate imaging [8,16,19]. These tumors have a significant size range, but typically they are large and extensive, with complex irregular shapes. Their growth rate is erratic and unpredictable. They can involve different parts of the body and may infiltrate, displace or compress the surrounding structures [8,19]. Thus, in addition to aesthetic disfiguration, they may lead to substantial morbidity [8,12]. In 10 % of cases, these lesions transform into malignant peripheral nerve sheath tumors (MPNSTs) [16,19]. This serious complication is associated with higher tumor burden and rapid tumor progression. This makes tumor burden assessment essential for detecting aggressive lesions at an early stage and for monitoring response to therapy [10]. The complex shape and the asymmetric growth of PNs make it impossible to assess tumor growth reliably with unior bi-dimensional measurements. Figure 1 shows a pelvis
short TI inversion recovery (STIR) MR scan with PN lesions. Volumetric measurements of PN lesions with a generalpurpose medical segmentation tool may take more than 30 min of user time. This is impractical if we take into account that radiologists typically have less than 5 min to read an MR scan, which consists of tens to hundreds of 2D slices [3]. Very few automatic and semiautomatic methods that address the challenging task of PN segmentation have been reported in the literature. Solomon et al. [13] describe an interactive method for volume quantification of PNs. Their algorithm is a histogram-based method that requires an initial user definition of the PNs in each slice without transfer of information from adjacent slices. Cai et al. [1] present a dynamic threshold (DT) level set-based method for 3D nerve sheath tumors segmentation in MR images. It requires user interaction for each lesion and expands automatically from the user’s predefined location to segment the entire PN lesion according to statistical measures in the 3D image space. Note that both methods assume that tumor pixels are brighter than the surrounding tissue in MRI STIR scans. In addition, they require intensive, continued and timeconsuming expert user interaction. In our recent publication [18], we presented a method for PN segmentation. It requires an initial delineation of the tumor area in a single slice and automatically segments the PN lesions in the entire scan based on connectivity analysis. However, that method’s performance is sub-optimal for very complex PN lesions in large images, as these may contain many disconnected PN lesions. In this paper, we present a user friendly method that locates the relatively bright PN tumor pixels in each MRI slice based on very simple and easily generated input. The simplicity of the required input makes our method unique among other state-of-the-art algorithms. Based on the method, we developed a tool, called PN Interactive Segmentation Tool (PNist), which is currently in use for PN segmentation in the clinical environment. Methods
Fig. 1 STIR coronal cross-section image with PN lesions contours marked in red
123
Our method consists of a training phase, performed once per tumor type, and an interactive phase in which a new scan is segmented. Figure 2 illustrates the flow of our method. The input to the training phase is a set of STIR MR scans and manual annotation of the PNs in those scans. These scans and annotations are used to create models of expected intensity distribution in the vicinity of tumors, which are then stored in a histogram database. Information from the database, together with a user scribble on a new scan, is then fed into the interactive phase of the method to segment PNs in the new STIR MRI scans.
Int J CARS Fig. 2 The flow chart of PNist. After a training phase based on samples of manually annotated PNs in STIR MRIs, a database of PN histogram models is established. PNs in a new STIR MR scan are segmented based on simple user scribble in conjunction with information from the database
Training phase The input for the training phase is a set of MRI STIR scans, each with a ground truth segmentation of tumor voxels, annotated by an expert. The output is a set of vectors that represent histograms of PN lesion environments. The training phase consists of two steps: (1) tumor sampling and (2) model histogram generation. 1. Tumor sampling In each image slice, every tumor component is fitted with two axis aligned frames: a large frame as the bounding box of the tumor region and a small frame with sides half the size of the large frame, located in its center. The majority of the small frame covers tumor pixels, and the large frame covers both tumor and normal surrounding tissue samples. These frames are used to model the gray-level intensities of tumor and normal regions. In the interactive segmentation phase (“Interactive segmentation phase” section), the user is required to draw a scribble over the tumor region to provide a seed for the segmentation process. In the training phase, we use both large and small frames to simulate two common scenarios of scribble drawing. The small frame histogram is used to simulate cases in which the user’s scribble is within the tumor area
and contains no healthy tissues. The large frame histogram is used to simulate cases in which the user’s scribble contains, besides tumor pixels, also healthy regions. These two histogram patterns are sufficient to cover the most common cases in scribble drawing. Figure 3a, b illustrates large and small frames overlaid on a STIR MRI slice. Each frame is then assigned with a threshold corresponding to 90 % on the cumulative histogram of the PN pixels in the frame. The value of 90 % was determined after a few trials with several other values on our datasets and was found to be the most appropriate to compensate for partial volume effect and minor manual ground truth segmentation errors. A 50-bin histogram is computed for each frame, by dividing the range between the minimum and the maximum gray levels in the frame into 50 bins. This results in two histograms per tumor component. While a higher number of bins is expected to provide better results, the value of 50 was found to be the lowest value that enables a short computation time and still provides reliable results. These histograms, along with their associated threshold bins, are the key for automatic threshold selection and tumor segmentation in the interactive segmentation phase. Figure 3c, d illustrates two histograms of the PN pixels in the frames (the regions that are marked red in Fig. 3a, b). Note that the threshold computation is based on
123
Int J CARS
(a)
(c)
(b)
Histogram of PN pixles in large frame
(d)
Histogram of PN pixles in small frame
0.06
Prevalence
Prevalence
0.06 0.04
0.02
0
0
10
20
30
40
0.04
0.02
0
50
0
10
20
Bin# Histogram of pixles in large frame
50
0.1 0.08
0.06
0.04
0.02
0
40
Histogram of pixles in small frame
(f)
0.08
Prevalence
Prevalence
(e)
30
Bin#
0.06 0.04 0.02
0
10
20
30
40
50
Bin#
0
0
10
20
30
40
50
Bin#
Fig. 3 MRI slice showing a large and b small frames around a tumor component. c, d The corresponding 50-bin histograms of the PN pixles only in the large (c) and small (d) frames. e, f The corresponding 50-bin histograms of the entire large (c) and small (d) frames. The histogram of the large frame shows two peaks (e), as expected from a sample that contains both tumor and background pixels. The histogram of the
small frame shows only one peak (f) since it is made almost entirely of tumor pixels. The threshold bin, computed as the brightest value below 90 % of the brightest PN pixels in the frame, is shown in red in c, d. This threshold is translated to the gray-level range of the entire frames’ histograms and is shown in red in e, f
these histograms. Figure 3e, f shows the histograms of the entire large and small frames.
normalizing the histogram vectors. This involves dividing all bins by the total number of pixels such that their sum is one, making the scale of the vectors independent of the number of pixels in the sample. Then, we cluster the normalized vectors into k groups with the k-means [15] algorithm where similarity between
2. Model histogram generation In this step, we create representative models of different tumor environments. First, we strip the histograms of all image-specific information by
123
Int J CARS
Fig. 4 Fifteen model histograms obtained after the training phase. These histograms and their associated thresholds (marked in red) serve as a database of model histograms which will be used to segment PNs in future STIR MR scans, based on a simple user scribble
two normalized histogram vectors is measured with Earth mover’s distance (EMD) [6] metric. We call the mean vector of each cluster a “model histogram.” Similarly, we name the mean threshold bin number associated with the histograms in the cluster “threshold bin number.” Note that the clustering process does not take into account the source of each normalized vector, in terms of small/large frame histogram. As a result, the information regarding small/large frame histogram does not exist in the averaged mean histogram. At the end of this process, we obtain k model histograms, each with an associated threshold bin number. These histograms and thresholds serve as a database of model histograms which will be used in the next phase. This learning process is performed only once, and the resulting model is saved and can be used henceforth for tumor segmentation on new MRI STIR scans. Figure 4 shows representative histogram models and their associated thresholds for k = 15.
models from the previous section; (2) an unseen MRI STIR with PN lesions; and (3) a user scribble consisting of one or more PN lesions in the scan. This interactive phase consists of five steps: (1) User input; (2) Classification and thresholding; (3) Connected component filtering; (4) Manual Error Correction and; (5) Multislice expansion in 3D. Figure 5 illustrates the stages of the interactive phase. These steps are repeated until the user determines that all tumor pixels in the scan are marked. 1. User input The input is a scribble of pixels on a STIR MRI slice with PN lesions (Fig. 5a). It is required that the sample includes at least some tumor pixels in addition to background. Note that a sample does not have to include the entire tumor and can include pixels from several lesions. At this point, the user should use his clinical knowledge to avoid including bright non-tumor structures in the sample—if they are sampled, they may be misclassified as tumors.
Interactive segmentation phase Once the training is complete, we segment tumors in scans of new patients based on a very simple user scribble. The inputs to this phase are as follows: (1) the database of histogram
2. Classification and thresholding We compute the 50-bin normalized histogram of the manually drawn pixels and search the database of histogram models for the closest matching model histogram with the EMD similarity metric
123
Int J CARS
(a)
(b)
(d)
(c)
(e)
Fig. 5 Illustration of the interactive segmentation phase stages. a A user scribble is used as an input to the segmentation algorithm; b Classification. First, the scribble histogram is calculated (top). Then, the histogram from the histograms database (shown in Fig. 4) that best matches the scribble histogram is found (middle). Finally, the associated threshold of the matching histogram is converted to the scribble’s gray-level dynamic range (bottom). c Thresholding. The converted threshold is
used to create a binary image of bright pixels in the slice. d Connected Component filtering. Only the bright pixels that are connected to the user’s scribble are filtered, to obtain the PN segmentation in the slice. This stage is followed by manual exclusion of bright non-tumor structures, which are connected to the PN lesions. e Multi-slice expansion in 3D—the segmentation from one slice can be copied to adjacent slices and used as input for the segmentation of those slices
(Fig. 5b). Once the most similar model histogram is found, we use it to compute a threshold that applies to the new sample as follows: First, we convert the model’s threshold bin number to an intensity value in the scan’s intensity values range. We then use this value to threshold the unseen MRI slice. This process results in a binary mask that indicates the pixels within the intensity range of the tumors (Fig. 5c). We then apply erosion followed by dilation operators [17], both with a structuring element matrix of 3 × 3, to eliminate noise and “holes” in areas that have the potential to be marked as tumors.
3. Connected component filtering The result of the previous step may misclassify healthy tissues with intensity range similar to that of tumors. To remove misclassified pixels, we filter the connected components in the mask result, adding only components that intersect with the original sample to the segmentation. This approach labels all the tumors that were sampled even if they belong to more than one mass. Note that this approach operates automatically to minimize the user interaction time of the method.
123
4. Manual error correction Inevitably, there are cases in which hyper-intensity healthy tissues are connected to the PN
Int J CARS
lesions and therefore can be wrongly segmented as PNs. To cope with this phenomenon, the user can manually perform error corrections to eliminate healthy, wrongly segmented areas in the slice. Figure 5d shows a PN segmentation result after this manual correction. 5. Multi-slice expansion The user scribble on a single slice can be used to segment tumor in adjacent slices. The user may choose to operate in a “3D mode” by defining the maximal slice distance to effect. Once the segmentation of a tumor has been determined in one slice, the result is copied onto the adjacent slices and used as input. The output result is then successively used as input for adjacent slices. The process stops when the maximal distance of slices has been reached or when the computed segmentation threshold is lower than the slice average, indicating the area does not contain tumors. While it is possible to perform the connected components filtering step in 3D mode to obviate the multi-slice expansion, we limit our method to a 2D slice-by-slice analysis due to the nature of the data. In many cases, PN STIR body scans consist of thick image slices, (our datasets include scans with slice thickness of up to 9 mm). In such cases, the scribble used for the segmentation of a single slice may not be appropriate for adjacent slices due to the large gap between the slices. In addition, tumor connectivity might be interrupted for the same reason. Therefore, we decided to perform the connectivity analysis and the segmentation in a 2D mode, while allowing the option to propagate the segmentation into adjacent slices. The choice of how many adjacent slices will be used in 3D mode is made by the user to avoid misclassification of healthy pixels as tumors, depending on the physical dimension of the data and the nature of the tumor in the case at hand.
Experimental results We developed a graphical user interface tool to allow convenient use of our method. It allows viewing and creation of new segmentation based on training, as previously described. A snapshot of PNist is shown in Fig. 6. We retrospectively obtained 32 MRI datasets from 16 NF1 patients with PN lesions. Four datasets with a total of 686 PN lesion regions were used for training, and 28 datasets were used to validate the method. The MR images were acquired by 1.5 Tesla MR system (GE Signa EXCITE HDx, Milwaukee, WI, USA) at the Tel Aviv Sourasky Medical Center. Images were acquired in coronal or axial planes. Slice resolution is 256 × 256, 512 × 512 or 576 × 576. The number of slices in a scan varies between 14 and 48. The voxel size in the datasets varies between 0.4 × 0.4 × 3.3 mm3 and 1.9×1.9×9 mm3 . TR values vary between 2020 and 6880 ms, TE values vary between 26 and 120 ms and TI values vary
between 140 and 150 ms. Tumor volume ranges from 5.11 to 713 cm3 . The scans show lesions in various locations: scalp, neck, shoulder, spine, abdomen, pelvis and calves. To characterize the diversity of PN appearance, tumors were partitioned into three categories according to their characteristics. The partitioning was performed visually, based on radiological assessment as follows: – Simple—lesion has not extended or infiltrated into adjacent tissue and its borders are well-defined. – Intermediate—extension or infiltration of adjacent tissue with clear borders. – Complex—extension or infiltration of adjacent tissue with unclear borders. An expert radiologist manually annotated PNs in the four training images. The histogram vectors of the training phase were clustered into k = 15 clusters. Twenty-eight scans from 12 patients were used to evaluate the method’s performance. These patients were not included in the training phase. Three means of evaluation were used to quantify the performance of the method: 1. Repeatability: We test how the tumor volume varies from one use of PNist to the next, on the same image. The intraobserver variability was measured, where each scan was repeatedly segmented by the same user two or three times on different days. 2. Reliability: We compared the volume of PNs segmentation generated with PNist to gold-standard segmentation obtained with a general-purpose volumetric segmentation tool (Analyze Direct Software, version 9.0, Mayo Clinic, Rochester, MN, USA), which is currently considered as the state-of-the-art medical image processing software. 3. Segmentation time: We compared the time required to complete segmentation with PNist and Analyze Direct Software. The metrics used for evaluation are absolute volume difference (AVD), volume overlap error (VOE) and average absolute surface distance (AASD), computed as follows: S (1) AVD = 100 × − 1 T 2 × |S T | VOE = 1 − × 100 (2) |S| + |T | where S and T refer to the two segmented volumes that are being compared. Gerig et al. [4] provide detailed explanation and computation methodology of AASD. The segmentation time is measured in minutes and consists of the time required to complete the segmentation, including error correction actions and computer computation time.
123
Int J CARS Fig. 6 Annotated screenshot of PN Interactive Segmentation Tool (PNist). The display area shows a slice of an abdominal STIR MRI scan with the tumors segmented in red. The program allows viewing, creating and editing medical segmentations. The function of each control is listed next to it
Maximal Intra−observer Absolute Volume Difference
Repeatability evaluation
Reliability evaluation The PN tumors in 28 data sets were segmented separately with PNist and Analyze. The results were compared to measure the accuracy of our tool. Figure 9 shows the correlation between volumes measured by the different tools (R 2 = 0.994). For small tumors (volume < 20 cm3 ), the correlation is lower, with R 2 = 0.8371. The average AVD is
123
Maximal Intra−observer AVD (%)
To quantify the intra-observer variability, 28 MRI scans were segmented using PNist two or three times on different days by the same radiologist. We consider the mean volume of all measurements as the best estimate for the real tumor volume. We define the maximal AVD for each scan as the highest AVD between each measurement and the mean volume for the scan. Figure 7 shows the maximal AVD for each scan. The average maximal AVD is 3 %. The VOE was computed for every pair of segmentations (three segmentations yield three different pairs). The average VOE for each scan is shown in Fig. 8. The standard deviation for each scan is less than 3.2 %. The average of all intraobserver VOE is 9.5 %. Note that the VOE is less than 10 % in all cases that were classified as simple and that the VOE for the complex cases is close to the average value. This indicates that the segmentation of the data set can be reproduced with low variability, even in complex cases.
12 Simple Intermediate Complex
10
8
6
4
2
0
0
200
400
600
800
3
Tumor volume (cm ) Fig. 7 Repeatability using PNist: for each dataset, the maximal AVD between segmentations performed on different days and their mean value is shown
6.8 %. The average of all VOE is 12.5 %. Note that complex cases also exhibit high correlation with the gold-standard segmentation, indicating the reliability of PNist in these cases. Figure 10 shows the average absolute surface distance (AASD) between the segmentation results of PNist and Analyze versus the tumor volume. The average of all AASD is 0.6451 mm. Note that as opposed to volume-based measures, AASD does not depend on the overall tumor volume.
Int J CARS
Average Absolute Surface Distance Average Absolute Surface Distance (mm)
Intra−observer Volume Overlap Error (%)
Average Intra−observer Volume Overlap Error 25 Simple Intermediate Complex
20
15
10
5
0
0
200
400
600
800
2 Simple Intermediate Complex 1.5
1
0.5
0
0
200
700 600
R2=0.994 500 400 300 Simple Intermediate Complex
0
100
200
300
400
500
600 3
Segmentation time: PNist vs. Analyze 40 35 30 25 20 15 Simple
10
Intermediate Complex
5
y=x 0
0
5
10
15
20
25
30
35
40
PNist segmentation time (minutes)
y=x 0
800
Fig. 10 Reliability using PNist: average absolute surface distance (AASD) between the segmentation results of PNist and Analyze
Analyze segmentation time (minutes)
Analyze segmentation volume (cm3)
PNist volume vs. Analyze volume
100
600
Tumor volume (cm )
Fig. 8 Repeatability using PNist: average VOE for each scan. Each data point represents the average overlap error between all possible pairs of measurements of the same data set. The simple group has significantly lower average overlap value than the other classes (6.0 % compared with 11.9 and 10.2 % for intermediate and complex, respectively)
200
400 3
Tumor volume (cm3)
700
PNist segmentation volume (cm )
Fig. 11 Segmentation times of PNist versus Analyze
Fig. 9 Correlation between PNist and Analyze measured volumes
Segmentation time Figure 11 shows the segmentation times of PNist versus Analyze. An improvement in segmentation time was achieved regardless of tumor volume or difficulty in all but two cases. Most cases show a significant improvement in segmentation time. The average user time of PNist is 7:12 versus 13:48 min with Analyze. Compared with the estimated segmentation time of Analyze, the required segmentation time using our method is on average 52 % less. The repeatability of PNist is similar to other PN segmentation methods in the literature [12,13,18]. However,
in terms of reliability and user time, PNist performs significantly better. For example, PNist provides an AVD of 6.8 % obtained with an average user time of less than 8 min. This is a significant improvement in comparison with the method of Solomon et al. [13] that exhibits an AVD of 37 %, as previously reported [18]. Note that both methods (ours and Solomon et al.) require pre- or post-processing to eliminate high intensity healthy regions which are connected to PNs. In addition, Solomon et al. report a lower AVD, 5.6 %, when their method is applied on a different dataset [13]. This can be explained by two reasons: (a) our database includes more tumor complexities than the cases tested in [13] and (b) AVD in [13] was calculated as the difference between segmentation results obtained by two observers, while we
123
Int J CARS
report difference between results of the automatic method and a gold standard obtained manually. PNist also outperforms our previously published method [18], which reports an AVD of 10 % obtained with an average user time of 13 min. The VOE results emphasize the superiority of PNist which exhibits a VOE of 12.5 %, in comparison with 43.4 % reported by using Solomon et al. and 26.9 % with our previous method [18].
Discussion We have presented a new method for the semiautomatic segmentation of PN lesions based on STIR MR images. In this section, we discuss two important issues regarding the training procedure required for reliable performance of the method and further approaches to automate the method. Training samples optimization Our method is based on model histograms, generated from training samples of PN lesions in STIR scans. The training samples should reflect the properties of the scans that will be segmented using PNist, in terms of gray-level distribution of lesions and other tissues in the image, the organ/body regions appearing in the scans and the complexity of the lesions. In our experiments, we used four training scans, which included 686 training lesion regions. While we did not analyze the minimal/maximal number of training scans required for adequate results, we believe that it is important to include various lesions at different locations and complexities and to include different acquisition parameters in the training database. Therefore, several training scans with small numbers of lesions each are preferred over a single training scan with multiple lesions, as a single scan may limit the method’s compatibility to a small range of acquisition parameters. Our database consists of scans acquired from multiple scanners of the same manufacturer. In cases where scans are acquired from scanners of multiple manufacturers, the training database may need to be extended to include the variability of STIR images provided by various manufacturers. Further approaches for automation PNist requires human interaction in drawing the scribble over the PN lesions in the scans and performing error corrections as necessary. This interaction can be reduced by considering the following approaches to extend the method: – Multi-resolution template: By defining moving window templates with varying sizes, a single image slice can be automatically searched for bright regions at various sizes.
123
This may obviate the user scribble, as bright PN lesions will be automatically detected and marked by the moving window. In this approach, however, a post-processing step will be required to exclude bright healthy regions in the image which may be wrongly segmented as tumors. – Histogram matching: The current method does not normalize the scans’ intensity values but rather normalizes the histogram vectors of the training samples and the scribble histogram only. Histogram normalization of the entire scans before the segmentation may improve the method’s results and reduce the human interaction in the error correction step. – FLAIR integration: Bright regions in STIR images include, besides PN lesions, also fluid and cyst areas. These areas, however, are dark in FLAIR images. Registration of the same scan’s STIR and FLAIR images can highly reduce the false-positive errors of our method, by eliminating fluid and cyst areas that appear bright in STIR images and wrongly classified as PN lesions.
Conclusions The assessment of changes in tumor volume is a key factor in treatment decisions. Accurate volumetric measurements are rarely performed for tumors of plexiform neurofibroma type, mainly due to complex shape, size and multiple lesion locations. These characteristics make current segmentation methods time consuming and prone to large variability. Lacking a better alternative, treatment decisions are largely based on a gross assessment of changes in tumor using MRI. PNist is a new method for the semiautomatic segmentation of PN lesions based on STIR MR images. Its simplicity and reliability make it unique among state-of-the-art methods. Existing semiautomatic segmentation methods usually rely on region of interest delineation which requires more time and a steady hand or seed initiation which can be tiresome when many small disconnected masses are involved. Our method relies on a user scribble over the object (or objects) of interest. It is performed in under a second, requiring less mental exertion and accuracy than other input types. This input contains sufficient information to locate the tumor and distinguish between the intensity levels of tumors and the surrounding tissues. The characterization of tumor gray levels based on similar normalized histograms is novel in this type of application. This method is not limited to PN tumor detection—the algorithm can be adapted to detect any type of anatomical structure with a distinct range of gray levels in relation to its background. All that is required is a different model histogram pool. This means that our method can be potentially used to generate fast and easy to use semiautomatic segmentation tools for other anatomies.
Int J CARS
Future work will contemplate including an existing baseline segmentation to improve the repeatability and user time of follow-up PN segmentation of the same patient. Acknowledgments The authors wish to thank the Gilbert Israeli Neurofibromatosis Center (GINFC) for their contribution of providing the real data and supporting the medical part of the paper. The authors also thank Vicki Myers for editorial assistance. All procedures followed were in accordance with the ethical standards of the responsible committee on human experimentation (institutional and national) and with the Helsinki Declaration of 1975, as revised in 2008 (5). Informed consent was obtained from all patients for being included in the study. Conflict of interest None.
References 1. Cai W, Kassarjian A, Bredella MA, Harris GJ, Yoshida H, Mautner VF, Wenzel R, Plotkin SR (2009) Tumor burden in patients with neurofibromatosis types 1 and 2 and schwannomatosis: determination on whole-body MR images1. Radiology 250(3):665–673 2. Eisenhauer E, Therasse P, Bogaerts J, Schwartz L, Sargent D, Ford R, Dancey J, Arbuck S, Gwyther S, Mooney M et al (2009) New response evaluation criteria in solid tumours: revised recist guideline (version 1.1). Eur J Cancer 45(2):228–247 3. Fleishon HB, Bhargavan M, Meghea C et al (2006) Radiologists’ reading times using pacs and using films: one practice’s experience. Acad Radiol 13(4):453 4. Gerig G, Jomier M, Chakos M (2001) Valmet: a new validation tool for assessing and improving 3d object segmentation. In: Medical image computing and computer-assisted intervention-MICCAI 2001, Springer, pp 516–523 5. Huson S, Hughes R (1994) The neurofibromatoses: a pathogenetic and clinical overview. Chapman & Hall, London 6. Levina E, Bickel P (2001) The earth mover’s distance is the mallows distance: some insights from statistics. In: Eighth IEEE international conference on computer vision, 2001. ICCV 2001. Proceedings. IEEE, vol 2, pp 251–256 7. Macdonald DR, Cascino TL, Schold S, Cairncross JG (1990) Response criteria for phase ii studies of supratentorial malignant glioma. J Clin Oncol 8(7):1277–1280
8. Mautner V, Hartmann M, Kluwe L, Friedrich R, Fünsterer C (2006) MRI growth patterns of plexiform neurofibromas in patients with neurofibromatosis type 1. Neuroradiology 48(3):160–165 9. Miller A, Hoogstraten B, Staquet M, Winkler A (1981) Reporting results of cancer treatment. Cancer 47(1):207–214 10. Nguyen R, Kluwe L, Fuensterer C, Kentsch M, Friedrich RE, Mautner VF (2011) Plexiform neurofibromas in children with neurofibromatosis type 1: frequency and associated clinical deficits. J Pediatr 159(4):652–655 11. Park JO, Lee SI, Song SY, Kim K, Kim WS, Jung CW, Park YS, Im YH, Kang WK, Lee MH et al (2003) Measuring response in solid tumors: comparison of recist and who response criteria. Jpn J Clin Oncol 33(10):533–537 12. Poussaint TY, Jaramillo D, Chang Y, Korf B (2003) Interobserver reproducibility of volumetric MR imaging measurements of plexiform neurofibromas. Am J Roentgenol 180(2):419–423 13. Solomon J, Warren K, Dombi E, Patronas N, Widemann B (2004) Automated detection and volume measurement of plexiform neurofibromas in neurofibromatosis 1 using magnetic resonance imaging. Comput Med Imaging Gr 28(5):257–265 14. Therasse P, Arbuck SG, Eisenhauer EA, Wanders J, Kaplan RS, Rubinstein L, Verweij J, Van Glabbeke M, van Oosterom AT, Christian MC et al (2000) New guidelines to evaluate the response to treatment in solid tumors. J Natl Cancer Inst 92(3):205–216 15. Tou J, Gonzalez RC (1974) Pattern recognition principles. Addison-Wesley Publishing Company, Reading, MA 16. Tucker T, Friedman J, Friedrich R, Wenzel R, Fnsterer C, Mautner V (2009) MRI growth patterns of plexiform neurofibromas in patients with neurofibromatosis type 1. J Med Genet 46(2):81–85 17. Van Den Boomgaard R, Van Balen R (1992) Methods for fast morphological image transforms using bitmapped binary images. CVGIP Gr Models Image Process 54(3):252–258 18. Weizman L, Hoch L, Bashat DB, Joskowicz L, Constantini S, Sira LB (2012) Interactive segmentation of plexiform neurofibroma tissue: method and preliminary performance evaluation. Med Biol Eng Comput 50(8):877–884 19. Williams VC, Lucas J, Babcock MA, Gutmann DH, Korf B, Maria BL (2009) Neurofibromatosis type 1 revisited. Pediatrics 123(1):124–133
123
ORIGINAL ARTICLE
PNist: interactive volumetric measurements of plexiform neurofibromas in MRI scans Lior Weizman · Dina Helfer · Dafna Ben Bashat · Li-tal Pratt · Leo Joskowicz · Shlomi Constantini · Ben Shofty · Liat Ben Sira
Received: 20 June 2013 / Accepted: 4 November 2013 © CARS 2013
Abstract Purpose Volumetric measurements of plexiform neurofibromas (PNs) are time consuming and error prone, as they require the delineation of the PN boundaries, which is mostly impractical in the daily clinical setup. Accurate volumetric measurements are seldom performed for these tumors mainly due to their great dispersion, size and multiple locations. This paper presents a semiautomatic method for segmentation of PN from STIR MRI scans. Methods Plexiform neurofibroma interactive segmentation tool (PNist) is a new tool to segment PNs in STIR MRI scans. Electronic supplementary material The online version of this article (doi:10.1007/s11548-013-0961-0) contains supplementary material, which is available to authorized users. L. Weizman (B) · D. Helfer · L. Joskowicz School of Engineering and Computer Science and The Edmond and Lily Safra Center for Brain Sciences (ELSC), The Hebrew University of Jerusalem, Jerusalem, Israel e-mail: [email protected] D. Helfer e-mail: [email protected] L. Joskowicz e-mail: [email protected] D. Ben Bashat · L. Pratt · S. Constantini · B. Shofty · L. Ben Sira Tel Aviv Medical Center, Tel Aviv, Israel e-mail: [email protected] L. Pratt e-mail: [email protected] S. Constantini e-mail: [email protected] B. Shofty e-mail: [email protected] L. Ben Sira e-mail: [email protected]
The method is based on histogram tumor models computed from a training set. Results Experimental results from 28 datasets show an average absolute volume difference of 6.8 % with an average user time of approximately 7 min versus more than 13 min with manual delineation. In complex cases, the PNist user time is less than half in compared to state-of-the-art tools. Conclusions PNist is a new method for the semiautomatic segmentation of PN lesions. Its simplicity and reliability make it unique among other state-of-the-art methods. It has the potential to become a clinical tool that allows the reliable evaluation of PN burden and progression.
Keywords Segmentation · Neurofibromatosis-1 · Plexiform neurofibromas · PNS tumors · Shwannoma
Introduction Tumor volume estimation in MRI is a key indicator for disease progression and response to treatment [9,14]. The most common measurement of radiological follow-up and response to treatment are based on linear measurements of the tumor defined by WHO, RECIST and MacDonald [2,7,14]. The disadvantages of these criteria are that they disregard the three-dimensional information available in modern imaging techniques, as well as suffering from limited reproducibility [11]. Furthermore, some tumors involve multiple parts of the body and have asymmetric growth, making them unsuitable for linear measurements. In these cases, the major factor for decision-making is the tumor burden, which is estimated based on volumetric quantification of the tumor. Manual delineation is too time consuming and labor intensive to be practical. As a result, the tumor burden is
123
Int J CARS
“eyeballed” for comparative purposes in most cases. This may compromise the assessment quality of clinical course and treatment decisions. Therefore, automatic and semiautomatic methods for tumor segmentation using medical image processing are essential to make volumetric measurements clinically feasible. In this paper, we focus on plexiform neurofibromas (PNs). These tumors are one of the hallmarks of neurofibromatosis1 (NF1) [5], a common genetic disorder that involves the central and peripheral nervous systems. The reported incidence of PN ranges between 30 and 40 % of NF1 patients. This is probably an underestimation since internal PNs are often undetected without appropriate imaging [8,16,19]. These tumors have a significant size range, but typically they are large and extensive, with complex irregular shapes. Their growth rate is erratic and unpredictable. They can involve different parts of the body and may infiltrate, displace or compress the surrounding structures [8,19]. Thus, in addition to aesthetic disfiguration, they may lead to substantial morbidity [8,12]. In 10 % of cases, these lesions transform into malignant peripheral nerve sheath tumors (MPNSTs) [16,19]. This serious complication is associated with higher tumor burden and rapid tumor progression. This makes tumor burden assessment essential for detecting aggressive lesions at an early stage and for monitoring response to therapy [10]. The complex shape and the asymmetric growth of PNs make it impossible to assess tumor growth reliably with unior bi-dimensional measurements. Figure 1 shows a pelvis
short TI inversion recovery (STIR) MR scan with PN lesions. Volumetric measurements of PN lesions with a generalpurpose medical segmentation tool may take more than 30 min of user time. This is impractical if we take into account that radiologists typically have less than 5 min to read an MR scan, which consists of tens to hundreds of 2D slices [3]. Very few automatic and semiautomatic methods that address the challenging task of PN segmentation have been reported in the literature. Solomon et al. [13] describe an interactive method for volume quantification of PNs. Their algorithm is a histogram-based method that requires an initial user definition of the PNs in each slice without transfer of information from adjacent slices. Cai et al. [1] present a dynamic threshold (DT) level set-based method for 3D nerve sheath tumors segmentation in MR images. It requires user interaction for each lesion and expands automatically from the user’s predefined location to segment the entire PN lesion according to statistical measures in the 3D image space. Note that both methods assume that tumor pixels are brighter than the surrounding tissue in MRI STIR scans. In addition, they require intensive, continued and timeconsuming expert user interaction. In our recent publication [18], we presented a method for PN segmentation. It requires an initial delineation of the tumor area in a single slice and automatically segments the PN lesions in the entire scan based on connectivity analysis. However, that method’s performance is sub-optimal for very complex PN lesions in large images, as these may contain many disconnected PN lesions. In this paper, we present a user friendly method that locates the relatively bright PN tumor pixels in each MRI slice based on very simple and easily generated input. The simplicity of the required input makes our method unique among other state-of-the-art algorithms. Based on the method, we developed a tool, called PN Interactive Segmentation Tool (PNist), which is currently in use for PN segmentation in the clinical environment. Methods
Fig. 1 STIR coronal cross-section image with PN lesions contours marked in red
123
Our method consists of a training phase, performed once per tumor type, and an interactive phase in which a new scan is segmented. Figure 2 illustrates the flow of our method. The input to the training phase is a set of STIR MR scans and manual annotation of the PNs in those scans. These scans and annotations are used to create models of expected intensity distribution in the vicinity of tumors, which are then stored in a histogram database. Information from the database, together with a user scribble on a new scan, is then fed into the interactive phase of the method to segment PNs in the new STIR MRI scans.
Int J CARS Fig. 2 The flow chart of PNist. After a training phase based on samples of manually annotated PNs in STIR MRIs, a database of PN histogram models is established. PNs in a new STIR MR scan are segmented based on simple user scribble in conjunction with information from the database
Training phase The input for the training phase is a set of MRI STIR scans, each with a ground truth segmentation of tumor voxels, annotated by an expert. The output is a set of vectors that represent histograms of PN lesion environments. The training phase consists of two steps: (1) tumor sampling and (2) model histogram generation. 1. Tumor sampling In each image slice, every tumor component is fitted with two axis aligned frames: a large frame as the bounding box of the tumor region and a small frame with sides half the size of the large frame, located in its center. The majority of the small frame covers tumor pixels, and the large frame covers both tumor and normal surrounding tissue samples. These frames are used to model the gray-level intensities of tumor and normal regions. In the interactive segmentation phase (“Interactive segmentation phase” section), the user is required to draw a scribble over the tumor region to provide a seed for the segmentation process. In the training phase, we use both large and small frames to simulate two common scenarios of scribble drawing. The small frame histogram is used to simulate cases in which the user’s scribble is within the tumor area
and contains no healthy tissues. The large frame histogram is used to simulate cases in which the user’s scribble contains, besides tumor pixels, also healthy regions. These two histogram patterns are sufficient to cover the most common cases in scribble drawing. Figure 3a, b illustrates large and small frames overlaid on a STIR MRI slice. Each frame is then assigned with a threshold corresponding to 90 % on the cumulative histogram of the PN pixels in the frame. The value of 90 % was determined after a few trials with several other values on our datasets and was found to be the most appropriate to compensate for partial volume effect and minor manual ground truth segmentation errors. A 50-bin histogram is computed for each frame, by dividing the range between the minimum and the maximum gray levels in the frame into 50 bins. This results in two histograms per tumor component. While a higher number of bins is expected to provide better results, the value of 50 was found to be the lowest value that enables a short computation time and still provides reliable results. These histograms, along with their associated threshold bins, are the key for automatic threshold selection and tumor segmentation in the interactive segmentation phase. Figure 3c, d illustrates two histograms of the PN pixels in the frames (the regions that are marked red in Fig. 3a, b). Note that the threshold computation is based on
123
Int J CARS
(a)
(c)
(b)
Histogram of PN pixles in large frame
(d)
Histogram of PN pixles in small frame
0.06
Prevalence
Prevalence
0.06 0.04
0.02
0
0
10
20
30
40
0.04
0.02
0
50
0
10
20
Bin# Histogram of pixles in large frame
50
0.1 0.08
0.06
0.04
0.02
0
40
Histogram of pixles in small frame
(f)
0.08
Prevalence
Prevalence
(e)
30
Bin#
0.06 0.04 0.02
0
10
20
30
40
50
Bin#
0
0
10
20
30
40
50
Bin#
Fig. 3 MRI slice showing a large and b small frames around a tumor component. c, d The corresponding 50-bin histograms of the PN pixles only in the large (c) and small (d) frames. e, f The corresponding 50-bin histograms of the entire large (c) and small (d) frames. The histogram of the large frame shows two peaks (e), as expected from a sample that contains both tumor and background pixels. The histogram of the
small frame shows only one peak (f) since it is made almost entirely of tumor pixels. The threshold bin, computed as the brightest value below 90 % of the brightest PN pixels in the frame, is shown in red in c, d. This threshold is translated to the gray-level range of the entire frames’ histograms and is shown in red in e, f
these histograms. Figure 3e, f shows the histograms of the entire large and small frames.
normalizing the histogram vectors. This involves dividing all bins by the total number of pixels such that their sum is one, making the scale of the vectors independent of the number of pixels in the sample. Then, we cluster the normalized vectors into k groups with the k-means [15] algorithm where similarity between
2. Model histogram generation In this step, we create representative models of different tumor environments. First, we strip the histograms of all image-specific information by
123
Int J CARS
Fig. 4 Fifteen model histograms obtained after the training phase. These histograms and their associated thresholds (marked in red) serve as a database of model histograms which will be used to segment PNs in future STIR MR scans, based on a simple user scribble
two normalized histogram vectors is measured with Earth mover’s distance (EMD) [6] metric. We call the mean vector of each cluster a “model histogram.” Similarly, we name the mean threshold bin number associated with the histograms in the cluster “threshold bin number.” Note that the clustering process does not take into account the source of each normalized vector, in terms of small/large frame histogram. As a result, the information regarding small/large frame histogram does not exist in the averaged mean histogram. At the end of this process, we obtain k model histograms, each with an associated threshold bin number. These histograms and thresholds serve as a database of model histograms which will be used in the next phase. This learning process is performed only once, and the resulting model is saved and can be used henceforth for tumor segmentation on new MRI STIR scans. Figure 4 shows representative histogram models and their associated thresholds for k = 15.
models from the previous section; (2) an unseen MRI STIR with PN lesions; and (3) a user scribble consisting of one or more PN lesions in the scan. This interactive phase consists of five steps: (1) User input; (2) Classification and thresholding; (3) Connected component filtering; (4) Manual Error Correction and; (5) Multislice expansion in 3D. Figure 5 illustrates the stages of the interactive phase. These steps are repeated until the user determines that all tumor pixels in the scan are marked. 1. User input The input is a scribble of pixels on a STIR MRI slice with PN lesions (Fig. 5a). It is required that the sample includes at least some tumor pixels in addition to background. Note that a sample does not have to include the entire tumor and can include pixels from several lesions. At this point, the user should use his clinical knowledge to avoid including bright non-tumor structures in the sample—if they are sampled, they may be misclassified as tumors.
Interactive segmentation phase Once the training is complete, we segment tumors in scans of new patients based on a very simple user scribble. The inputs to this phase are as follows: (1) the database of histogram
2. Classification and thresholding We compute the 50-bin normalized histogram of the manually drawn pixels and search the database of histogram models for the closest matching model histogram with the EMD similarity metric
123
Int J CARS
(a)
(b)
(d)
(c)
(e)
Fig. 5 Illustration of the interactive segmentation phase stages. a A user scribble is used as an input to the segmentation algorithm; b Classification. First, the scribble histogram is calculated (top). Then, the histogram from the histograms database (shown in Fig. 4) that best matches the scribble histogram is found (middle). Finally, the associated threshold of the matching histogram is converted to the scribble’s gray-level dynamic range (bottom). c Thresholding. The converted threshold is
used to create a binary image of bright pixels in the slice. d Connected Component filtering. Only the bright pixels that are connected to the user’s scribble are filtered, to obtain the PN segmentation in the slice. This stage is followed by manual exclusion of bright non-tumor structures, which are connected to the PN lesions. e Multi-slice expansion in 3D—the segmentation from one slice can be copied to adjacent slices and used as input for the segmentation of those slices
(Fig. 5b). Once the most similar model histogram is found, we use it to compute a threshold that applies to the new sample as follows: First, we convert the model’s threshold bin number to an intensity value in the scan’s intensity values range. We then use this value to threshold the unseen MRI slice. This process results in a binary mask that indicates the pixels within the intensity range of the tumors (Fig. 5c). We then apply erosion followed by dilation operators [17], both with a structuring element matrix of 3 × 3, to eliminate noise and “holes” in areas that have the potential to be marked as tumors.
3. Connected component filtering The result of the previous step may misclassify healthy tissues with intensity range similar to that of tumors. To remove misclassified pixels, we filter the connected components in the mask result, adding only components that intersect with the original sample to the segmentation. This approach labels all the tumors that were sampled even if they belong to more than one mass. Note that this approach operates automatically to minimize the user interaction time of the method.
123
4. Manual error correction Inevitably, there are cases in which hyper-intensity healthy tissues are connected to the PN
Int J CARS
lesions and therefore can be wrongly segmented as PNs. To cope with this phenomenon, the user can manually perform error corrections to eliminate healthy, wrongly segmented areas in the slice. Figure 5d shows a PN segmentation result after this manual correction. 5. Multi-slice expansion The user scribble on a single slice can be used to segment tumor in adjacent slices. The user may choose to operate in a “3D mode” by defining the maximal slice distance to effect. Once the segmentation of a tumor has been determined in one slice, the result is copied onto the adjacent slices and used as input. The output result is then successively used as input for adjacent slices. The process stops when the maximal distance of slices has been reached or when the computed segmentation threshold is lower than the slice average, indicating the area does not contain tumors. While it is possible to perform the connected components filtering step in 3D mode to obviate the multi-slice expansion, we limit our method to a 2D slice-by-slice analysis due to the nature of the data. In many cases, PN STIR body scans consist of thick image slices, (our datasets include scans with slice thickness of up to 9 mm). In such cases, the scribble used for the segmentation of a single slice may not be appropriate for adjacent slices due to the large gap between the slices. In addition, tumor connectivity might be interrupted for the same reason. Therefore, we decided to perform the connectivity analysis and the segmentation in a 2D mode, while allowing the option to propagate the segmentation into adjacent slices. The choice of how many adjacent slices will be used in 3D mode is made by the user to avoid misclassification of healthy pixels as tumors, depending on the physical dimension of the data and the nature of the tumor in the case at hand.
Experimental results We developed a graphical user interface tool to allow convenient use of our method. It allows viewing and creation of new segmentation based on training, as previously described. A snapshot of PNist is shown in Fig. 6. We retrospectively obtained 32 MRI datasets from 16 NF1 patients with PN lesions. Four datasets with a total of 686 PN lesion regions were used for training, and 28 datasets were used to validate the method. The MR images were acquired by 1.5 Tesla MR system (GE Signa EXCITE HDx, Milwaukee, WI, USA) at the Tel Aviv Sourasky Medical Center. Images were acquired in coronal or axial planes. Slice resolution is 256 × 256, 512 × 512 or 576 × 576. The number of slices in a scan varies between 14 and 48. The voxel size in the datasets varies between 0.4 × 0.4 × 3.3 mm3 and 1.9×1.9×9 mm3 . TR values vary between 2020 and 6880 ms, TE values vary between 26 and 120 ms and TI values vary
between 140 and 150 ms. Tumor volume ranges from 5.11 to 713 cm3 . The scans show lesions in various locations: scalp, neck, shoulder, spine, abdomen, pelvis and calves. To characterize the diversity of PN appearance, tumors were partitioned into three categories according to their characteristics. The partitioning was performed visually, based on radiological assessment as follows: – Simple—lesion has not extended or infiltrated into adjacent tissue and its borders are well-defined. – Intermediate—extension or infiltration of adjacent tissue with clear borders. – Complex—extension or infiltration of adjacent tissue with unclear borders. An expert radiologist manually annotated PNs in the four training images. The histogram vectors of the training phase were clustered into k = 15 clusters. Twenty-eight scans from 12 patients were used to evaluate the method’s performance. These patients were not included in the training phase. Three means of evaluation were used to quantify the performance of the method: 1. Repeatability: We test how the tumor volume varies from one use of PNist to the next, on the same image. The intraobserver variability was measured, where each scan was repeatedly segmented by the same user two or three times on different days. 2. Reliability: We compared the volume of PNs segmentation generated with PNist to gold-standard segmentation obtained with a general-purpose volumetric segmentation tool (Analyze Direct Software, version 9.0, Mayo Clinic, Rochester, MN, USA), which is currently considered as the state-of-the-art medical image processing software. 3. Segmentation time: We compared the time required to complete segmentation with PNist and Analyze Direct Software. The metrics used for evaluation are absolute volume difference (AVD), volume overlap error (VOE) and average absolute surface distance (AASD), computed as follows: S (1) AVD = 100 × − 1 T 2 × |S T | VOE = 1 − × 100 (2) |S| + |T | where S and T refer to the two segmented volumes that are being compared. Gerig et al. [4] provide detailed explanation and computation methodology of AASD. The segmentation time is measured in minutes and consists of the time required to complete the segmentation, including error correction actions and computer computation time.
123
Int J CARS Fig. 6 Annotated screenshot of PN Interactive Segmentation Tool (PNist). The display area shows a slice of an abdominal STIR MRI scan with the tumors segmented in red. The program allows viewing, creating and editing medical segmentations. The function of each control is listed next to it
Maximal Intra−observer Absolute Volume Difference
Repeatability evaluation
Reliability evaluation The PN tumors in 28 data sets were segmented separately with PNist and Analyze. The results were compared to measure the accuracy of our tool. Figure 9 shows the correlation between volumes measured by the different tools (R 2 = 0.994). For small tumors (volume < 20 cm3 ), the correlation is lower, with R 2 = 0.8371. The average AVD is
123
Maximal Intra−observer AVD (%)
To quantify the intra-observer variability, 28 MRI scans were segmented using PNist two or three times on different days by the same radiologist. We consider the mean volume of all measurements as the best estimate for the real tumor volume. We define the maximal AVD for each scan as the highest AVD between each measurement and the mean volume for the scan. Figure 7 shows the maximal AVD for each scan. The average maximal AVD is 3 %. The VOE was computed for every pair of segmentations (three segmentations yield three different pairs). The average VOE for each scan is shown in Fig. 8. The standard deviation for each scan is less than 3.2 %. The average of all intraobserver VOE is 9.5 %. Note that the VOE is less than 10 % in all cases that were classified as simple and that the VOE for the complex cases is close to the average value. This indicates that the segmentation of the data set can be reproduced with low variability, even in complex cases.
12 Simple Intermediate Complex
10
8
6
4
2
0
0
200
400
600
800
3
Tumor volume (cm ) Fig. 7 Repeatability using PNist: for each dataset, the maximal AVD between segmentations performed on different days and their mean value is shown
6.8 %. The average of all VOE is 12.5 %. Note that complex cases also exhibit high correlation with the gold-standard segmentation, indicating the reliability of PNist in these cases. Figure 10 shows the average absolute surface distance (AASD) between the segmentation results of PNist and Analyze versus the tumor volume. The average of all AASD is 0.6451 mm. Note that as opposed to volume-based measures, AASD does not depend on the overall tumor volume.
Int J CARS
Average Absolute Surface Distance Average Absolute Surface Distance (mm)
Intra−observer Volume Overlap Error (%)
Average Intra−observer Volume Overlap Error 25 Simple Intermediate Complex
20
15
10
5
0
0
200
400
600
800
2 Simple Intermediate Complex 1.5
1
0.5
0
0
200
700 600
R2=0.994 500 400 300 Simple Intermediate Complex
0
100
200
300
400
500
600 3
Segmentation time: PNist vs. Analyze 40 35 30 25 20 15 Simple
10
Intermediate Complex
5
y=x 0
0
5
10
15
20
25
30
35
40
PNist segmentation time (minutes)
y=x 0
800
Fig. 10 Reliability using PNist: average absolute surface distance (AASD) between the segmentation results of PNist and Analyze
Analyze segmentation time (minutes)
Analyze segmentation volume (cm3)
PNist volume vs. Analyze volume
100
600
Tumor volume (cm )
Fig. 8 Repeatability using PNist: average VOE for each scan. Each data point represents the average overlap error between all possible pairs of measurements of the same data set. The simple group has significantly lower average overlap value than the other classes (6.0 % compared with 11.9 and 10.2 % for intermediate and complex, respectively)
200
400 3
Tumor volume (cm3)
700
PNist segmentation volume (cm )
Fig. 11 Segmentation times of PNist versus Analyze
Fig. 9 Correlation between PNist and Analyze measured volumes
Segmentation time Figure 11 shows the segmentation times of PNist versus Analyze. An improvement in segmentation time was achieved regardless of tumor volume or difficulty in all but two cases. Most cases show a significant improvement in segmentation time. The average user time of PNist is 7:12 versus 13:48 min with Analyze. Compared with the estimated segmentation time of Analyze, the required segmentation time using our method is on average 52 % less. The repeatability of PNist is similar to other PN segmentation methods in the literature [12,13,18]. However,
in terms of reliability and user time, PNist performs significantly better. For example, PNist provides an AVD of 6.8 % obtained with an average user time of less than 8 min. This is a significant improvement in comparison with the method of Solomon et al. [13] that exhibits an AVD of 37 %, as previously reported [18]. Note that both methods (ours and Solomon et al.) require pre- or post-processing to eliminate high intensity healthy regions which are connected to PNs. In addition, Solomon et al. report a lower AVD, 5.6 %, when their method is applied on a different dataset [13]. This can be explained by two reasons: (a) our database includes more tumor complexities than the cases tested in [13] and (b) AVD in [13] was calculated as the difference between segmentation results obtained by two observers, while we
123
Int J CARS
report difference between results of the automatic method and a gold standard obtained manually. PNist also outperforms our previously published method [18], which reports an AVD of 10 % obtained with an average user time of 13 min. The VOE results emphasize the superiority of PNist which exhibits a VOE of 12.5 %, in comparison with 43.4 % reported by using Solomon et al. and 26.9 % with our previous method [18].
Discussion We have presented a new method for the semiautomatic segmentation of PN lesions based on STIR MR images. In this section, we discuss two important issues regarding the training procedure required for reliable performance of the method and further approaches to automate the method. Training samples optimization Our method is based on model histograms, generated from training samples of PN lesions in STIR scans. The training samples should reflect the properties of the scans that will be segmented using PNist, in terms of gray-level distribution of lesions and other tissues in the image, the organ/body regions appearing in the scans and the complexity of the lesions. In our experiments, we used four training scans, which included 686 training lesion regions. While we did not analyze the minimal/maximal number of training scans required for adequate results, we believe that it is important to include various lesions at different locations and complexities and to include different acquisition parameters in the training database. Therefore, several training scans with small numbers of lesions each are preferred over a single training scan with multiple lesions, as a single scan may limit the method’s compatibility to a small range of acquisition parameters. Our database consists of scans acquired from multiple scanners of the same manufacturer. In cases where scans are acquired from scanners of multiple manufacturers, the training database may need to be extended to include the variability of STIR images provided by various manufacturers. Further approaches for automation PNist requires human interaction in drawing the scribble over the PN lesions in the scans and performing error corrections as necessary. This interaction can be reduced by considering the following approaches to extend the method: – Multi-resolution template: By defining moving window templates with varying sizes, a single image slice can be automatically searched for bright regions at various sizes.
123
This may obviate the user scribble, as bright PN lesions will be automatically detected and marked by the moving window. In this approach, however, a post-processing step will be required to exclude bright healthy regions in the image which may be wrongly segmented as tumors. – Histogram matching: The current method does not normalize the scans’ intensity values but rather normalizes the histogram vectors of the training samples and the scribble histogram only. Histogram normalization of the entire scans before the segmentation may improve the method’s results and reduce the human interaction in the error correction step. – FLAIR integration: Bright regions in STIR images include, besides PN lesions, also fluid and cyst areas. These areas, however, are dark in FLAIR images. Registration of the same scan’s STIR and FLAIR images can highly reduce the false-positive errors of our method, by eliminating fluid and cyst areas that appear bright in STIR images and wrongly classified as PN lesions.
Conclusions The assessment of changes in tumor volume is a key factor in treatment decisions. Accurate volumetric measurements are rarely performed for tumors of plexiform neurofibroma type, mainly due to complex shape, size and multiple lesion locations. These characteristics make current segmentation methods time consuming and prone to large variability. Lacking a better alternative, treatment decisions are largely based on a gross assessment of changes in tumor using MRI. PNist is a new method for the semiautomatic segmentation of PN lesions based on STIR MR images. Its simplicity and reliability make it unique among state-of-the-art methods. Existing semiautomatic segmentation methods usually rely on region of interest delineation which requires more time and a steady hand or seed initiation which can be tiresome when many small disconnected masses are involved. Our method relies on a user scribble over the object (or objects) of interest. It is performed in under a second, requiring less mental exertion and accuracy than other input types. This input contains sufficient information to locate the tumor and distinguish between the intensity levels of tumors and the surrounding tissues. The characterization of tumor gray levels based on similar normalized histograms is novel in this type of application. This method is not limited to PN tumor detection—the algorithm can be adapted to detect any type of anatomical structure with a distinct range of gray levels in relation to its background. All that is required is a different model histogram pool. This means that our method can be potentially used to generate fast and easy to use semiautomatic segmentation tools for other anatomies.
Int J CARS
Future work will contemplate including an existing baseline segmentation to improve the repeatability and user time of follow-up PN segmentation of the same patient. Acknowledgments The authors wish to thank the Gilbert Israeli Neurofibromatosis Center (GINFC) for their contribution of providing the real data and supporting the medical part of the paper. The authors also thank Vicki Myers for editorial assistance. All procedures followed were in accordance with the ethical standards of the responsible committee on human experimentation (institutional and national) and with the Helsinki Declaration of 1975, as revised in 2008 (5). Informed consent was obtained from all patients for being included in the study. Conflict of interest None.
References 1. Cai W, Kassarjian A, Bredella MA, Harris GJ, Yoshida H, Mautner VF, Wenzel R, Plotkin SR (2009) Tumor burden in patients with neurofibromatosis types 1 and 2 and schwannomatosis: determination on whole-body MR images1. Radiology 250(3):665–673 2. Eisenhauer E, Therasse P, Bogaerts J, Schwartz L, Sargent D, Ford R, Dancey J, Arbuck S, Gwyther S, Mooney M et al (2009) New response evaluation criteria in solid tumours: revised recist guideline (version 1.1). Eur J Cancer 45(2):228–247 3. Fleishon HB, Bhargavan M, Meghea C et al (2006) Radiologists’ reading times using pacs and using films: one practice’s experience. Acad Radiol 13(4):453 4. Gerig G, Jomier M, Chakos M (2001) Valmet: a new validation tool for assessing and improving 3d object segmentation. In: Medical image computing and computer-assisted intervention-MICCAI 2001, Springer, pp 516–523 5. Huson S, Hughes R (1994) The neurofibromatoses: a pathogenetic and clinical overview. Chapman & Hall, London 6. Levina E, Bickel P (2001) The earth mover’s distance is the mallows distance: some insights from statistics. In: Eighth IEEE international conference on computer vision, 2001. ICCV 2001. Proceedings. IEEE, vol 2, pp 251–256 7. Macdonald DR, Cascino TL, Schold S, Cairncross JG (1990) Response criteria for phase ii studies of supratentorial malignant glioma. J Clin Oncol 8(7):1277–1280
8. Mautner V, Hartmann M, Kluwe L, Friedrich R, Fünsterer C (2006) MRI growth patterns of plexiform neurofibromas in patients with neurofibromatosis type 1. Neuroradiology 48(3):160–165 9. Miller A, Hoogstraten B, Staquet M, Winkler A (1981) Reporting results of cancer treatment. Cancer 47(1):207–214 10. Nguyen R, Kluwe L, Fuensterer C, Kentsch M, Friedrich RE, Mautner VF (2011) Plexiform neurofibromas in children with neurofibromatosis type 1: frequency and associated clinical deficits. J Pediatr 159(4):652–655 11. Park JO, Lee SI, Song SY, Kim K, Kim WS, Jung CW, Park YS, Im YH, Kang WK, Lee MH et al (2003) Measuring response in solid tumors: comparison of recist and who response criteria. Jpn J Clin Oncol 33(10):533–537 12. Poussaint TY, Jaramillo D, Chang Y, Korf B (2003) Interobserver reproducibility of volumetric MR imaging measurements of plexiform neurofibromas. Am J Roentgenol 180(2):419–423 13. Solomon J, Warren K, Dombi E, Patronas N, Widemann B (2004) Automated detection and volume measurement of plexiform neurofibromas in neurofibromatosis 1 using magnetic resonance imaging. Comput Med Imaging Gr 28(5):257–265 14. Therasse P, Arbuck SG, Eisenhauer EA, Wanders J, Kaplan RS, Rubinstein L, Verweij J, Van Glabbeke M, van Oosterom AT, Christian MC et al (2000) New guidelines to evaluate the response to treatment in solid tumors. J Natl Cancer Inst 92(3):205–216 15. Tou J, Gonzalez RC (1974) Pattern recognition principles. Addison-Wesley Publishing Company, Reading, MA 16. Tucker T, Friedman J, Friedrich R, Wenzel R, Fnsterer C, Mautner V (2009) MRI growth patterns of plexiform neurofibromas in patients with neurofibromatosis type 1. J Med Genet 46(2):81–85 17. Van Den Boomgaard R, Van Balen R (1992) Methods for fast morphological image transforms using bitmapped binary images. CVGIP Gr Models Image Process 54(3):252–258 18. Weizman L, Hoch L, Bashat DB, Joskowicz L, Constantini S, Sira LB (2012) Interactive segmentation of plexiform neurofibroma tissue: method and preliminary performance evaluation. Med Biol Eng Comput 50(8):877–884 19. Williams VC, Lucas J, Babcock MA, Gutmann DH, Korf B, Maria BL (2009) Neurofibromatosis type 1 revisited. Pediatrics 123(1):124–133
123
