1.

Introduction

Clinical research on dry-eye syndrome is currently receiving much attention.^1–12 In particular, in the study of meibomian gland dysfunction, which is a major cause of dry-eye,¹ gland dysfunction can be diagnosed by directly observing the morphology of the meibomian glands using meibography. Recently,^2–4 a noncontact, patient-friendly infrared meibography method was developed that is fast and greatly reduced discomfort to the patient. This advance allows meibography to be used routinely in clinical research.

When analyzing meibography images, it is important to quantify the area of meibomian gland loss. Grades are then assigned to the image based on the area of gland loss, which in turn serve as an indicator of meibomian gland dysfunction. The loss of meibomian glands in the upper and lower eyelids have been studied by Arita et al.,² Pult et al.,⁵and Srinivasan et al.⁹ Pult and Reide-Pult have also assessed the subjective and objective grading of meibography images,⁸ and Pult and co-workers have also demonstrated that precise calculation of the area of meibomian gland loss using a computer software can lead to greatly improved repeatability in the application of grading scales.^5–7 In addition, other meibomian gland features such as the width and tortuosity of the glands can be accurately computed and analyzed. Similar progress have also been reported by Srinivasan et al.⁹ These recent works demonstrated the advantages of a computerized approach to analyzing meibomian glands in meibography images.

In these recent works on computerized classification of meibomian gland loss, the images were analyzed using the image editing software ImageJ (National Institute of Health; http://imagej.nih.gov/ij). The user needs to be involved in identifying the gland region on the image. In other words, the user needs to tell the software where the glands are. If a large number of images is involved, this process is tedious and time consuming. Also, different examiners may draw the gland region differently, leading to inter-observer variability. Hence, it is desirable to have a method of processing meibography images that is fast, less laborious for the clinician, and is less dependent on the subjectiveness particular examiner.

A fully automated, computational approach of analyzing meibography images can address many of the difficulties mentioned, and has the potential to aid ophthalmologists in the study of meibomian glands using meibography. The purpose of this paper is to take a first step in this direction by presenting algorithms that can detect meibomian glands and classify meibography images with minimal user input. Image editing software such as ImageJ and Photoshop that require input from the user are not used. Instead, algorithms that cater specifically to the detection of meibomian glands are specially developed.

Meibography images exhibit a wide range of gland morphologies. To be of practical help to ophthalmology, one should try to sample from a wide range of images when developing the detection and classification algorithm. In particular, of most concern to clinicians is identification of images with winding glands or other complex gland patterns, because these represent cases which are intermediate between the healthy and unhealthy, and hence require the most clinical attention. However, although identification of complex gland pattern is easy for trained experts, it is a challenging task for an algorithm without help from a human user. Hence, this paper focuses on the easiest healthy and unhealthy (criteria defined in Sec. 2.1 below) cases, two examples of which are shown in Fig. 1(d) and 1(e). There are many images of clinical interest that are neither healthy nor unhealthy according to our criteria, and we exclude these from our computational analysis in this paper.

Fig. 1

Extracting gland and inter-gland lines. (a) Lines along the centers of the glands can be located by computing the local maxima. The maxima are separated into the shown colored groups using a clustering algorithm. (b) Pixel intensity profile (after Gaussian smoothing) for one row of the image shown in (a). Gland centers are located at the local maxima. (c) Steps to process a cluster of local maxima points into a continuous curve. (d) and (e) The gland (red) and inter-glands (green) lines obtained for a healthy and an unhealthy image.

As shown in Fig. 1(d) and 1(e), glands form a zebra-strip pattern in a healthy eye, whereas this pattern is absent in an unhealthy eye. The approach is to detect the lines along the center of and between the glands (which will be called the gland and inter-gland lines), and also the width of the glands, and then define features based on them to train a classifier to differentiate the images.

2.

Methods

2.2.

Detection of Gland Length

Meibography images present difficult challenges for image processing. Apart from low contrast, nonuniformly illuminated, and out-of-focus images, frequently there are artifacts such as specular reflections and intruding eyelashes that can interfere with the detection of the glands. In this paper, a pre-processing step was taken to first manually edit away the artifacts. Then, the following computational methods are applied to detect the gland lines and widths.

2.2.1.

Normalization of nonuniform illumination

Depending on the direction and focus of the infrared light source, meibography images may be nonuniformly illuminated. Figure 2(a) and 2(b) shows examples where these artifacts can be observed. If the raw image Fig. 2(a) is just enhanced using histogram equalization,¹³ the nonuniform illumination inherent in the raw image results in an image Fig. 2(b) that is bright at the center but dark at the edges. Although the stationary points (cf. Sec. 2.2.2 below) can be located in the bright region, this is difficult in the dark regions. Hence, a preprocessing step to normalize the nonuniform image intensity is needed so as to extract the stationary points in the darker regions as well. The result of normalization followed by enhancement is shown in Fig. 2(c). The stationary points obtained from Fig. 2(c) are added to those obtained from Fig. 2(b) (cf. the full algorithm in Fig. 3).

Fig. 2

(a) Raw meibography images have low contrast, making it difficult to visually recognize the meibomian gland structures. Specular reflections (red arrow) interfere with conventional enhancement methods. (b) Effects of nonuniform illumination on histogram equalization. The resultant image is brighter at the center [green arrow (1)] and darker at the sides [purple arrows (2) and (3)], making it difficult to detect glands at different locations. (c) Enhanced image after compensating for nonuniform illumination.

Fig. 3

Flowchart outlining steps for extracting gland and inter-gland lines. Reference for specific algorithm: FIFO, dilation,¹³ floodfill,¹³ skeletonization,¹⁴ and pruning.¹⁵

To normalize the image intensity, the nonuniform light illumination is modeled as

Eq. (1)

$$I(x)=b(x)J(x)+n(x),$$

where $I(x)$ is the observed image intensity at pixel $x$, $J(x)$ is the original signal to be restored, $b(x)$ is the unknown nonuniform illumination field, and $n(x)$ is the noise. The goal is to recover the original signal $J(x)$ by estimating $b(x)$. $b(x)$ is estimated using fuzzy c-means clustering method¹⁶ and the original signal is estimated as $J(x)=I(x)/b(x)$. After getting $J(x)$, a 2D contrast enhancement method¹⁷ is applied to improve the contrast. The result is shown in Fig. 2(c).

2.2.2.

Extracting gland and inter-gland lines

To detect the gland and inter-gland lines, the enhanced images [e.g., Fig. 2(b) and 2(c)] are smoothed using a Gaussian kernel ($45\times 45$ window, ${\sigma }_x={\sigma }_y=7.25$).^* Consider the intensity profile of the smoothed image in the horizontal direction [Fig. 1(b)]. The gland centers and inter-gland points are located at the local maxima and minima of this profile, where the gradient of the pixel intensity vanishes. The colored pixels of Fig. 1(a) show that the maxima points lie along the centers of the glands. To separate the maxima points into different groups [as shown in Fig. 1(a), where pixels belonging to the same group are represented using the same color], two maximum points are considered as belonging to the same group if they are separated by a Euclidean distance of lesser than 10 pixels.^† After grouping the points, to transform each group of points into one continuous line, the morphological processing steps in Fig. 1(c) are applied. The idea is to first dilate all the pixels so as to merge them into one connected component [Fig. 1(c) ii, merging], fill up any inner holes [Fig. 1(c) iii, filling], thin it to one pixel thick [Fig. 1(c) iv, thinning], and finally remove side branches [Fig. 1(c) v, pruning]. The full algorithm is summarized in Fig. 3.

Figure 1(d) and 1(e) shows the gland (maxima, red) and inter-gland (minima, green) lines obtained for a sample healthy and unhealthy image. By visual inspection, it is observed that there is a tendency for the lines in the healthy image to be longer than those in the unhealthy one. The arclength of a gland line is used as an approximation for the length of the corresponding gland, and the average arclength of all the lines (both gland and inter-gland lines) in an image is used as a feature of that image for classification.^‡

2.3.

Detection of Gland Width Using SIFT-Shannon Entropy

Scale invariant feature transform (SIFT)¹⁸ is an algorithm to detect and describe local features in images. For meibography images, the scale feature in SIFT is used to represent the thickness or width of the glands. Figure 4(a) shows an example of the so-called SIFT key points (red circles) computed by applying SIFT on a histogram equalized image (no normalization). The size of the key point (its scale) detects the width of the gland and the inter-gland distance [blue box in Fig. 4(a)].

Fig. 4

SIFT key points (red circles) of a healthy and an unhealthy image. Radius indicates the scale of the key point. Blue box: The scales detect the gland widths and inter-gland distances.

Figure 4(a) and 4(b) compares the key points of a healthy and an unhealthy image. For an unhealthy image, the scales are generally smaller than those of a healthy image. The average width of all the glands in an image is approximated by the average of the scales of the key points. Let $s^{\sigma }$ be the scale of the key point $\sigma$ on the image. Then the average scale $\overline{s}$ is defined as

Eq. (2)

$$\overline{s}=\frac{1}{N}\sum _{\sigma }{s}^{\sigma },$$

where $N$ is the total number of SIFT key points and $\sum _{\sigma }$ means sum over all the key points.

It was also observed that for the healthy glands, because of the zebra-strip gland patterns, neighboring key points have the same scale, and are located in an orderly manner along the center of and between the glands [blue box in Fig. 4(a)]. For an unhealthy image, because there are no gland patterns, the key points are randomly scattered and have nonuniform scales. This difference in local scale distribution can be captured using Shannon entropy, which is a measure of the uniformity of a distribution. Let $s_i^{\sigma }$ be the scale of the $i$’th nearest key point to the key point $\sigma$. Define the normalized scale $p_i^{\sigma }$ for the $i$th nearest key point as

Eq. (3)

$$p_i^{\sigma }=\frac{{(s_i^{\sigma })}^2}{\sum _{i=1}^n{(s_i^{\sigma })}^2},$$

where $n$ is the number of nearest neighbor SIFT points considered (we use $n=20$). The Shannon entropy of the key point $\sigma$, $S(\sigma )$, is defined as,

Eq. (4)

$$S(\sigma )=-\sum _{i=1}^n{p}_i^{\sigma }\ln p_i^{\sigma }.$$

$S(\sigma )$ is maximized only when all the $p_i^{\sigma }$ are equal (i.e., uniform). The average entropy of the image $\overline{S}$ is

Eq. (5)

$$\overline{S}=\frac{1}{N}\sum _{\sigma }S(\sigma ).$$

$\overline{S}$ should be high for a healthy image because of its uniform local distribution of scales, and low for an unhealthy one because its scales are randomly and nonuniformly distributed.

The average arclength, $\overline{s}$, and $\overline{S}$ will be used in the next section as features for classification.

3.

Results

3.1.

Separation Healthy and Unhealthy Images

Figure 5(a) shows how the healthy and unhealthy images are distributed according to their average arclength. The $x$-axis is the average arclength of the lines in an image and, to aid visualization, the $y$-axis shows the coefficient of variation of the lines [i.e., (standard deviation)/(mean)]. Each image is represented as a point. The healthy images (blue points) are well-separated from the unhealthy images (red points) along the $x$ direction, meaning that the average arclength is a good feature for separating the healthy and unhealthy class. Figure 5(b) shows the distribution of the healthy and unhealthy images according to their $\overline{s}$ ($y$-axis, average scale) and $\overline{S}$ ($x$-axis, average entropy). Once again, the healthy and unhealthy points are well-separated, meaning that they are also good features.

Fig. 5

Separation of healthy and unhealthy images in feature space. Each point represents an image.

4.

Discussion

In this paper, meibography images are classified using criteria based on gland lines and width, which are different from the well-established criteria based on the area of meibomian gland loss.⁷ In order to compute the area of meibomian gland loss, the gland region needs to be segmented from the background nongland region. Without any input from the user, this is computationally challenging because the image pixel intensity changes gradually between the borders and edges of the glands, making it difficult for the algorithm to decide the appropriate bounderies between glands. Gland lines and width, on the other hand, are easier to detect without the user, and the results presented in this paper can be considered as successful first steps based on these simpler features. In future work, these two features will be combined to segment the individual glands. This will allow the gland area to be computed, which will then allow images to be graded according to the established criteria.

This study focused on images of the upper eyelids because it was easier to obtain a uniformly focused image of the tarsal plate. Clinically, however, it is sometimes easier to work with the lower lids because they are less uncomfortable for the patients. The methods presented here have also been tested on 10 healthy lower lid images. It was found that the gland and inter-gland detection algorithm is effective for the lower lid images. For the detection of gland widths, the bright and dark strips of the zebra-strip pattern are different from the upper lids, with the dark strips being very narrow and the bright strips being much wider and frequently touching. Width detection for the bright and dark strips should be done separately, and using different scales. The calculated features for the lower lids are also different compared to the upper lids, so classification of upper and lower lids should be performed separately. Overall, the techniques here are applicable to both upper and lower lids.

In addition to gland length and width, meibomian glands can also be characterized by other morphological features such as shape, contour, and tortuosity.^5–7,9 Such features are important when grading images that are intermediate between the healthy and unhealthy ones. A limitation of the current study is that these features were not assessed. The healthy and unhealthy glands images considered here are extreme cases, and, although easy for any trained examiner, is challenging for an automated algorithm. This paper focuses on them because any algorithm working without user input should first be able to classify these simplest cases before addressing more difficult intermediate cases. It is shown here that based on the detected gland lines and width, the classification performed by the algorithm agrees closely with that done by human experts.

When applied to the classification of intermediate images, the current algorithm utilizing only gland length and width failed to produce satisfactory result. Instead of being distinctly separated in feature space like the healthy and unhealthy cases (Fig. 5), the intermediate cases overlaps with the healthy and unhealthy points. This makes it difficult for the SVM to classify. However, this is expected because studies have already shown that the meibomian gland structure of intermediate cases must be described by additional morphological features such as tortuosity,^5–7,9 whereas in the current computational method, only the length and width of the gland are considered. To treat the intermediate cases, additional morphological features can be computed based on the gland lines already obtained, and be used to expand the dimension of the feature for the SVM. The assessment of these morphological features and their application to the grading of intermediate images will be reported in future work.

An important difficulty encountered in this work concerns how images are taken. Slightly different ways of taking an image will not influence the gland structure. But because it introduces complication during the computational analysis, the final conclusions might be affected. This problem is important in this paper because the objective is to have a user-free assessment of the gland structure, so any error caused by detection of errorneous artifacts present in the image will not have the chance to be corrected by the user. This ultimately introduces errors into the final results. To be specific, there may be artifacts such as intruding eyelashes, specular reflections from the tear film, and misalignment of the gland region in the image. In this paper, the primary concern is to demonstrate the ability of the proposed algorithms in detecting meibomian glands, and so the computational difficulty of handling artifacts was dealt with heuristically. An image editing software was used to edit away parts of the image not directly related to the gland region by manually drawing out the region of interest in an image. The gland detection algorithm are then applied to the identified region. Although this process of drawing out gland region manually may appear similar to what has already been done in previous studies,^5–7,9 the important difference is that here the region inside the identified area is subjected to gland detection by the algorithm, whereas in previous studies the area of interest are drawn by the user with the intention of locating the glands.

As the ultimate objective is to achieve completely automated detection of gland regions without any user input, it is important to also develop an algorithm for locating the area of interest. It is found that this can be achieved by a simple change in the imaging protocol. The magnification when taking an image should be lowered such that the upper eyelid margin and the edge of the upper tarsal plate are both visible on the image; the margin and edge can then be used as references lines to define the gland region. This two lines can be detected using image processing techniques,¹³ and then the area of interest is obtained. Such a user-free way of locating the area of interest is desirable because it eliminates the subjectiveness introduced by clinicians when drawing the area of interest. Different examiners may draw the area differently; an algorithmic method, on the other hand, will always produce the same region everytime. In the current study, images analyzed were collected from an earlier clinical study where the need to include the upper eyelid margin and upper tarsal plate in the images was not foreseen. In subsequent studies, this precaution will be taken and the area of interest will be automatically identified. This work will be reported in a future paper.

5.

Conclusion

In this paper, a use-free computational approach to detecting meibomian glands and classifying meibography images was presented. The gland and inter-gland lines, and gland width were detected algorithmically, and then used as features for classification of the images. The classification results by the computational approach agrees closely with that done by human experts, with a specificity of 96% for healthy images, and sensitivity of 98% for unhealthy ones.