CN101667289A - Retinal image segmentation method based on NSCT feature extraction and supervised classification - Google Patents

Abstract

The invention discloses a retinal image segmentation method based on NSCT feature extraction and supervised classification, which relates to medical image processing. The steps are: (1) Use the red component of the retinal training image and the retinal image to be segmented to obtain the region of interest; (2) For the retinal training image and the green component of the retinal image to be segmented, iterate the edge of the region of interest Extend; (3) perform NSCT transformation on the expanded image respectively, and decompose it into i layers; (4) extract one-dimensional features by using the j direction subband coefficients of each layer, and extract the feature vector layer by layer, and normalize (5) establish training samples for the eigenvectors of the normalized retinal training images; (6) select a classifier, and use the training samples to train the classifier, and use the normalized eigenvectors of the retinal images to be segmented Input the classifier to segment the retinal image to be segmented. The invention has the advantages of clear image segmentation edge and high accuracy, and is used for retinal segmentation of medical images.

Description

Retinal Image Segmentation Method Based on NSCT Feature Extraction and Supervised Classification

technical field

The invention belongs to the technical field of image processing, relates to the application of retinal detection, and can be used to extract retinal blood vessels from retinal images in medicine.

Background technique

The development of medicine is closely related to human health, so digital image processing technology has aroused great interest in biomedical circles from the very beginning. As early as the end of the 1970s, literature statistics pointed out that a very wide application of image processing is medical image processing. In medicine, whether it is in basic disciplines or clinical applications, it is a field with many types of image processing. However, due to the technical difficulty of medical image processing, it is difficult to achieve clinical practicality in many processes. In recent years, with the reduction of the cost of digital image processing equipment, the use of digital image processing technology to improve the quality of various medical images has reached the practical stage.

Cardiovascular and cerebrovascular diseases such as hypertension, cerebrovascular sclerosis, and coronary arteriosclerosis are the main causes of death and disability in the elderly in my country. The tissue level of damage caused by such diseases is firstly the changes in the microcirculation and microvascular levels. Fundus retinal microvessels are the only deep microvessels in the human body that can be directly observed non-invasively. Its changes are closely related to the course, severity and prognosis of diseases such as hypertension. Diseases such as high blood pressure, diabetes, and arteriosclerosis can be found through the examination of the retinal vascular system. In retinal images, vessels are dominant and structurally stable, thus reliable vessel extraction is a prerequisite for retinal image analysis and processing.

There are two main types of existing retinal segmentation, one is to use Gaussian filter, Hessian matrix, etc. to enhance the contrast between blood vessels and the background, and then perform further operations, or perform threshold processing, or adopt region growth. In the article Adaptive local thresholding by verification-basedmultithreshold probing with application to vessel detection in retinal images in 2003, Jiang et al. proposed an adaptive local threshold segmentation based on a verification-based multi-threshold detection scheme, that is, first through the assumed threshold. Binarization yields a hypothetical target, and a validation procedure decides whether to accept or reject that target.

The other is to use a classifier to segment blood vessels. In this method, the key issue is how to build the eigenvectors. In the article Ridge based vessel segmentation in color images of theretina in 2004, J.Staal et al. proposed a method of first using ridge detection to obtain feature vectors, and then classifying and segmenting blood vessels. In the article Retinal Vessel Segmentation Using the 2-D Gabor Wavelet and Supervised Classification in 2006, Soares et al. et al. proposed a method of feature extraction based on the multi-scale and multi-directional characteristics of Gabor wavelet, and then segmented by supervised classification.

In retinal image blood vessel segmentation, the most important thing is to increase the detection rate of blood vessels. However, due to the global unbalanced gray scale of retinal images, poor contrast between blood vessels and background, noise and the influence of various lesion areas, the above methods are not effective for blood vessels. There are certain errors in the segmentation results, and the accuracy needs to be further improved.

Contents of the invention

The object of the invention is to provide a retinal image segmentation method based on Nonsubsampled Contourlet Transform (NSCT) and supervised classification. It overcomes the shortcomings of the existing technology and further improves the accuracy of segmentation.

To achieve the above object, the technical solution of the present invention comprises the following steps:

1. Feature extraction steps

(1) For the retinal training image and the retinal image to be segmented, use its red component to obtain the region of interest ROI;

(2) For the retinal training image and the green component of the retinal image to be segmented, iteratively expand the edge of the region of interest;

(3) Perform NSCT transformation on the expanded retinal training image and the retinal image to be segmented respectively, decompose it into i layers, and each layer has subband coefficients in j directions;

(4) extract one-dimensional features by using j direction sub-band coefficients in each layer, and extract features layer by layer to form feature vectors and normalize them;

2. Classifier training and segmentation steps

1) Establish training samples for the feature vectors of the normalized retinal training images;

2) Select a classifier and use the training samples to train the classifier, input the normalized feature vector of the retinal image to be segmented into the classifier, and segment the retinal image to be segmented.

In the above-mentioned retinal image segmentation method, the use of j direction subband coefficients in each layer to extract one-dimensional features described in step (3), and extracting the feature vector layer by layer, is carried out as follows:

(3a) For the j direction subband coefficients decomposed in each layer, the comparison result of the blood vessel gray level in the retinal image and the background gray level is used as the selected feature object. If the blood vessel gray level in the retinal image is smaller than the background gray level, select the smallest one The coefficient is used as a feature. If the grayscale of blood vessels in the retinal image is greater than the background grayscale, the largest coefficient is selected as a feature;

(3b) According to step (3a), the same feature extraction operation is performed layer by layer, and the obtained features are formed into feature vectors;

(3c) Add the gray value of the retinal training image and the green component of the retinal image to be segmented into the feature vector as a one-dimensional feature to obtain the final feature vector v.

Based on the multi-scale, multi-direction, and translation invariance of NSCT, the present invention proposes a new feature extraction method in combination with the specific expression form of NSCT transformation coefficients of retinal images, that is, using feature descriptions extracted on different scales Retinal blood vessels of different widths are combined with features obtained on several scales to form feature vectors, so the segmentation of retinal image blood vessel edges is more accurate; at the same time, because the present invention adopts the supervised classification method of training first and then classifying, in the case of supervision The retinal image is segmented, so the error of the retinal image segmentation result is small. Simulation results show that the present invention has better segmentation accuracy than existing retinal image segmentation methods.

Description of drawings

Fig. 1 is the realization flow schematic diagram of the present invention;

Fig. 2 is a schematic diagram of the results of various intermediate steps for obtaining the ROI region of the retinal image in the present invention;

Fig. 3 is a schematic diagram of the result of ROI edge expansion on the retinal green component;

Fig. 4 is a schematic diagram of feature extraction principle;

Figure 5 is two images selected from the results of experiments on the DRIVE database;

Figure 6 is two images selected from the results of experiments on the STARE database;

Figure 7 is a ROC curve diagram showing the results of testing the DRIVE database;

Fig. 8 is a comparison diagram of the simulation results of the method of the present invention and the existing Soares et al. method.

Specific implementation method

With reference to Fig. 1, the concrete realization process of the present invention is as follows:

Step 1, for the retinal training image and the retinal image to be segmented, use its red component to obtain the ROI of the region of interest. The specific steps refer to Figure 2 as follows:

(1.1) Divide the red component of the retinal image shown in Figure 2(a) by 255, as shown in Figure 2(b), and perform edge detection on Figure 2(b) through the Gaussian filter LOG to obtain Figure 2(c );

(1.2) Figure 2(c) is first expanded and then corroded, so that the edge fractures are connected, and Figure 2(d) is obtained;

(1.3) In Figure 2(d), add a contour along the edge of the image;

(1.4) Determine the external area through the threshold. In Figure 2(b), find the maximum value of its grayscale max red, and mark the point with a grayscale smaller than max red×0.15 as 1, thus obtaining a binary image. Then remove from it the part marked as 1 with an area smaller than 10 pixels, resulting in Figure 2(e);

(1.5) In Figure 2(e), remove the contour added in step (3), and then fill Figure 2(c), the filling starting point is the point marked 1 in Figure 2(e), as shown in Figure 2 (f);

(1.6) Perform the inversion operation on Figure 2(f), and then corrode it to obtain Figure 2(g);

(1.7) Perform the inversion operation on Figure 2(g), then remove the object marked as 1 with an area smaller than 5000 pixels, and then invert it to achieve the purpose of filling the missing area, the result is shown in Figure 2(h);

(1.8) Perform an open operation on Figure 2(h) to remove false contours, and then remove objects marked 1 with an area smaller than 5000 pixels, and finally obtain the ROI of the region of interest, as shown in Figure 2(i).

Step 2: Extend the ROI edge of the retinal training image and the green component of the retinal image to be segmented to eliminate the strong contrast between the retinal base and the area outside the aperture, and avoid a large number of misjudgments at the edge of the aperture. Specific steps are as follows:

(2.1) Find the pixels of the outer boundary of the ROI, these pixels are located outside the ROI, but are 4 neighbors with the pixels in the region;

(2.2) For the outer boundary pixel obtained in step (2.1), its gray value is set as the mean value of the smaller area that is centered on the pixel and belongs to the ROI, and the area size can be selected as 5 × 5;

(2.3) Add the outer boundary pixels obtained in step (2.2) to the ROI, perform a certain number of iterations, usually 80 iterations.

Fig. 3 is a schematic diagram of the result of ROI edge expansion on the green component of the retinal image. Figure 3(a) is the green component of the retinal image shown in Figure 2(a), and Figure 3(b) is the result of ROI edge extension on Figure 3(a).

Step 3: Perform Nonsubsampled Contourlet Transfrom (NSCT) transformation on the expanded retinal training image and the retinal image to be segmented respectively.

NSCT is a multi-resolution, local, and multi-directional image representation method. It not only inherits the multi-resolution time-frequency analysis characteristics of wavelet transform, but also has good anisotropy characteristics, and can express smooth curves with fewer coefficients than wavelet. NSCT consists of a non-subsampled multi-stage decomposition and a non-subsampled multi-stage directional filter with translation invariance. The non-subsampling multilevel decomposition is implemented by non-subsampling Laplacian decomposition, and the non-subsampling multilevel directional filter is realized by non-subsampling directional filter bank. NSCT decomposition is performed on the expanded retinal training image and the retinal image to be segmented. The present invention is divided into 4 layers, and each layer has 8 directions.

Step 4, perform feature extraction on the retinal image. The principle of the feature is: in the NSCT domain of the image, select the coefficient map in any direction on the third layer, as shown in Figure 4(a). In this coefficient map, the edge of the blood vessel is the zero crossing point, that is, the place where the sign of the coefficient changes. Assuming that the width of the blood vessel is very large, for example, hundreds of pixels wide, the shape of the coefficient at the blood vessel in the direction perpendicular to the blood vessel is shown in Figure 4(b). At this time, if the width of the blood vessel keeps decreasing and the two zero-crossing points keep approaching, the final figure 4(b) tends to figure 4(c). In the retinal image, the width of the blood vessel is only a dozen pixels, so its NSCT coefficient tends to be in the form of Figure 4(c).

The process of feature extraction is as follows:

(4.1) Select the second layer of retinal image NSCT decomposition, and the third layer and the fourth layer generate feature vectors;

(4.2) For the 8 direction subband coefficients decomposed in each layer, the comparison result of the blood vessel gray level in the retinal image and the background gray level is used as the selected feature object. If the blood vessel gray level in the retinal image is smaller than the background gray level, select the smallest one Coefficients are used as features, as shown in Figure 4(d); if the grayscale of blood vessels in the retinal image is greater than the grayscale of the background, the largest coefficient among them is selected as a feature.

(4.3) According to step (4.2), the same feature extraction operation is carried out layer by layer, and the obtained features are formed into feature vectors;

(4.4) Add the gray value of the retinal training image and the green component of the retinal image to be segmented into the feature vector as a one-dimensional feature to obtain the final feature vector v={v _i |i=1, 2, 3, 4}.

Step 5, select a classifier, use training samples to train the classifier, input the normalized feature vector of the retinal image to be segmented into the classifier, and segment the retinal image to be segmented.

The present invention uses a Bayesian classifier for classification, and its conditional probability density function is represented by a Gaussian mixture model, and the pixels in the retinal image are divided into two categories: C ₁ = {vascular pixel}, C ₂ = {background pixel};

The Bayesian decision rule is as follows:

if p(v|C ₁ )P(C ₁ )＞p(v|C ₂ )P(C ₂ )

then decide C ₁ ; otherwise C ₂ ;

Among them, p(v|C _i ) is the class conditional probability density function, also known as the likelihood, p(C _i ) is the prior probability of C _i , and v is the feature vector;

Estimated p(C _i )=N _i /N, that is, the sample ratio of C _i in the training set, the class conditional probability density p(v|C _i ) is represented by a Gaussian mixture model, which is obtained by linear combination of several Gaussian functions, namely

$\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; v</span>}<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; C</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; j</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\overset{{\mathrm{<span\; class="notranslate">\; k</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; v</span>}<span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; j</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; C</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; P</span>}}_{\mathrm{<span\; class="notranslate">\; ij</span>}}$

Among them, k _i is the number of Gaussian functions representing p(v|C _i ), p(v|j, C _i ) is a multidimensional Gaussian distribution, and P _ij is a weight;

For each category C _i , given k _i Gaussian distributions, estimate the parameters and weights of each Gaussian distribution through the expectation-maximization EM algorithm;

When using a classifier to segment an image to be segmented, first obtain a probability map, where the value of each pixel is

P _v =p(C ₁ |v)/(p(C ₂ |v)+p(C ₁ |v));

Select the threshold p=0.5 to divide the probability map pixels into two categories, and finally complete the classification.

The effect of the present invention is further illustrated by the following simulation images and data.

1. Simulation image

The retinal images used in the present invention come from two public color image databases, DRIVE and STARE. The DRIVE database consists of 40 color images, of which 7 have lesions and have their manual segmentation results. The 40 images were divided into two groups: a training set and a test set, each containing 20 images, and the test set contained 3 lesion images. The images were manually segmented by personnel trained by professional ophthalmologists. The results of the training set graph being split by the first group of people are saved in setA. The results of the segmentation of the test set image by the first group of personnel are stored in setA, and at the same time, the test set image is also segmented by the second group of personnel and stored in set B. In set A, 12.7% of the pixels are labeled as blood vessels, while in set B, 12.3% of the pixels are labeled as blood vessels.

The STARE database consists of 20 digitized slides. There are lesions in ten of them. These images are segmented separately by two observers and stored in ah set and vk set respectively. In set A (ah set), that is, the results of the first observer, 10.4% of the pixels were labeled as blood vessels, while in set B (vk set), that is, the results of the second observer, 14.9% of the pixels were labeled as blood vessels. The difference between the results of these two observer segmentations is quite large. The second observer marked more extremely thin blood vessels than the first observer. This shows that the first observer is more conservative than the second observer.

For the DRIVE database, the training samples are generated from 20 labeled training images in the training set, and then the classifier is trained to classify the 20 images to be segmented in the test set.

For the 20 pictures in the STARE database, the present invention uses the leave-one-out method to test, that is, one image is used as the image to be segmented, and the other images are used as training images.

When these two databases take training sample sets, they all use set A. Since the training sample size is quite large, in all experiments, 1 million samples are randomly selected to train the classifier. In the test, for the class conditional probability density of the GMM classifier, the number of Gaussian models of the blood vessel and the background is taken as the same value k=k ₁ =k ₂ .

Second, the objective evaluation index

Two quantitative metrics are given: the area under the ROC curve Az and the accuracy rate. The ROC curve is the ratio of the correct detection rate to the false detection rate when the threshold p changes from 0 to 1 on the probability map P _v obtained by the classifier classification. The correct detection rate is the ratio of the number of pixels belonging to real blood vessels among the blood vessel pixels separated by the classifier to the total number of actual blood vessel pixels. The false detection rate is the ratio of the number of pixels that classify non-vessel pixels as vessels to the total number of actual non-vessel pixels. For the ROC curve, the closer the curve is to the upper left corner, the better the performance of the method. In other words, the closer the value of the area under the ROC curve Az is to 1, the better the performance of the method. Accuracy is the ratio of the total number of correctly classified pixels to the total number of retinal image pixels regardless of whether the pixel is vascular or non-vascular.

3. Simulation results and analysis

Fig. 5 is the results of two images selected from the pair DRIVE database when k=20, and their manual segmentation images. Figure 5(a) is the probability map, Figure 5(b) is the segmentation result, Figure 5(c) is the manual segmentation map in Set A, and Figure 5(d) is the manual segmentation map in Set B.

Fig. 6 is the experimental result of two images selected from the STARE database when k=20, and their artificially segmented images. The image in the second row is derived from an image of a diseased retina, and the image in the first row is derived from an image of a normal retina. Figure 6(a) is the probability map, Figure 6(b) is the segmentation result, Figure 6(c) is the manual segmentation map in set A, and Figure 6(d) is the manual segmentation map in Set B.

Figure 7(a) shows the ROC curve of the test results on the DRIVE database. Figure 7(b) shows the ROC curve of the test results on the STARE database. In order to compare the difference between the two different artificial segmentation results of set A and set B, the correct detection rate and error rate obtained when the probability map is set A and set B as the standard segmentation results respectively when p=0.5 are given. ratio of detection rates.

It can be seen intuitively from Fig. 6 and Fig. 7 that the segmentation method of the present invention can achieve better results.

Table 1 shows the performance comparison between the method of the present invention and the existing methods.

Experiments were carried out under the conditions of k=10, 15, and 20, and the experimental results are shown in Table 1. At the same time, the results of the present invention are compared with the methods extracted by Jiang et al., Staal et al.

Table 1 Performance comparison of different methods

As can be seen from the results in Table 1, the method of the present invention has achieved very good results for both DRIVE and STARE databases. By comparison, we find out that the effect of the inventive method is equivalent to that of Soares et al., and is better than the method of Jiang et al.Staal et al. This illustrates the effectiveness of the feature extracted by NSCT given by the present invention, and a new retinal feature is given.

Fig. 8 compares the result of the segmentation of an image in the DRIVE database by the present invention and the result of Soares et al. Fig. 8 (a) is the artificial segmentation result of set A, Fig. 8 (b) is the blood vessel segmentation result of the present invention, Fig. 8 (c) the result of Soares et al., Fig. 8 (d) is the result of Soares et al. minus The results of blood vessel segmentation in this paper, Figure 8(e) is the result of Soares et al. minus the manual segmentation result of set A, and Figure 8(f) is the segmentation result of this paper minus the manual segmentation result of set A. In Figure 8(d), the gray area represents the same part of the results of Soares et al. and the thick blood vessel detection results in this paper, the white area represents the area detected by Soares et al. method but not detected by the thick blood vessel method in this paper, and the black area Indicates the area detected by the method of this paper but not detected by the method of Soares et al..

As can be seen from Fig. 8(e) and Fig. 8(f), the result of segmentation by the existing Soares et al. method is generally thicker than the blood vessels in the real result, while the result of segmentation by the method of the present invention has an edge positioning effect Better than the method of Soares et al. This shows that the positioning of the edges of thick vessels in the Soares et al. method is not accurate. It can also be seen from Fig. 8(d) that the difference between the segmentation result of the method of the present invention and the result of Soares et al. is mainly concentrated at the edge of the thick blood vessel.

Comparing the results of the red oval area in Figure 8(a) in Figure 8(b)(c), it can be found that in the results of Soares et al., the space between two adjacent thick blood vessels is also judged as a blood vessel, while This problem does not exist in the results of the method of the present invention.

Claims (2)

1. A retinal image segmentation method based on NSCT feature extraction and supervised classification, comprising the steps of: 1. Feature extraction steps (1) For the retinal training image and the retinal image to be segmented, use its red component to obtain its region of interest; (2) For the retinal training image and the green component of the retinal image to be segmented, iteratively expand the edge of the region of interest; (3) Perform NSCT transformation on the expanded retinal training image and the retinal image to be segmented respectively, decompose it into i layers, and each layer has subband coefficients in j directions; (4) extract one-dimensional features by using j direction sub-band coefficients in each layer, and extract features layer by layer to form feature vectors and normalize them; 2. Classifier training and segmentation steps 1) Establish training samples for the feature vectors of the normalized retinal training images; 2) Select a classifier and use the training samples to train the classifier, input the normalized feature vector of the retinal image to be segmented into the classifier, and segment the retinal image to be segmented. 2. The method according to claim 1, wherein said step (3) utilizes j direction subband coefficients of each layer to extract one-dimensional features, and extracts layer by layer to form a feature vector, as follows: (3a) For the j direction subband coefficients decomposed in each layer, the comparison result of the blood vessel gray level in the retinal image and the background gray level is used as the selected feature object. If the blood vessel gray level in the retinal image is smaller than the background gray level, select the smallest one The coefficient is used as a feature. If the grayscale of blood vessels in the retinal image is greater than the background grayscale, the largest coefficient is selected as a feature; (3b) According to step (3a), the same feature extraction operation is performed layer by layer, and the obtained features are formed into feature vectors; (3c) Add the gray value of the retinal training image and the green component of the retinal image to be segmented into the feature vector as a one-dimensional feature to obtain the final feature vector v.

Images