Abstract:
A method is proposed for binarising an image by deriving an intensity threshold and classifying pixels according to whether their intensity is below or above the threshold. In the derivation of the threshold, prior konwledge is used to define a region of interest (ROI) in the image. Furthermore, prior knowledge is used to select a range in the frequency distribution of the intensities of the pixels in the ROI, and that only data within this frequency range is used to derive the threshold. These techniques provide a highly effective mechanism for incorporating prior knowledge into the threshold selection which is critical whether the image is a medical image or not. In particular, a threshold can be found to binarise images which exhibits high robustness to imaging artefacts such a gray level inhomogeneity and noise.

Description:
FIELD OF THE INVENTION  
       [0001]     The present invention relates to methods for processing an image so as to classify pixels of the image based on an intensity threshold. In particular, the invention relates to such a method having an improved process for selection of the threshold. The invention is applicable to both medical and non-medical images.  
       BACKGROUND OF INVENTION  
       [0002]     Binarisation is a well-known technique for image segmentation—that is classifying pixels of the image into two classes. Binarisation performs this classification based on whether a given pixel of the image has an intensity (gray-level) above or below a threshold. Binarisation has been widely applied to a number of image processing and computer vision applications, as a preliminary segmentation step. It makes an implicit assumption that an object of interest in the image has different intensity values from other (background) portions of the image.  
         [0003]     Many techniques exist for selection of the threshold. For example, in some such processes, the threshold can be selected in a process involving user interaction, while in other processes the threshold is selected entirely automatically. In some such processes the threshold is selected locally (i.e. such that the threshold varies from one pixel to another), while in other processes the threshold is the same over the whole image.  
         [0004]     Most automatic threshold selection methods employ a histogram of the gray levels in the image. For example, Otsu [1] proposed a selection of the threshold to maximise the separability of the resultant classes in gray levels, which is performed by minimising the within-class variance. Li and Lee [2] selected the threshold by minimising the cross entropy between the image and its segmented version. Kittler and lllingworth [3] selected the threshold by minimising the Bayes errors under the assumption that the object and pixel gray level values are normally distributed. Kapur et al [4] provided a maximum entropy approach. Wong and Sahoo [5] maximised the entropy with constraints on the region homogeneity and object boundary. Saha and Udupa [6] proposed a technique which maximised class uncertainty and homogeneity of the regions. Cheng et al [7] used the concept of fuzzy c-partition and the maximum fuzzy entropy principle to select a threshold.  
         [0005]     Cheung at al (U.S. Pat. No. 5,231,580A, 1993) disclosed an automatic method to characterise nerve fibres using local thresholds. It first partitions the entire image into sub-images and finds the threshold for each sub-image using a histogram-based thresholding method. Then, the pixel-wise threshold is approximated by interpolating the thresholds of neighbouring subimages.  
       SUMMARY OF THE INVENTION  
       [0006]     It is observed that the existing methods for selecting a threshold described above lack a mechanism for incorporating prior knowledge about the images to be binarised.  
         [0007]     Thus, the present invention aims to provide a new and useful technique for selecting a threshold for binarising an image, and in particular one which enables prior knowledge to be explicitly incorporated.  
         [0008]     In general terms, the invention proposes firstly that this prior knowledge is used to define a region of interest (ROI) in the image, such that the analysis of frequency distribution of pixel intensities (represented by a frequency histogram) is performed only for pixels in the ROI. Secondly, the invention proposes that the prior knowledge is used to select an intensity range, and that only pixels within this intensity range are used to generate the frequency distribution from which the threshold is selected.  
         [0009]     These two ideas are in principle separate, but in combination they provide a highly effective mechanism for incorporating prior knowledge into the threshold selection. The advantage is critical whether the image is a medical one or not. In particular, a threshold can be found to binarise images which exhibits high robustness to imaging artefacts such as gray level inhomogeneity and noise.  
         [0010]     Specifically, one expression of the invention is a method of binarising an image composed of pixels having respective intensity values, the method comprising: 
        (i) using prior knowledge about the image to derive a region of interest within it;     (ii) using prior knowledge about the image to derive an intensity range of pixels in the said region of interest;     (iii) obtaining a frequency distribution of the intensities within the said intensity range of pixels within the said region of interest;     (iv) using the said frequency distribution to derive an intensity threshold; and     (v) binarising the image by classifying pixels in the said region of interest according to whether their intensities are above or below the said intensity threshold.        
 
         [0016]     The invention may alternatively be expressed as a computer system which is set up to perform such a method. Alternatively, it can be expressed as software for performing the method. 
     
    
     BRIEF DESCRIPTION OF THE FIGURES  
       [0017]     Preferred features of the invention will now be described, for the sake of illustration only, with reference to the following figures in which:  
         [0018]      FIG. 1  shows the steps in a method which is an embodiment of the invention;  
         [0019]      FIG. 2  shows an MR SPGR intercommissural axial slice of a brain, which is a suitable subject for the method of  FIG. 1 ;  
         [0020]      FIG. 3  shows a region of interest within the image of  FIG. 2  derived by a first step of the method of  FIG. 1 ;  
         [0021]      FIG. 4  is a gray-level histogram of the ROI shown in  FIG. 3 , and a threshold selected in one form of a step of the method of  FIG. 1 ; and  
         [0022]      FIG. 5  shows the binarised image using the threshold selected in the method of  FIG. 1 . 
     
    
     DETAILED DESCRIPTION OF THE EMBODIMENTS  
       [0023]     Referring firstly to  FIG. 1 , the overall steps of a method which is an embodiment of the invention are shown.  
         [0024]     In step  1 , an image is input.  
         [0025]     In step  2 , prior knowledge of the image is used to define a region of interest (ROI) which is a subset of the image. This process can be done by whatever means, either automatic, semi-automatic, or even manual.  
         [0026]     In step  3  an analysis is performed on the frequency of occurrence of intensities within the ROI, and a range of frequencies is defined, again using prior knowledge.  
         [0027]     For example, without losing generality, we denote the image to be processed as f(x), where f(x) is the gray level at a pixel labelled x. It is further supposed that the processed image has L gray levels denoted by r i  where i is an integer in the range 0 to L- 1  and r 0 &lt;r 1 &lt; . . . r L-1 . It is also assumed that the object of interest has higher intensity values than the background. Suppose that due to prior knowledge or test we know that the proportion of the region of interest which is occupied by the object is in the percentage range per 0  to per 1 .  
         [0028]     Let h(i) denote the frequency of gray level r i , and let H(i) denote the cumulative frequency which is  
           ∑       i   ′     =   0     i     ⁢     h   ⁡     (     i   ′     )         ,       
 
 where i′ is an integer dummy index. Considering two values of i written as m and n, the frequency of intensities in the range r m  to r n  is  
         ∑       i   ′     =   n     n     ⁢       h   ⁡     (     i   ′     )       .         
 
 Thus, we can use per 0  to calculate a gray level r low , such that we are sure that all the pixels having lower intensity represent background. r low  can be written as:  
               r   low     =       min   i     ⁢       {     i   ❘       H   ⁡     (   i   )       ≥     per   0         }     .               (   1   )             
 
         [0029]     Similarly, we can use per 1  to calculate a gray level r high  such that we are sure that all pixels having higher intensity represent the object:  
               r   high     =       min   i     ⁢       {     i   ❘       H   ⁡     (   i   )       ≥     per   1         }     .               (   2   )             
 
         [0030]     In a step  4  of the method of  FIG. 1 , the threshold is selected using an algorithm which operates on the frequencies within the selected range from r low  to r high . The details of several ways in which this can be carried out within the scope of the invention are given below. Thus, a selected threshold is output in step  5 .  
         [0031]     Image binarisation is then performed using this threshold, to create an image in which all pixels (at least in the ROI) are classified into two classes. Further image processing steps may optionally be performed at this stage.  
         [0032]     We now turn to a discussion of three techniques by which step  4  can be carried out.  
         [0033]     1. Range-constrained Least Valley Detection Method (RCLVD)  
         [0034]     If the frequency range derived in step  3  is correctly estimated then it will include a valley in the frequency distribution of intensities. This valley separates the background and the object. Thus, valley detection can be exploited to select the threshold. This has the following steps:  
         [0035]     1) A frequency interval δh is specified.  
         [0036]     2) The gray level range [r low , r high ] is partitioned into K+1 intervals with an equal frequency range δh. For an interval labelled by integer index j, the lower end of its intensity range is denoted r 1   j  and the upper end is denoted r 2   j . Thus:  
                 r   1   0     =     r   low       ,       r   2   0     =       min   i     ⁢     {     i   ❘       H   ⁡     (   i   )       ≥     (       per   0     +     δ   ⁢           ⁢   h       )         }         ,                   r   1   1     =     r   2   0       ,       r   2   1     =       min   i     ⁢     {     i   ❘       H   ⁡     (   i   )       ≥     (       H   ⁡     (     r   1   1     )       +     δ   ⁢           ⁢   h       )         }         ,             …                 r   1   K     =     r   2     K   -   1         ,       r   2   K     =       min   i     ⁢     {     i   ❘       H   ⁡     (   i   )       ≥       (       H   ⁡     (     r   1   K     )       +     δ   ⁢           ⁢   h       )     .                             H   ⁢     (       r   1   K     +     δ   ⁢           ⁢   h       )       ≥       per   1     ⁢           ⁢   and   ⁢           ⁢     H   ⁡     (     r   1   K     )         &lt;       per   1     .               
 
         [0037]     3) The average frequency  h   j  for each of the intervals j is calculated given by  h   j =(H(r 2   j )−H(r 1   j))/(r   2   j −r 1   j )  
         [0038]     4) Let J denote the interval for which  h   j  is a minimum. The threshold of this RCVLD method, which is denoted θ RCVLD , may be selected to be any value in the range r 1   j  to r 2   j , such as θ RCVLD =(r 2   j +r 1   j )/2.  
         [0039]     2. Range-constrained Weighted Variance Method (RCWV)  
         [0040]     Let r k  fall within the range r low  to r high , and suppose that the pixels of the ROI are in two classes C 1  and C 2 , where C 1  is pixels of the background class and consists of pixels with gray levels r low  to r k , and C 2  is pixels of the object class and is composed of pixels with gray levels r k +1 to r high . The range-constrained weighted variance method maximises the “weighted between-class variance” defined as:  
             θ   RCWV     ⁡     (       W   1     ,     W   2       )       =       max     r   k       ⁢     (         Pr   ⁡     (     C   1     )       ⁢     D   ⁡     (     C   1     )       ⁢     W   1       +       Pr   ⁡     (     C   2     )       ⁢     D   ⁡     (     C   2     )       ⁢     W   2         )         ,       
 
 where W 1  and W 2  are two positive constants selected by the user and representing the weights of the two respective class variances, Pr(.) denotes the class probability, i.e.  
           Pr   ⁡     (     C   1     )       =       ∑     i   =     r   low         r   k       ⁢     h   ⁡     (   i   )           ,     
     ⁢       Pr   ⁡     (     C   2     )       =       ∑     i   =       r   k     +   1         r   high       ⁢     h   ⁡     (   i   )           ,       
 
 and D(C 1 ) and D(C 2 ) are given by:  
           D   ⁡     (     C   1     )       =           (       μ   0     -     μ   T       )     2     ⁢           ⁢   and   ⁢           ⁢     D   ⁡     (     C   2     )         =       (       μ   1     -     μ   T       )     2         ,       where   ⁢           ⁢     μ   T       =       ∑     i   =     r   low         r   high       ⁢     i   ×     h   ⁡     (   i   )                 
         μ   0     =         ∑     i   =     r   low         r   k       ⁢     i   ×     h   ⁡     (   i   )       ⁢           ⁢   and   ⁢           ⁢     μ   1         =       ∑     i   =       r   k     +   1         r   high       ⁢     i   ×       h   ⁡     (   i   )       .               
 
 When W 1  is bigger than W 2 , background homogeneity is emphasised. 
 
         [0041]     3. Range-constrained Fuzzy C-partition Thresholding Method (RCFCP)  
         [0042]     This third method is related to the technique used in [7], and the justification for it is as given there. In general terms, let A b /A 0  be the fuzzy sets of fuzzy events “background/object” (which denotes a fuzzy partition of the set {r low , . . . , r high  with a membership function μ A     b   /μ A     0    respectively). The probability of these fuzzy events are given by:  
           P   ⁡     (     A   i     )       =       ∑     j   =     r   low         r   high       ⁢         μ     A   i       ⁡     (   j   )       ×     h   j           ,       
 
 where A i ∈{A b ,A 0 , and the weighted entropy with this fuzzy partition can be calculated as:
 
 S ( W   1   ,W   2 )= W   1   ×P ( A   b )×log  P ( A   b )+ W   2   ×P ( A   0 )×log  P ( A   0 )
 
 where W 1  and W 2  are two positive constants, and log(.) is the natural logarithm. 
 
         [0043]     Let r low ≦a&lt;c≦r high . The membership functions can be defined as follows:  
           μ     A   b       ⁡     (   x   )       =     {               1   ,             r   low     ≤   x   ≤   a                 (     x   -   c     )     /     (     a   -   c     )             a   &lt;   x   &lt;   c             0         c   &lt;   x   ≤     r   high             ⁢     
     ⁢   and   ⁢     
     ⁢       μ     A   0       ⁡     (   x   )         =     {             1   ,             r   low     ≤   x   ≤   a                 (     x   -   a     )     /     (     c   -   a     )             a   &lt;   x   &lt;   c             0         c   &lt;   x   ≤     r   high             .               
 
         [0044]     The optimum parameters a* and c* are chosen to maximise the entropy S(W 1 , W 2 ), and the optimum threshold is θ RCFCP =(a*+c*)/2.  
         [0045]     Having now presented the steps of the embodiment in principle, we turn to an example of the embodiment in operation. This example uses the form of step  4  referred to above as RCLVD.  
         [0046]     The starting point of the method is the image shown in  FIG. 2 , an MR (Magnetic Resonance) image which is a T 1 -weighted or SPGR (spoiled gradient recalled acquisition) axial slice around the intercommissural plane. This image is input in step  1  of the method.  
         [0047]     In step  2  of the method, we calculate the pixels enclosed by the skull (i.e. find the ROI) using the following steps: the usual histogram-based thresholding method is used to binarise the axial slice; a morphological closing operation is used to connect small gaps; the largest connected component is identified; and the holes within the component are filled. The resulting ROI (the pixels enclosed by the skull) is shown in  FIG. 3 .  
         [0048]     In step  3 , the two percentages pero and per, are set as 14% and 28%. This selection is based on previous experiments and/or other prior knowledge.  
         [0049]     In step  4  of the method (RCLVD), we select the δh to be 1% (alternatively any value in the range 1% to 5% would be suitable).  FIG. 4  shows the histogram of frequencies in the ROI, and the calculated threshold θ RCLvD  is shown as the line indicated. This completes the procedure of the embodiment.  
         [0050]     The output threshold of the method is used as in conventional techniques to binarise the image. The binarised image is shown in  FIG. 5 .  
         [0051]     Although only a single embodiment of the invention has been described, many variations are possible within the scope of the invention as will be clear to a skilled reader.  
       REFERENCES  
       [0052]     The disclosure of the following references is incorporated herein by reference in their entirety: 
    [1] Otsu N., “A threshold selection method from gray-level histograms”, IEEE Transactions on Systems, Man and Cybernetics, 1979; 9: p62-66.     [2] Li C. H., Lee C. K., “Minimum cross entropy thresholding”, Pattern Recognition 1993; 26: p617-625.     [3] Kittler J., Illingworth J., “Minimum error thresholding”, Pattern Recognition 1986; 19: p41-47.     [4] Kapur J. N., Sahoo P. K., Wong A. K. C., “A new method for gray-level picture thresholding using the entropy of the histogram”, Computer Vision Graphics and Image Processing, 1985, 29; 273-285.     [5] Wong A. K. C. and Sahoo P. K., “A gray-levl threshold selection method based on maximum entropy principle”, IEEE Transactions on Systems, Manand Cybernetics, 1989; 19: p866-871.     [6] Saha P. K. and Udupa J. K., “Optimum image thresholding via class uncertainty and region homogeneity”, IEEE Transactions on Pattern Analysis and Machine Intelligence, 2001; 23: p689-706.     [7] Cheng H. D., Chen J., and Li J., “Threshold selection based on fuzzy c-partition entropy approach”, Pattern Recognition 1998; 31: p857-870.