WO2010151229A1 - Procédé et système de segmentation d'une image du cerveau - Google Patents
Procédé et système de segmentation d'une image du cerveau Download PDFInfo
- Publication number
- WO2010151229A1 WO2010151229A1 PCT/SG2010/000228 SG2010000228W WO2010151229A1 WO 2010151229 A1 WO2010151229 A1 WO 2010151229A1 SG 2010000228 W SG2010000228 W SG 2010000228W WO 2010151229 A1 WO2010151229 A1 WO 2010151229A1
- Authority
- WO
- WIPO (PCT)
- Prior art keywords
- region
- points
- estimated
- regions
- image
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Ceased
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/143—Segmentation; Edge detection involving probabilistic approaches, e.g. Markov random field [MRF] modelling
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
- G06V10/00—Arrangements for image or video recognition or understanding
- G06V10/20—Image preprocessing
- G06V10/28—Quantising the image, e.g. histogram thresholding for discrimination between background and foreground patterns
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30016—Brain
Definitions
- the present invention relates to a method and system for segmenting a brain image.
- the segmented brain image may in turn be used for identifying stroke in a patient.
- Accurate gray and white matter segmentation of CT (Computer Tomography) images is important in the estimation of cerebrospinal fluid, white and gray matter regions, which in turn is useful in the prediction of clinical outcomes.
- CT Computer Tomography
- CSF cerebrospinal fluid
- WM white matter
- the infract intensity range is close to the intensity ranges of the hypo-dense WM (i.e. lower values of intensity amongst WM intensities) and the hyper-dense CSF (i.e. higher values of intensity amongst CSF intensities).
- Hu et al. 2008 proposed an approach of segmenting CT brain images using K-means and EM (Expectation-Maximization) clustering.
- the approach in Lee et al. 2008 is limited to the segmentation of the CSF (Cerebrospinal fluid), skull, parenchyma and calcification regions, and does not consider the segmentation of gray and white matter.
- Hu et al 2007 proposed an automatic segmentation method based on pattern recognition, fuzzy theory, anatomy, fractal dimensions etc.
- the method in Hu et al 2007 uses arbitrary parameters and has no quantitative measure of performance.
- a shape or texture based analysis can perform well only when the image quality is good. Qian et al.
- segmentation algorithms are based on the adjustment of several parameters [Heydarian et al 2009] and this requires experience and expertise. Moreover, such segmentation algorithms have their limitations when used on an organ composed of different tissue types (e.g. the brain). In addition, although genetic algorithms have been proposed for segmentation purposes, there are often issues with the convergence criteria and the huge computation time of these algorithms [Maulik 2009]. Also, most commercially available algorithms are not developed with CT images in mind due to the poor image quality of CT images. For example, some algorithms (such as the SPM5 which works for MR images only) require registration to an initial template if they are to be used for CT images. Summary of the invention
- the present invention aims to provide a new and useful method and system for identifying a CSF region in a brain image. This may be part of a process of segmenting the brain image into a CSF region, a WM region and a GM region.
- the present invention proposes, as part of a method for segmenting a brain image composed of intensity values at multiple points, estimating an upper limit for the intensity values of a CSF region in the image, by selecting a value for this upper limit such that the points of the image having an intensity value less than this upper limit and greater than zero include a subset of the points which form a spatially connected group of maximum size and which have a peaked intensity distribution.
- the invention exploits both the expected spatial distribution and expected intensity distribution of the CSF region. This makes it possible for the method to provide reliable discrimination of the CSF region even in CT images with poor image quality.
- the upper limit value is preferably varied such that the intensity distribution of the corresponding maximum size connected group of points gives rise to an optimum value of a predetermined significance criterion.
- the upper limit value is varied such as to maximize a predefined significance parameter.
- the upper limit estimated in this way may be further refined in further process steps, for example before it is used for segmentation of the CSF region of the image.
- Points in the image having intensity values above the estimated upper limit may be used in a process to obtain WM and GM regions.
- the method may be performed automatically (which is here used to mean that, although human interaction may initiate the algorithm, human interaction is not required while the algorithm is carried out).
- the algorithm may alternatively be performed semi-automatically (in which case there is human interaction with the computer during the processing).
- the algorithm is fully automatic, and there are no parameters to be set or adjusted by the user. No skull stripping is required, and there is no registration or machine learning.
- the invention may alternatively be expressed as a computer system for performing such a method.
- This computer system may be integrated with a device for capturing brain images.
- the invention may also be expressed as a computer program product, such as one recorded on a tangible computer medium, containing program instructions operable by a computer system to perform the steps of the method.
- Fig. 1 illustrates a flow diagram of a method for segmenting a brain image in an embodiment of the present invention
- Fig. 2 illustrates a flow diagram of a step for estimating initial parameters in the method of Fig. 1 ;
- Fig. 3 illustrates example histograms generated for the largest spatially connected component in different intensity ranges
- Fig. 4 illustrates a flow diagram of a step for optimizing initial parameters in the method of Fig. 1 ;
- Fig. 5 illustrates example Gaussian distributions corresponding to white and gray matter regions;
- Fig. 6 illustrates a flow diagram of a step for setting threshold values in the method of Fig. 1 ;
- Figs. 7(a) - (c) respectively illustrate initial estimated WM, GM and CSF regions in selected slices of an initial 3D segmented image obtained in the method of Fig. 1 and Figs. 7(d) - 7(f) respectively illustrate the estimated regions in Figs. 7(a) - 7(c) after a smoothing function is applied to these regions;
- Fig. 8 illustrates an example of a mask obtained for removing errors due to skull and dura mater in the method of Fig. 1 ;
- Fig. 9(a) illustrates the overlap regions between the initial estimated WM and GM regions in Figs. 7(d) and 7(e) and Fig. 9(b) illustrates the overlap regions between the initial estimated WM and CSF regions in Figs. 7(d) and 7(f);
- Figs. 10(a) and 10(b) respectively illustrate examples of unclassified points and small error regions in the image after obtaining the initial estimated WM, GM and CSF regions;
- Figs. 11 (a) - (c) illustrate a correction of a set of error points in the initial estimated WM 1 GM and CSF regions
- Fig. 12 illustrates the results produced by the method of Fig. 1 and the selection of different types of Ground Truth regions
- Fig. 13 illustrates a ROC curve obtained using ground truth GM and WM regions; and Fig. 14 illustrates a curve plotting median values of points in ground truth
- a method 100 which is an embodiment of the present invention, and which segments a brain image into a CSF region, a WM region and a GM region.
- the brain image comprises intensity values at respective points and is a 3D image since method 100 exploits the fact that in three-dimension and within a parenchyma region in a brain, the ventricle region is the largest connected component.
- the intensity values may be in Hounsfield units, or may be converted into Hounsfield units as a first step of the method.
- the points may be voxels of the 3D image.
- the brain image which is the input to method 100 may be a CT brain image. However, the brain image may be of a different image modality (e.g. an MRI image).
- a set of initial parameters are estimated and in step 104, some of the estimated initial parameters are optimized.
- threshold values are set based on the estimated initial parameters and in step 108, the image is segmented into estimated CSF, WM and GM regions using the threshold values.
- steps 112 and 114 overlap regions between the estimated regions, unclassified points and small error regions are resolved.
- step 116 further error regions are identified and resolved using known spatial characteristics of regions in the image.
- a final segmented image is obtained.
- Step 102 Estimate initial parameters
- step 102 an upper limit for the intensity values in the CSF region is first estimated.
- Initial parameters are then estimated using the estimated upper limit (i.e. CSFU).
- the initial parameters are the means and standard deviations of the intensity values in the CSF, WM and GM regions.
- Fig. 2 illustrates the sub-steps of step 102.
- an intensity range [0, X] comprising intensity values in the CSF region is first defined.
- X is an upper end of the intensity range and is variable.
- X is initially set as 10 so that for most of the acquired CT images, there is a high possibility that the intensity range lies within the range of intensity values in the CSF region of the image [Hacker and Artmann, 1976].
- points with intensity values within the intensity range are extracted from the image.
- a smoothing function is then applied on these extracted points. This smoothing is advantageous as it removes spatial error regions in the parenchyma. It is preferable to use an optimum smoothing function so that the shapes of the boundary regions can be preserved while smoothing the core regions in the image. For example, a median filter of window size [5,5] may be used on an image with a resolution of 1mm x 1mm to achieve optimum results.
- a histogram corresponding to the maximum size spatially connected three dimensional region in the extracted points from sub-step 204 is generated.
- This subset comprises those extracted points which constitute the largest spatially connected region in the image.
- This largest spatially connected region represents a three-dimensional core of ventricular regions in the image. Note that conceivably there may be cases in which all the extracted points are connected to each other, and in this case all those points would be in the "subset"; i.e. the term "subset" is used in this document to embrace the set of all the extracted points.
- a histogram of the intensity values of the subset of the extracted points is generated.
- a smoothing function is then applied to the generated histogram.
- the histogram is a plot where the x-axis indicates the intensity values of the points of the connected region, and the y- axis shows the frequencies with which these intensity values occur within the connected region.
- the method includes at least one iteration of (i) expanding the intensity range, and (ii) repeating sub-steps 202 - 206.
- the iterations continue until a difference between a currently generated histogram and a previously generated histogram indicates that the intensity range has been raised so high that it has begun to include a significant amount of the WM region. This is elaborated as follows.
- a significance parameter is calculated for the generated histogram as (Hmax - Hendy ⁇ Hmax - Hend) whereby Hmax is a maximum frequency (or voxel count) in the generated histogram and Hend is a frequency with which the maximum intensity value (i.e. the last bin) in the generated histogram occurs [Mandel 1984, Bevington 1969].
- Hmax is a maximum frequency (or voxel count) in the generated histogram
- Hend is a frequency with which the maximum intensity value (i.e. the last bin) in the generated histogram occurs [Mandel 1984, Bevington 1969].
- the difference between Hend and Hmax in the generated histogram is smaller when the upper end of the intensity range is very large or very small.
- the significance parameter is maximum when the upper end of the intensity range is at a critical intermediate value.
- Fig. 3 illustrates example histograms generated for the largest spatially connected component in different intensity ranges.
- the curves 302, 304, 306, 308 and 310 are respectively the example histograms generated for the intensity ranges [-5 to 10], [-5 to 15], [-5 to 20], [-5 to 25] and [-5 to 30].
- method 100 relies on the "flattening" of the post-peak slope of the histogram due to the contribution of WM intensity values in the post-peak slope region (as shown in the following sub-steps). Note that the histogram generated for the initial intensity range (defined in sub-step 202) need not be a monotonous function.
- sub-step 210 a check is performed to determine if the significance parameter for the currently generated histogram is less than the significance parameter for the previously generated histogram.
- sub-step 212 is performed in which the intensity range is expanded by increasing the upper end of the intensity range and sub-steps 204 - 210 are repeated.
- sub-step 210 determines that the difference between the currently generated histogram and the previously generated histogram indicates that the significance parameter has fallen, i.e. the connected component is starting to include significant numbers of points which are actually in the WM region. In this case, sub-step
- This estimated upper limit (i.e. CSFU) is used in the following steps for segmenting the image.
- a mean and a standard deviation of the intensity values in each of the CSF, WM and GM regions are estimated using CSFU as follows.
- CSFmean is estimated as an intensity value with a maximum frequency (i.e. Hmax) in the histogram generated when estimating the upper limit for the intensity values in the CSF region.
- CSFstd is estimated as a function of the intensity value at the left-side half width half maximum (HWHM) of the same histogram.
- HWHM left-side half width half maximum
- the "left-side" HWHM means the HWHM in the histogram which is lower than Hmax.
- a first range of intensity values is defined wherein the lower end of this first range is the estimated upper limit, CSFU.
- this first range is [CSFU, CSFU + L] whereby L is set as 50. From experimentation, it was found that with L set as 50, the first range comprises intensity values of most forms of parenchyma regions in different types of data (i.e. different image modalities). A first set of points with intensity values within this first range are then identified within the image and a histogram is formed using this first set of points. The histogram is a plot where the x-axis indicates the intensity values of the first set of points, and the y-axis shows the frequencies with which these intensity values occur.
- WMmean, GMmean, WMstd and GMstd are then derived using the histogram.
- an intensity boundary to partition points in the image into an estimated WM region and an estimated GM region (or in other words, an estimated upper limit for the intensity values in the WM region, WMU) is obtained.
- the intensity boundary is defined as an intensity value with a maximum frequency in the formed histogram.
- Second and third ranges of values are defined to meet at the value WMU.
- the widths of the second and third ranges are equal, and the second and third ranges are [CSFU, WMU] and [WMU, WMU+(WMU - CSFU)] respectively.
- a second and third set of points with intensity values lying within the second and third ranges respectively are identified in the image.
- the second set of points is taken as an estimated WM region
- the third set of points is taken as an estimated GM region.
- the WMmean and WMstd are then estimated as a mean and a standard deviation of the intensity values of the second set of points
- GMmean and GMstd are estimated as a mean and a standard deviation of the intensity values of the third set of points.
- Step 104 Optimize initial parameters
- step 104 a least squares fit algorithm is used to optimize (in other words, improve) the estimates of WMmean, WMstd, GMmean and GMstd.
- Fig. 4 illustrates the sub-steps of step 104.
- sub-step 402 a distribution function is calculated for each of the estimated WM and GM regions based on WMmean, WMstd and GMmean, GMstd. This is done by setting the mean and standard deviation of the distribution function for the WM region as WMmean and WMstd respectively, and the mean and standard deviation of the distribution function for the GM region as GMmean and GMstd respectively.
- the distribution function may be any distribution as long as it is a function of a mean and a standard deviation.
- the distribution function may be a Gaussian distribution or a lognormal distribution.
- the distribution function is a Gaussian distribution and is estimated by normalizing at the value of the respective mean (WMmean or GMmean). The advantage of using a Gaussian distribution is that this is one of the simplest functions available.
- the statistical parameters (i.e. the mean and standard deviation) of the distribution functions for the WM and GM regions are varied locally to fit the distributions of the intensity values in the estimated WM and GM regions respectively as follows.
- the statistical parameters are varied within a tight range of ⁇ 5% of their initial values as set above. For example, if WMmean is 35, the mean of the distribution function for the WM region is varied between 33.25 to 36.75. Similarly, if WMstd or GMstd is 3, the standard deviations of the distribution functions for the WM and GM regions are varied between 2.85 to 3.15.
- the final statistical parameters of the distribution function for the WM region are selected as the parameters which achieve a minimum sum of squared differences between the distribution function and the distribution of intensity values in the estimated WM region.
- the final statistical parameters are obtained using the least squares fitting algorithm in Matlab. To avoid spurious fits, the parameters are varied only locally in the tight range as mentioned above.
- WMmean and WMstd are then set as the final mean and standard deviation of the distribution function for the WM region. These are the improved estimates of WMmean and WMstd.
- the final statistical parameters of the distribution function for the GM region are selected as the parameters which achieve a minimum sum of squared differences between the distribution function and the distribution of the intensity values in the estimated GM region.
- the final statistical parameters are obtained using the least squares fitting algorithm in Matlab.
- GMmean and GMstd are then set as the final mean and standard deviation of the distribution function for the GM region.
- Fig. 5 illustrates example Gaussian distributions corresponding to the white and gray matter regions.
- curves 502, 504 and 506 respectively plot the parenchyma histogram, an example Gaussian distribution corresponding to the white matter region (WM gaussian fit) and an example Gaussian distribution corresponding to the gray matter region (GM gaussian fit).
- Step 106 Set threshold values based on the initial parameters by minimizing overlap regions
- a lower threshold value and an upper threshold value are set for each of the CSF, WM and GM regions based on the initial parameters, CSFmean, CSFstd, WMmean, WMstd, GMmean, GMstd.
- the lower threshold value for the WM region is set as WMmean - 1.96 * WMstd
- the upper threshold value for the GM region is set as GMmean + 1.96 * GMstd.
- Fig. 6 illustrates the sub-steps for setting the upper threshold value for the WM region and the lower threshold value for the GM region. These sub-steps aim to achieve a minimum overlap area between the estimated WM and GM regions.
- the upper threshold value for the WM region and the lower threshold value for the GM region are first set as WMmean + ss*WMstd and GMmean - ss*GMstd respectively whereby ss is a variable.
- a first plurality of points with intensity values lying between the lower and upper threshold values for the WM region and a second plurality of points with intensity values lying between the lower and upper threshold values for the GM region are extracted from the image.
- the first and second plurality of points form trial WM and GM regions respectively.
- the initial value of ss may be set as a very small value, for example 0.5 and in this case, the overlap area is likely to be zero.
- a smoothing function is then applied to the first and second plurality of points.
- an overlap area between the first and second plurality of points is determined whereby the overlap area is the number of points belonging to both the first and second plurality of points.
- a check is then performed to determine if the overlap area is larger than an error value representative of an error in the total number of points in the first and second plurality of points.
- the error value is calculated as a function of ⁇ j ( Area W M +AreaGM) whereby AreaWM is the number of points in the first plurality of points and AreaGM is the number of points in the second plurality of points. The experimental results discussed below are obtained by setting the error value equal to ⁇ (AreaWM +AreaGM) .
- the upper threshold value for the WM region and the lower threshold value for the GM region are set as their current values. If not, the upper threshold value for the WM region and the lower threshold value for the GM region are adjusted and sub-steps 604 - 608 are repeated until the check in sub-step 608 returns a positive result.
- ss is selected to be the value which maximizes the sizes of the trial WM and GM regions with the constraint that the overlap between them is bounded by the error value.
- the adjustment of the thresholds is performed by increasing the value of ss. In one example, the value of ss is increased in steps of 0.1.
- a similar process can be used to obtain the cut-off between the CSF region and the WM region.
- CSFmean and CSFstd estimates have been obtained, and these are used to derive a Gaussian function.
- the lower threshold for the WM region is then redefined as WMmean-sscsf * WMstd.
- the following sub-steps 604a - 608a similar to sub-steps 604 - 608 are then performed.
- a third plurality of points with intensity values lying between the lower and upper threshold values for the CSF region and a fourth plurality of points with intensity values lying between the lower and upper threshold values for the WM region are extracted from the image.
- the third and fourth plurality of points form trial CSF and WM regions respectively.
- the initial value of sscsf may be set as a very small value, for example 0.5 and in this case, the overlap area is likely to be zero.
- a smoothing function is then applied to the third and fourth plurality of points.
- an overlap area between the third and fourth plurality of points is determined whereby the overlap area is the number of points belonging to both the third and fourth plurality of points.
- a check is then performed to determine if the overlap area is larger than an error value representative of an error in the total number of points in the third and fourth plurality of points.
- the error value is calculated as a function of whereby AreaCSF is the number of points in the third plurality of points and AreaWM is the number of points in the fourth plurality of points. The experimental results discussed below are obtained by setting the error value equal to ⁇ (AreaCSF +AreaWM) .
- the upper threshold value for the CSF region and the lower threshold value for the WM region are set as their current values. If not, the upper threshold value for the CSF region and the lower threshold value for the WM region are adjusted and sub-steps 604a - 608a are repeated until the check in sub-step 608a returns a positive result.
- sscsf is selected to be the value which maximizes the sizes of the trial CSF and WM regions with the constraint that the overlap between them is bounded by the error value.
- the adjustment of the thresholds is performed by increasing the value of sscsf. In one example, the value of sscsf is increased in steps of 0.1.
- Step 108 Segment image using the threshold values
- the image is segmented using the threshold values obtained in step 106.
- an estimated CSF region, WM region and GM region comprising intensity values lying between the upper and lower threshold values for the CSF region, WM region and GM region respectively are extracted from the image.
- Figs. 7(a) - 7(c) respectively illustrate the initial estimated WM, GM and CSF regions in selected slices of the 3D segmented image obtained in step 108 whereas Figs. 7(d) - 7(f) respectively illustrate the estimated regions in Figs. 7(a) - 7(c) after a smoothing function is applied to these regions.
- the slices illustrated in Figs. 7(a) - 7(c) are axial slices selected at random from slices identified using domain knowledge regarding the distribution of CSF, WM, and GM in different planes.
- Step 110 Remove errors due to skull and dura mater
- a fourth range of intensity values is defined based on the lower threshold value for the CSF region and the upper threshold value for the GM region.
- the fourth range of intensity values is [CSFIower, GMupper] wherein CSFIower is the lower threshold value for the CSF region and GMupper is the upper threshold value for the GM region.
- a fourth set of points with intensity values lying within the fourth range is extracted from the image.
- a mask comprising non-zero values in an area within an inner skull boundary is created.
- the inner skull boundary is a boundary of a largest connected component formed by the fourth set of points in the image and the mask is created by filling the volume within the inner skull boundary with non-zero values.
- Fig. 8 illustrates an example of a mask obtained in step 110.
- Step 112 Resolve overlap regions between the estimated CSF and WM regions and between the estimated WM and GM regions
- step 112 overlap regions between the estimated CSF and WM regions, and overlap regions between the estimated WM and GM regions are resolved.
- Fig. 9(a) illustrates the overlap regions between the estimated WM and GM regions in Figs. 7(d) and 7(e) respectively whereas
- Fig. 9(b) illustrates the overlap regions between the estimated WM and CSF regions in Figs. 7(d) and 7(f) respectively.
- a first set of overlap points belonging to both the estimated CSF and WM regions is extracted and a mean of the intensity values of the first set of overlap points is then calculated. For each overlap point in the first set, a check is then performed to determine if the intensity value of the overlap point is greater than the mean of the intensity values of the first set of overlap points. If so, the estimated WM region is updated to comprise the overlap point. If not, the estimated CSF region is updated to comprise the overlap point. This is because points with higher intensity values are more likely points in the WM region whereas points with lower intensity values are more likely points in the CSF region.
- a second set of overlap points belonging to both the estimated WM and GM regions is extracted and a mean of the intensity values of the second set of overlap points is then calculated. For each overlap point in the second set, a check is then performed to determine if the intensity value of the overlap point is greater than the mean of the intensity values of the second set of overlap points. If so, the estimated GM region is updated to comprise the overlap point. If not, the estimated WM region is updated to comprise the overlap point. This is because points with higher intensity values are more likely points in the GM region whereas points with lower intensity values are more likely points in the WM region.
- a median filter is applied to the image prior to extracting the estimated CSF, WM and GM regions.
- a median filter with a narrow window for example, [3,3]
- Step 114 Resolve unclassified points and small error regions
- step 114 a set of unclassified points which do not belong to any of the estimated regions is detected from the image and resolved via voting. Small error regions are also resolved via voting in step 114.
- Figs. 10(a) and 10(b) illustrate examples of unclassified points 1002 and small error regions 1004 respectively.
- the set of unclassified points is resolved as follows. For each unclassified point, it is determined which estimated region comprises a maximum number of votes, in other words, comprises the most number of points in a neighbourhood region around the unclassified point. In one example, the neighbourhood region is a square region. The estimated region comprising the maximum number of votes (i.e. a majority region) is then updated to comprise the unclassified point.
- a size of the neighbourhood region around the unclassified point is increased and this is repeatedly performed until a majority region for the unclassified point can be determined.
- the size of the largest neighbourhood region in step 114 depends on the size of the largest region formed by the unclassified points in the image.
- the detected set of unclassified points form regions of limited size in the image.
- the number of points in the regions of limited size is less than or equal to 0.01% of the number of points in a parenchyma region in the image (in one example, the number of points in the parenchyma region is estimated as the total number of points in the estimated CSF, WM and GM regions). This is because the unclassified points forming the larger regions in the image are more likely points in the vessels, calcification segments and sulci etc.
- Points in small error regions are resolved in the same way as the unclassified points.
- small error regions are defined as regions having a number of points less than or equal to 0.01% of the number of points in a parenchyma region in the image.
- Step 116 Identify and resolve further error regions using known spatial characteristics of regions in the image
- step 116 known spatial distributions of points in the WM, GM and CSF regions in orthogonal planes (axial, coronal and/or sagittal plane) of typical brain images are used to identify the most (or the least probable tissue type) in a particular spatial location. This information is then used to reallocate error points in the estimated CSF, WM and GM regions.
- step 116 To the brain skull boundary of the parenchyma region in a typical brain image, the probability that a point belongs to the CSF region or the GM region is much higher than the probability that the point belongs to the WM region e.g. see reference for distribution functions [14]. This observation is applied in step 116 to remove errors created due to the incorrect identification of the points in the "hyper-dense CSF regions" and the "hypo-dense GM regions" as points in the WM region of the cortical region.
- a first set of error points is extracted.
- the first set of error points are points which form small regions in the cortical region of the image and which belong to the estimated WM region.
- the small regions are defined as connected regions having a number of points less than a fraction (for example, 0.05%) of the number of points in a parenchyma region in the image.
- the number of points in the parenchyma region is estimated as the total number of points in the estimated CSF, WM and GM regions. These small regions are generally attached to the contour of the parenchyma region.
- the first set of error points may be located as follows. A perimeter of the brain is first determined in a slice of the image. The perimeter of the brain is then thickened by approximately 2mm on each side. Next, small regions in the estimated WM region connected to the perimeter of the brain are located. The points in these small regions form the first set of error points.
- a check is performed to determine if the intensity value of the error point is greater than WMmean and if so, the estimated GM region is updated to comprise the error point. If not, the estimated CSF region is updated to comprise the error point.
- Figs. 11(a) - (c) illustrate the correction of the first set of error points.
- Fig. 11(a) illustrates a selected slice of the input image
- Fig. 11(b) illustrates the estimated GM, WM and CSF regions with the first set of error points
- Fig. 11(c) illustrates the estimated GM, WM and CSF regions with the first set of error points corrected.
- step 116 Points in regions above the ventricular system in the image are extracted.
- the ventricular system is taken as a maximum connected region in the estimated CSF region.
- the image is filtered prior to extracting these points.
- the second set of error points comprise the points in the regions above the ventricular system which belong to both the estimated WM region and regions connected to or in a vicinity of the MSP in the image.
- the regions connected to or in a vicinity of the MSP in the image are located as follows. The MSP is first located on a slice in the image (the MSP is a line on the slice). The MSP is then thickened by approximately 2mm on each side. Regions connected to the MSP are then located and this helps to localize the small WM regions not connected to the MSP in the vicinity of the MSP.
- the number of points in each region connected to or in a vicinity of the MSP may be less than 0.05% of the number of points in the parenchyma region in the image.
- a point in a vicinity of the MSP is defined as a point which is at most 3 points away from the MSP.
- a check is performed to determine if the intensity value of the error point is greater than WMmean and if so, the estimated GM region is updated to comprise the error point. If not, the estimated CSF region is updated to comprise the error point.
- step 116 In typical brain scans, it is observed that points in the WM region are connected in three-dimension. In other words, a small isolated group of points in three- dimension identified as points belonging to the WM region is most likely an error. This observation is applied in step 116 as follows.
- a third set of error points belonging to the estimated WM region and forming small connected 3D regions in the image is extracted.
- a small connected 3D region is defined as a connected 3D region with a number of points less than 0.05% of the number of points in a parenchyma region in the image.
- the number of points in the parenchyma region is estimated as the total number of points in the estimated CSF, WM and GM regions.
- a check is performed to determine if the intensity value of the error point in the third set is greater than WMmean and if so, the estimated GM region is updated to comprise the error point. If not, the estimated CSF region is updated to comprise the error point.
- Step 118 Obtain a final segmented image
- step 118 a final segmented image is obtained by combining the estimated CSF, WM and GM regions.
- Fig. 12 illustrates the results obtained using method 100.
- A1 , B1 , C1 , D1 , E1 , F1 illustrate slices of the input brain images (these are the original unenhanced CT images) whereas A2, B2, C2, D2, E2, F2 illustrate the corresponding slices of the final segmented images obtained using method 100 overlaid with corresponding slices of ground truth images in which the locations of the CSF, WM and GM regions are marked by clinical experts.
- the slices illustrated in Fig. 12 are axial slices obtained at different levels in the axial plane of the brain image. These levels range from the inferior portion of the brain to the superior portion of the brain.
- Areas 1202, 1204, 1206, 1210, 1212 and 1214 are examples of correctly segmented areas whereas area 1208 is an example of an incorrectly segmented area.
- the results in Fig. 12 show that the number of correctly segmented areas is significantly higher than the number of incorrectly segmented areas.
- Fig. 13 illustrates a ROC curve obtained by evaluating the ground truth GM and WM regions against estimated GM and WM regions obtained using the intensity boundary (WMU). These estimated GM and WM regions comprise intensity values within the ranges [WMU, WMU + (WMU - CSFU)] and [CSFU, WMU] respectively. As shown in Fig. 13, choosing WMU to be near an intensity value corresponding to a peak of the parenchyma histogram (i.e.
- the histogram of the intensity values in the parenchyma region of the image results in a sensitivity of 85.6 and a specificity of 86.1.
- the WMU chosen in this manner provides a good initial estimate of the boundary between the WM region and the GM region. Further fine tuning of this initial estimate is then performed in the subsequent steps as described above.
- Fig. 14 illustrates a curve plotting the median values of the points in the ground truth WM and GM regions against the median values of the points in the estimated WM and GM regions from method 100.
- the square markers denote the points in the WM regions and the circle markers denote the points in the GM regions.
- the median values of the points in the estimated regions from method 100 are approximately the same as the median values of the points in the ground truth regions.
- Method 100 uses an adaptive approach which combines intensity, spatial and statistical properties of the image to optimize parameters automatically whereby these parameters may be used in the setting of threshold values. Thus, method 100 does not require a user to set or adjust parameters. Neither does method 100 require the resetting of a large number of parameters (unlike, in level set algorithms). Furthermore, method 100 does not employ skull stripping steps, registration steps or a machine learning algorithm. Thus, method 100 is less prone to errors arising due to the need for multiple parameter adjustments, initial guesses or convergence criteria, registration, machine learning or skull stripping. This allows method 100 to be more robust to intensity variability and poor image contrast.
- Method 100 calculates threshold values using the shapes of a plurality of histograms generated using intensity values in the image. Thus, it is robust to instrumentation parameters that may affect the image quality of the image.
- spatial smoothing is performed in method 100, hence reducing the amount of overlap between the intensity distributions of different tissue regions.
- spatial smoothing is performed in step 102, hence minimizing the contribution of intensity values in the WM region when calculating the estimated upper limit for the intensity values in the CSF region.
- spatial filtering is performed and this helps to remove error points so that the threshold values can be estimated from the less contaminated regions in the image.
- method 100 only utilizes the histograms corresponding to certain tissue regions to estimate the upper limit for the intensity values in the CSF region (only the histograms of points in a largest connected component are generated in one example of step 102).
- the shapes of the histograms are mainly influenced by the intensity values in the CSF region and not the WM region.
- steps to reduce false positive and false negative errors are performed. This allows the method 100 to produce more accurate results.
- steps are performed to examine the boundaries of the estimated regions which may be predominant areas in which errors exist. For example, the higher intensity points in the boundaries are classified into the region known to have points of higher intensity values whereas the lower intensity points in the boundaries are classified into the region known to have points with lower intensity values.
- the unclassified points and small error regions i.e. small contaminated regions
- the contaminations in the WM region are also resolved using known spatial characteristics of regions in the image.
- points in the estimated WM region form an isolated small region, they are considered points of a contamination region because it is expected that WM regions exist as large 3D connected regions.
- method 100 is also versatile. It can be used to process input images of different image modalities.
- method 100 is tested using a median filter and an anisotropic diffusion filter [Perna and Malik] for performing its smoothing operations and it is shown that method 100 is robust to the type of filter used for the smoothing operations.
- Torbey MT Selim M
- Knorr J Bigelow C
- Mah L Quantitative Analysis of the Loss of Distinction Between Gray and White Matter in Comatose Patients After Cardiac Arrest. Stroke 31 :2163-2167.
- Nowinski WL Gupta V
- Wai Yen Chan Wai Yen Chan
- YY Sitoh Kang Sim. Use of a local gray to white matter relationship to study the brain in health and disease. Submitted to Academic Radiology.
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Multimedia (AREA)
- Probability & Statistics with Applications (AREA)
- Software Systems (AREA)
- Apparatus For Radiation Diagnosis (AREA)
- Image Analysis (AREA)
Abstract
L'invention concerne un procédé de segmentation d'une image du cerveau en une région CSF (fluide cérébro-spinal), une région WM (matière blanche) et une région GM (matière grise). Une limite supérieure des valeurs d'intensité d'une région CSF dans l'image est estimée afin que les points de l'image qui ont une intensité inférieure à cette limite supérieure comprennent un sous-ensemble des points qui forment un groupe connecté spatialement et qui ont une distribution avec pic d'intensité. En d'autres termes, l'invention exploite à la fois la distribution spatiale attendue et la distribution d'intensité attendue de la région CSF. Cela permet au procédé d'assurer une discrimination fiable de la région CSF, même sur des images CT avec une qualité de l'image médiocre. Différents procédés sont proposés pour utiliser la limite supérieure et pour améliorer la précision de la segmentation.
Priority Applications (3)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| EP10792425.0A EP2446418A4 (fr) | 2009-06-23 | 2010-06-17 | Procédé et système de segmentation d'une image du cerveau |
| SG2011092822A SG176860A1 (en) | 2009-06-23 | 2010-06-17 | A method and system for segmenting a brain image |
| US13/379,649 US8831328B2 (en) | 2009-06-23 | 2010-06-17 | Method and system for segmenting a brain image |
Applications Claiming Priority (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| SG200904356 | 2009-06-23 | ||
| SG200904356-3 | 2009-06-23 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| WO2010151229A1 true WO2010151229A1 (fr) | 2010-12-29 |
Family
ID=43386783
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| PCT/SG2010/000228 Ceased WO2010151229A1 (fr) | 2009-06-23 | 2010-06-17 | Procédé et système de segmentation d'une image du cerveau |
Country Status (4)
| Country | Link |
|---|---|
| US (1) | US8831328B2 (fr) |
| EP (1) | EP2446418A4 (fr) |
| SG (1) | SG176860A1 (fr) |
| WO (1) | WO2010151229A1 (fr) |
Families Citing this family (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2012056362A1 (fr) * | 2010-10-25 | 2012-05-03 | Koninklijke Philips Electronics N.V. | Système de segmentation d'une image médicale |
| WO2015013719A1 (fr) | 2013-07-26 | 2015-01-29 | Li-Cor, Inc. | Filtre de bruit adaptatif |
| US10395350B2 (en) | 2013-07-26 | 2019-08-27 | Li-Cor, Inc. | Adaptive background detection and signal quantification systems and methods |
| WO2015013722A1 (fr) | 2013-07-26 | 2015-01-29 | Li-Cor, Inc. | Systèmes et procédés pour régler des paramètres initiaux d'affichage |
| US9530206B2 (en) * | 2015-03-31 | 2016-12-27 | Sony Corporation | Automatic 3D segmentation and cortical surfaces reconstruction from T1 MRI |
| EP3289563B1 (fr) | 2015-04-30 | 2020-08-05 | Koninklijke Philips N.V. | Classification des tissus cérébraux |
| RU2713707C2 (ru) * | 2015-04-30 | 2020-02-06 | Конинклейке Филипс Н.В. | Классификация тканей головного мозга |
| KR101995900B1 (ko) | 2017-09-11 | 2019-07-04 | 뉴로핏 주식회사 | 3차원 뇌지도 생성 방법 및 프로그램 |
| CN109872328B (zh) * | 2019-01-25 | 2021-05-07 | 腾讯科技(深圳)有限公司 | 一种脑部图像分割方法、装置和存储介质 |
| CN111160260B (zh) * | 2019-12-30 | 2023-07-04 | 湖北航天技术研究院总体设计所 | 一种sar图像目标检测方法及系统 |
| EP3996036A1 (fr) | 2020-11-09 | 2022-05-11 | Koninklijke Philips N.V. | Angiographie par résonance magnétique à l'aide d'images pondérées en t2 |
| CN116129107B (zh) * | 2022-11-17 | 2025-10-03 | 华中科技大学 | 基于长短期记忆自注意力模型的三维医学图像分割方法及系统 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6430430B1 (en) * | 1999-04-29 | 2002-08-06 | University Of South Florida | Method and system for knowledge guided hyperintensity detection and volumetric measurement |
| WO2006114003A1 (fr) * | 2005-04-27 | 2006-11-02 | The Governors Of The University Of Alberta | Methode et systeme de detection et de segmentation automatiques de tumeurs et d'oedemes associes (tumefaction) dans des images a resonance magnetique (irm) |
| WO2008069762A1 (fr) * | 2006-12-06 | 2008-06-12 | Agency For Science, Technology And Research | Procédé pour identifier une zone pathologique d'un balayage, comme une zone d'accident ischémique cérébral d'un balayage irm |
| US7428323B2 (en) * | 2000-11-14 | 2008-09-23 | Yitzchak Hillman | Method and system for automatic diagnosis of possible brain disease |
Family Cites Families (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6658080B1 (en) * | 2002-08-05 | 2003-12-02 | Voxar Limited | Displaying image data using automatic presets |
| WO2004077359A1 (fr) | 2003-02-27 | 2004-09-10 | Agency For Science, Technology And Research | Procede et appareil d'extraction d'un systeme cerebro-ventriculaire a partir d'images |
| WO2006121410A1 (fr) * | 2005-05-11 | 2006-11-16 | Agency For Science, Technology And Research | Procede, dispositif et logiciel permettant de segmenter le cerveau a partir de donnees d'imagerie par resonance magnetique |
-
2010
- 2010-06-17 SG SG2011092822A patent/SG176860A1/en unknown
- 2010-06-17 US US13/379,649 patent/US8831328B2/en not_active Expired - Fee Related
- 2010-06-17 EP EP10792425.0A patent/EP2446418A4/fr not_active Withdrawn
- 2010-06-17 WO PCT/SG2010/000228 patent/WO2010151229A1/fr not_active Ceased
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6430430B1 (en) * | 1999-04-29 | 2002-08-06 | University Of South Florida | Method and system for knowledge guided hyperintensity detection and volumetric measurement |
| US7428323B2 (en) * | 2000-11-14 | 2008-09-23 | Yitzchak Hillman | Method and system for automatic diagnosis of possible brain disease |
| WO2006114003A1 (fr) * | 2005-04-27 | 2006-11-02 | The Governors Of The University Of Alberta | Methode et systeme de detection et de segmentation automatiques de tumeurs et d'oedemes associes (tumefaction) dans des images a resonance magnetique (irm) |
| WO2008069762A1 (fr) * | 2006-12-06 | 2008-06-12 | Agency For Science, Technology And Research | Procédé pour identifier une zone pathologique d'un balayage, comme une zone d'accident ischémique cérébral d'un balayage irm |
Non-Patent Citations (4)
| Title |
|---|
| CLINE H.E ET AL.: "Three dimensional Segmentation of MR images of the head using Probability and Connectivity", JOURNAL OF COMPUTER ASSISTED TOMOGRAPHY, vol. 14, no. 6, 1990, pages 1037 - 1045, XP000601360 * |
| DZUNG L. PHAM ET AL.: "Current Methods In Medical Image Segmentation", ANNU. REV. BIOMED ENG 2000., vol. 02, 2000, pages 315 - 37 7, XP009062827 * |
| See also references of EP2446418A4 * |
| WANG ET AL.: "Correction for variations in MRI scanner sensitivity in brain studies with histogram matching", MRM, vol. 39, 1998, pages 322 - 327, XP008148820 * |
Also Published As
| Publication number | Publication date |
|---|---|
| SG176860A1 (en) | 2012-01-30 |
| EP2446418A1 (fr) | 2012-05-02 |
| EP2446418A4 (fr) | 2013-11-13 |
| US20120099779A1 (en) | 2012-04-26 |
| US8831328B2 (en) | 2014-09-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US8831328B2 (en) | Method and system for segmenting a brain image | |
| Park et al. | Skull stripping based on region growing for magnetic resonance brain images | |
| CN105389811B (zh) | 一种基于多级阈值分割的多模态医学图像处理方法 | |
| EP2070045B1 (fr) | Diagnostic assisté par ordinateur avancé de nodules du poumon | |
| Balan et al. | Smart histogram analysis applied to the skull-stripping problem in T1-weighted MRI | |
| Chiverton et al. | Statistical morphological skull stripping of adult and infant MRI data | |
| Hu et al. | Segmentation of brain from computed tomography head images | |
| Capelle et al. | Unsupervised segmentation for automatic detection of brain tumors in MRI | |
| CN108885786B (zh) | 医学图像处理 | |
| Gupta et al. | Automatic segmentation of cerebrospinal fluid, white and gray matter in unenhanced computed tomography images | |
| Goceri | Automatic kidney segmentation using Gaussian mixture model on MRI sequences | |
| US20100049035A1 (en) | Brain image segmentation from ct data | |
| KR101474162B1 (ko) | 흉부 컴퓨터 단층 영상의 간유리음영 결절 자동 분할 시스템 및 그 분할 방법 | |
| Shelke et al. | Automated segmentation and detection of brain tumor from MRI | |
| Silvoster M et al. | Efficient segmentation of lumbar intervertebral disc from MR images | |
| Somasundaram et al. | Skull stripping of MRI head scans based on chan-vese active contour model | |
| Resmi et al. | A novel automatic method for extraction of glioma tumour, white matter and grey matter from brain magnetic resonance images | |
| Myint et al. | Effective kidney segmentation using gradient based approach in abdominal CT images | |
| Rao et al. | Implementation of clustering techniques for brain tumor detection | |
| Samundeeswari et al. | A novel multilevel hybrid segmentation and refinement method for automatic heterogeneous true NSCLC nodules extraction | |
| Banik et al. | Automatic detection, extraction and mapping of brain tumor from MRI scanned images using frequency emphasis homomorphic and cascaded hybrid filtering techniques | |
| Abirami et al. | Kidney segmentation for finding its abnormalities in abdominal CT images | |
| WO2006121410A1 (fr) | Procede, dispositif et logiciel permettant de segmenter le cerveau a partir de donnees d'imagerie par resonance magnetique | |
| Hajiesmaeili et al. | A new approach to locate the hippocampus nest in brain MR images | |
| Sumir et al. | Segmentation of brain tumor from MRI images using fast marching method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| 121 | Ep: the epo has been informed by wipo that ep was designated in this application |
Ref document number: 10792425 Country of ref document: EP Kind code of ref document: A1 |
|
| DPE2 | Request for preliminary examination filed before expiration of 19th month from priority date (pct application filed from 20040101) | ||
| WWE | Wipo information: entry into national phase |
Ref document number: 13379649 Country of ref document: US |
|
| NENP | Non-entry into the national phase |
Ref country code: DE |
|
| WWE | Wipo information: entry into national phase |
Ref document number: 2010792425 Country of ref document: EP |