Optic disc localization in retinal images using histogram matching

In this article, we propose a new method for localizing optic disc in retinal images. Localizing the optic disc and its center is the first step of most vessel segmentation, disease diagnostic, and retinal recognition algorithms. We use optic disc of the first four retinal images in DRIVE dataset to extract the histograms of each color component. Then, we calculate the average of histograms for each color as template for localizing the center of optic disc. The DRIVE, STARE, and a local dataset including 273 retinal images are used to evaluate the proposed algorithm. The success rate was 100, 91.36, and 98.9%, respectively.


Introduction
Retina is the innermost layer of the eye which can be visualized using adequate apparatus such as fundus camera. The two main structures used in retinal image analysis are blood vessels and optic disc. Optic disc is the brightest region in the retinal image and the blood vessels originate from its center [1]. Optic disc is a key reference for recognition algorithms [2,3], blood vessels segmentation [4], and diagnosing some diseases such as diabetes [5]. Histogram is the main character of each image and histogrambased methods are used as the first step of most preprocessing methods to improve the contrast and illumination of retina images. One of the main drawbacks of uneven illumination in retina images and their poor quality is the inability to analyze the optic disc. Applying illumination equalization (histogram equalization, histogram specification, and other normalization methods) as preprocessing methods to retina images considerably improves the contrast, and illumination for further analysis tasks such as optic disc localization and vessel segmentation [6,7]. In this article, we propose a new method based on the histograms of some optic discs extracted from retinal images. For this purpose, we extract the optic disc of the first four retinal images in DRIVE dataset. Then, we calculate the average of histograms for each color component as template to localize the center of optic disc.
The rest of this article is organized as follows. "Review of previous methods" section is devoted to review the latest proposed methods for optic disc localization. In "Anatomy of the retina" section, we briefly review the anatomy of retina. "Method" section presents the proposed method for optic disc localization. Experimental results are given in "Results" section. Finally, "Conclusion and future work" section is devoted to concluding remarks.

Review of previous methods
Osareh [8] proposed a method based on template matching for localizing the center of optic disc. In this algorithm, some of retinal images in dataset were used to create a template and the correlation between each image and template is computed. The point which has the maximum correlation value is selected as the center of optic disc.
Youssif et al. [9] used directional pattern of the retinal blood vessels to localize the center of optic disc. Hence, a simple matched filter was proposed to match the direction of the vessels at the optic disc vicinity. The retinal vessels were segmented using a simple and standard 2D Gaussian matched filter. Consequently, vessels' direction map of the segmented retinal vessels was obtained using the same segmentation algorithm. Then, the segmented vessels were thinned and filtered using local intensity to represent the optic disc center candidates. The Gaussian matched filter was resized in four different sizes, and the difference between the output of the matched filter and the vessels' directions was measured. The minimum difference provided an estimate of the optic disc-center coordinates.
Li and Chutatape [10] proposed a new method to localize optic disc center. The candidate regions were first determined by clustering the brightest pixels in retinal images. This strategy can only work when there is no abnormality in the retina image. Principal component analysis was applied to these candidate regions. The minimum distance between the original retinal image and its projection onto disk space was located as the center of optic disc.
Rangayyan et al. [11,12] proposed two different methods. In the first method, optic disc center was localized based on the property that it appears as the focal point of the blood vessels in retina mage. The method includes detection of the blood vessels using Gabor filters and detection of peaks in the node map via phase portrait analysis. In the second method, edge detection using the Sobel operators and detection of circles using the Hough transform were employed to localize optic disc and its center.
Aquino et al. [13] used two independent methodologies to detect optic disc in retina images. Location methodology obtains a pixel that belongs to the optic disc using image contrast analysis and structural filtering techniques. Then, a boundary segmentation methodology estimates a circular approximation of the optic disc boundary by applying mathematical morphology, edge detection techniques, and the circular Hough transform.
Siddalingaswamy and Gopalakrishna Prabhu [14] proposed a new approach for the automatic localization and accurate boundary detection of the optic disc. Iterative thresholding method followed by connected component analysis was employed to localize the approximate center of the optic disc. Then, geometric model based on implicit active contour model was applied to find the exact boundary of the optic disc.
Foracchia et al. [15] presented a new technique for localizing the optic disc center in retinal images. The method was based on the preliminary detection of the main retinal vessels. All retinal vessels originate from the optic disc and their path follows a similar directional pattern (parabolic course) in all images. To describe the general direction of retinal vessels at any given position in the image, a geometrical parametric model was proposed, where two of the model parameters are the coordinates of the optic disc.
Carmona et al. [16] used genetic algorithm method to obtain an ellipse approximating the optic disc in retinal images. A set of hypothesis points were initially obtained that exhibited geometric properties and intensity levels similar to the optic disc contour pixels. Then, a genetic algorithm was used to find an ellipse containing the maximum number of hypothesis points in an offset of its perimeter, considering some constraints.
A number of other interesting algorithms can be found in the literature that used vessel segmentation results for optic disc localization [17][18][19][20].
Some methods are based on the Hough transform which is capable of finding geometric shapes. Therefore, the circular shape of optic disc was detected using Hough transform and other algorithms such as thresholding and morphological operations [21][22][23][24][25]. 1-Superior temporal blood vessels 2-Superior nasal blood vessels 3-Fovea 4-Optic disc 5-Inferior temporal blood vessels 6-Inferior nasal blood vessels

Method
Most of the methods for localizing optic disc fail when pathological regions exist in retina images [8,10]. Some other algorithms suffer from high computational cost [9,[11][12][13][14][15][16][17][18][19][20]. Here, a new robust method for localizing the center of optic disc in presence of pathological regions is proposed. Since in this method preprocessing algorithms such as segmentation are not used, the computational cost is drastically reduced with respect to some counterparts.  In this method, similar to [8], we use a number of retinal images to create a template for optic disc. However, instead of creating an image as template, we construct three histograms as template, each corresponding to one color component. At the first step to decrease the effect of noise, we apply an average filter with the size of 6 × 6 pixels to retina images. Then, we use a window with the typical size of the optic disc (80 × 80 pixels) to extract the optic disc of each retinal image. In the next step, we separate color components (red, blue, and green) of each optic disc to obtain the histogram of each color component. Finally, the mean histogram of each color component for all retinal image samples is calculated as template. Histogram is a graph showing the number of pixels at each different intensity value found in an image. As illustrated before, we use the histogram of each three channels (red, green, and blue) as template for optic disc localization. Then, to decrease the effect of pathological regions and exudates that are high-bright regions like optic disc, we use the histogram of pixels which has the intensity value lower than 200. Therefore, we decrease the effect of high intensity regions that are common in optic disc, pathological regions and exudates and the role of vessels for optic disc localization will increase.

Template matching
Up to now, we determined three histograms as template for localizing the center of optic disc. For localizing the center of optic disc, at first step to decrease the effect of noise an average filter with the size of 6 × 6 pixels is applied to retina image. Then, an 80 × 80 pixels window is moved through retinal image. In each moving window, we separate the channels (red, blue, and green) and obtain the histogram of each channel. Then, we calculate the correlation between the histogram of each channel in the moving window and the histograms of its corresponding channel in template. For this purpose, we can use correlation or cross-correlation function to obtain the similarity of the two histograms; however, the optic disc centers obtained using these methods are not accurate. The function used for correlation between two histograms is expressed in the following equation: where a and b are two histograms that we want to calculate their correlation and c is the result of the correlation. Therefore, if the two histograms (a and b) are similar Σ i (a i − b i ) 2 ≈ 0, and c ≈ 1, else and c<<1. Therefore, using Equation (1) we can calculate the correlation between two histograms and the result of correlation is in the range of [0 1].
For each moving window, we obtain three values as the results of correlation between the histograms. The result of histograms matching is computed as the weighted sum of the three obtained values: where (i,j) is the center of moving window. c r , c g , and c b are the results of correlation for three channels (red, green, and blue) and t r , t g , and t b are weights used for each channel. In Equation (2), we can use different weights for c r , c g , and c b . The green channel has the highest weight because the contrast of the green channel is higher than red and blue channels [27]. In some retinal images, blue channel is noisy; therefore, to decrease the effect of blue channel on our localizing method, we determine the lowest weight for blue channel. The best weights that result high accuracy rate for optic disc localizing method are t r =0.5, t g =2, and t b =1. To localize the center of optic disc, we apply thresholding on the correlation function C(i,j). For finding the best threshold, we did a global scanning of different values and the best equation to determine the threshold (Th) was obtained as follows.
where max(C) is the element of C with the maximum value. Therefore, the threshold value for each image is half of the maximum value of the correlation function. The center of gravity of the binary image obtained from thresholding is considered as the center of optic disc.

Results
We used the optic disc of the first four retinal images in DRIVE dataset [28] to obtain their histograms as template. The first four retinal images from DRIVE dataset that have been used to extract the histograms of their optic disc are shown in Figure 3. The mean histogram of each color component for the optic discs of these four retina images is calculated as template.
In Figure 4, we can see three histograms obtained as template for localizing the center of optic disc.
The proposed method was applied on a dataset including 40 retina images from DRIVE dataset (565 × 584 pixels), 81 retinal images from STARE dataset (605 × 700 pixels) [29], and 273 retinal images from a local dataset (720 × 576 pixels). Retina images in local dataset were captured by a Canon CR5 in Razi clinic a from normal and abnormal eyes. The success rate was 100, 91.36, and 98.9% for these three datasets, respectively. In Figure 5, some retina images in datasets and the results of applying threshold before determining the center of optic disc are shown.
In Figures 6, 7, and 8, the results of the proposed method for some retinal images in our datasets are shown. In presence of abnormality in the eye (pathological regions or exudates), using the histogram analysis of optic disc is more effective. Pathological regions, exudates, and optic disc are bright regions in the retina images. Therefore, methods such as template matching [6] or methods which are based on the segmentation results of blood vessels [9,10,[15][16][17] fail to localize the center of optic disc in presence of pathological regions and exudates in retina image. Figures 6, 7, and 8 show the result of the proposed method on normal retina image and retina images with pathological regions and exudates. Despite the existence of dark hemorrhages or bright exudates and pathological regions, the results of the proposed method are satisfactory and it shows the effectiveness of the proposed method for localizing the center of optic disc. In Figure 9, some retinal images with incorrectly detected optic disc center are shown.
In Figure 9a, there is not any vessel in vicinity of optic disc and the characteristic of optic disc-like brightness and high number of vessels in vicinity of optic disc cannot be seen; therefore, our proposed method failed to localize the optic disc center. For the retinal image in Figure 9b, optic disc is in the corner of image and there is really no vessel in optic disc. Therefore, our proposed method failed to localize the optic disc center. Therefore, in situation like Figure 9a that there are not any vessels in optic disc vicinity or in situation that we have pathological region with high number of vessels, our proposed method failed to localize optic disc center. Comparing the proposed method and its counterparts, Table 1 shows the effectiveness of the proposed method. The datasets used in the proposed method are larger than datasets used in other methods and thanks to avoiding pre-processing algorithms such as segmentation, the proposed method takes less computation time in comparison to the counterpart methods. In Table 1, the running times of some methods are indicated. The system's configuration used in each algorithm was different. Therefore, simple comparison of the running time of different methods does not appear to be correct. A parameter which determines the running time of different methods is computational complexity. In methods that use segmentation, a large number of operations per pixel are needed and consequently these methods require high order of computational cost. In our proposed method, no preprocessing such as segmentation is required; therefore, it takes less computational time. To clarify that system's configuration does not have high effect we used a computer with system's configuration of Intel Core 2 Duo, 1.7 GHz, and 512 MB RAM and the average running time for this configuration was 32.5 s.
In Figure 10, we can see the sensitivity of the proposed method to the threshold (Th) applied in the final step.  Therefore, the average distance between the estimated and manually identified optic disc centers based on the different amount of thresholdings is plotted in Figure 10. From Figure 10, we can understand the effect of thresholding on the average distance between the estimated and manually identified optic disc centers. Therefore, the best threshold value is half of the maximum value of the correlation function obtained before applying threshold.

Conclusion and future work
In this article, we presented a new method for localizing the center of optic disc. The average distance between the estimated and the manually identified optic disc centers is 17 and 26 pixels in [9] and 23.2 and 119 pixels in [11] for DRIVE and STARE datasets, respectively. These values in the proposed method are 15.9, 11.4, and 8.9 pixels for DRIVE, STARE, and local datasets, respectively. Therefore, the estimated optic disc centers obtained using the proposed method are more accurate in comparison to other algorithms such as methods introduced in [9,11,12]. In this article, we used the histograms of some optic discs and in presence of pathological regions and exudates in retinal images, and we could determine the center of optic disc correctly. Most of the counterpart methods perform well when there are no pathological regions or exudates in retinal images. In this article, the first four retinal images in DRIVE dataset were used to obtain their histograms as template, using more retinal images such as retinal images in STARE and local datasets may improve the effectiveness of our proposed method. Also to decrease the running time of our proposed method, we can combine our proposed method with other methods. For example, as we know template matching method proposed in [8] fails in situation like pathological regions and exudates exit and also the accuracy of template matching method for localizing optic disc center of retina images without any pathological regions and exudates is low. Therefore, we can use template matching method for retina images to obtain candidate regions that probability of existing optic disc in them is more than other regions in retina images. Then, instead of applying our proposed method on the whole of retina images, we apply it to candidate regions to obtain optic disc center. Therefore, the running time of our proposed method will considerably decrease. In future work, we use optic disc center obtained as the first step for localizing the boundary of optic disc and also we can use the optic disc center for recognition algorithm in our future research for human recognition based on the retinal images. Distance between the estimated and manually optic disc centers (Pixel) Figure 10 Sensitivity of the proposed method to thresholding based on average distance between the estimated and manually identified optic disc centers for all retina images in datasets.