Automatic Detection of Subretinal Fluid and Cyst in Retinal Images

Age-related macular degeneration is one of the most frequent reasons of acquired blindness in economically advanced countries. The constant growing of AMD patient population is more and more challenging. AMD means degeneration of the macula which is the region of the retina responsible for central vision. Since AMD affects only this specific part of the retina, unattended patients lose their fine shape- and face recognition, reading ability, and central vision [9].

Basically, AMD has two forms: dry and wet form, and the latter causes rapid and serious visual impairment in 10% of the cases [13]. In this type of the disease, abnormal angiogenesis starts from the choroid under the macula. Fluid and blood leak out of the neovascularized membrane into retina layers that ruins the photoreceptors.

Experiments have demonstrated that the vascular endothelial growth factor (VEGF) plays a vital role in the formation of choroidal neovascularization [4]. Currently, the most common and effective clinical treatment for wet AMD is anti-VEGF therapy, which is a periodic intravitreal (into the eye) injection [11].

In the last decade, Optical Coherence Tomography (OCT) has been widely used in the diagnosis of AMD and follow-up therapy. Spectral domain OCT (SD-OCT) produces 3D volumes of data, which have been useful in clinical practice. Existing OCT systems are partially suited to monitoring the progress of the disease, but OCT shows many features about AMD such as hyper-reflective dots (HRD), subretinal fluid and cysts. Figure 1 illustrates an SD-OCT B-scan with biomarkers of AMD.

Original Optical Coherence Tomography (SD-OCT) image with biomarkers of AMD in inverted display (as our medical experts use it in daily routine).

Full size image

A large number of publications in the scientific literature deal with the problem of detecting retinal layers based on various techniques. One approach is the automatic segmentation procedure using graph theory [1, 2, 6]. In this approach, the graph nodes usually relate to image pixels, the graph edges are assigned to pairs of pixels, the edge weights depend on the intensity differences between the node pixels, and also may depend on the spatial distance between the pixels. Image segmentation then becomes a graph cutting problem, which can be solved via dynamic programming. These approaches are less tolerant to noise, that is a disadvantage, because real images are often very noisy. Another idea relies on the well-known energy-minimizing active contour method which, unfortunately, also has problems in handling low contrast and noise. Yazdanpanah et al. [18] suggested a multi-phase framework with a circular shape prior in order to model the boundaries of retinal layers and estimate the shape parameters. They used a contextual scheme to balance the weight of different terms in the energy functional. Machine learning is widely used in recent years, also for retinal image analysis. Lang et al. [12] used random forest classifier to segment retinal layers. The random forest classifier learns the boundary pixels between layers and produces an accurate probability map for each boundary, which is further processed to finalize boundaries. Procedures based on active contour or machine learning provide effective solution, but these methods are too time-consuming. Hassan et al. [8] used a structure tensor approach combined with a nonlinear diffusion process for layer detection. A structure tensor is a second-moment matrix that shows similarities and prominent orientations of the image gradient. Some other approaches use optimized boundary tracking [5] or polynomial smoothing [14]. These algorithms are rather complex. We have developed an algorithm which uses simple operations to localize subretinal areas. It is based on vertical profile analysis [10].

Relatively few publications deal with the problem of automatically detecting cysts. Gonzalez et al. [7] described a method based on watershed segmentation and different machine learning classifiers. They focused on feature extraction which can help to eliminate false regions. Other approaches can also be found in the literature, two of which are discussed in Sect. 2.1.

In this study, we deal with the automatic localization of subretinal fluid areas and cysts and also analyze major retinal layers, since layer information can help localizing and distinguishing fluid and cyst regions. We present an algorithm that automatically delineates the ILM (inner-limiting membrane) and RPE (retinal pigment epithelium) retinal layers. We also describe a method to detect subretinal fluid and cyst segments and distinguish them from each other. We compare our results with some other algorithms from literature.

In this section, we present several algorithms that use different approaches to determine cyst and subretinal fluid. First, we briefly describe two existing approaches that we re-implemented according to the original papers for comparison. Then, we describe our novel approach in more detail. The procedures first delineate the inner and outer boundary layers (ILM and RPE, resp.) for easier determination of the important areas.

2.1 Literature Procedures

Firstly, we describe Wieclawek’s [16] algorithm. First, the input image is normalized to the [0 1] interval, because images can be made with different settings, so that their intensity range may vary. The OCT images are affected by distortions like noise, so the authors used a non-linear filtering to reduce this effect. Next, they applied a spatial averaging filtering technique which is based on the real product of complex diffusion. The tools of mathematical morphology was used to delineate specific cystic areas, based on the observation that cysts appear as darker segments in the images. Among other operations, they used H-minima transform to highlight important regions. The single control parameter is a threshold value. This value has been fixed experimentally to 30% of the maximum brightness in the image. The next step was the binarization of the obtained image with a given threshold value. Since the result may still contain false regions, they filtered all objects that are above ILM and below RPE layers.

Wilkins et al. [17] investigated the problem from another point of view. They discarded color information in the first step of the algorithm and determined the major layer boundaries. To improve image quality and filter noise, they used the combination of a median filter and a bilateral filter during preprocessing. After binarization, the method determined the boundaries of the remaining possible cyst segments and they defined three conditions based on empirical studies. They investigated the extent of the objects, the degree of scattering between the intensity of the pixels in the segment, and whether the object is located between the ILM and RPE layers.

Sample OCT image before (a) and after (b) applying the function \(\kappa ^{*}_{\nu }\). It filters out the noise outside the retina and highlights boundary layers, so it is easier to delineate them.

Full size image

2.2 Proposed Method

OCT images are affected by distortions like “shadowing” by blood vessels, that may yield to false detections. In the first step of our proposed algorithm, we improve the image quality by noise filtering and contrast enhancement using a fuzzy operator [3]. This operation can highlight major retinal layers. We analyze vertical profiles of the filtered image and large intensity steps in pixel density are assumed to correspond to the change of tissue. The used fuzzy function is defined as

where \(\nu \) is a threshold, x is pixel intensity and \(\lambda \) denotes the sharpness of the filtering. The function \(\kappa ^{*}_{\nu }\) can highlight boundary layers while suppressing noise. We determine dynamically the input parameter \(\nu \) in a simple way. We sample from the top range of the image and calculated average intensity for this ROI. The parameter \(\lambda \) was set to 3 empirically. Figure 2 shows an example of applying the function \(\kappa ^{*}_{\nu }\).

After filtering, we divide the image into bars with fixed width. A bar consists of 10 contiguous pixel columns and we calculate horizontal projections of each bar to determine boundaries. One of the main steps of our proposed method is to analyze the vertical profiles. This signal is usually noisy, so there is a need for filtering the data. We use the Savitzky-Golay filter [15] which is a smoothing digital filter. This filter is effective in preserving the relevant high frequency components of the signal, which is an important aspect for our detection method.

Determining the outer layer boundary is harder than that of the inner boundary, because Choriocapillaris and Chorodoidal vessels are located under the RPE layer. The intensity of these regions vary, so several peaks appear in the projections. Fortunately, in most cases, these minimum points are not prominent, and do not cause problem in choosing the right locations. The algorithm chooses the most important local minimum from the projected data to identify the possible inner and outer layer. In the next step, we filter out the outliers and we fit a curve to the remaining points.

After we determined the boundary layers, the next step is the segmentation of fluid and cyst areas. It can be observed that these regions appear as spots with brighter intensity in the image. For processing, we use the inverse of the signals, because our medical colleagues used the inverted presentation of images for visual assessment and also exported the image data for us in this format. The zones of the disease and the intensity of the vitreous body of the eye are almost within the same range (if distortions are not considered). Anisotropic diffusion is used to eliminate various errors from imaging or blood shadows. Using the filter, it is more apparent that some parts of the retina are within a given intensity range, so we quantize the grayscale image into five intensity levels. Our observations showed that the layers of the retina are only in some intensity ranges, so this operation facilitates the separation of the 8–10 main retinal layers. During the binarization, we keep the brightest points, because we know that the reflectivity of the cyst similar to that of the vitreous body. By this step, we create a mask for the active contour process. To achieve the appropriate segmentation result, the input parameters of the model were given based on empirical studies. After that, there may be holes in some objects, so we use hole filling. Figure 3 illustrates these steps.

Intermediate stages of the proposed algorithm. (a) Anisotropic diffusion, (b) Quantization, (c) Binarization, (d) Result of active contour.

Full size image

So far we have identified possible important segments. Next, we need to separate cyst, fluid and false segments from each other. This step may be omitted in some cases, when there is no object in the picture that would be detected as a pathological mutation. We developed a condition set for filtering cysts and fluid regions. We test the fulfillment of the four criteria for classification at the object level.

The greater distortion of layers was observed in those parts where these symptoms appear, so firstly, we examine whether the actual layer is distorted or not. This plays an important role in distinguishing between cysts and fluid areas. To determine whether the layers are creased or not, we use the top layer boundary points. We calculate the minimal y coordinate of the top points. For this, the method does not take into account the left and right 25% of the image. The sides of the image may not contain information, because of the image registration, and large distortions also may appear in these parts of the images. In Fig. 4 we illustrate these mentioned effects. Next, we search the minimal y coordinate of the top points, divide this image into two parts along the established peak, and we investigate the maximal y point on both sides. We estimate the degree of creasing. Sometimes there is no change in the middle of the image, so we determined a threshold to decide if there is any crease in the slice or not. The threshold for the minimum y point and the given maximum y point was defined in 5 pixels experimentally.

Example of distortion in consecutive slices.

Full size image

Various a priori information can be used to distinguish between cysts and fluid areas, and to filter out the false segments. We investigated a condition system and we considered where the object is located within the retina, what is the extent and the shape of the object, and whether the layer is distorted or not. Fluid areas have larger extent and they are located close to the bottom layer boundary or may appear in the distorted area. In the case of higher distortion of the retina, they may also appear on the left or right side. When determining fluid regions, it is also important to examine the cases where there is no creasing in the layer. In the case of cysts, we need to find objects with oval shapes and the observations show that these segments are in the increased zone. Contrary to fluid regions, cysts are found in higher layers. We distinguish the symptoms based on these characteristics. False segments may also appear in the image, but these are small objects, so we can remove them easily with an area-based filtering. The key stages of the procedure are summarized in Fig. 5.

Flowchart of the proposed algorithm.

Full size image

3.1 Image Data

Our evaluation dataset contained 11 Heidelberg Spectralis OCT scans of wet age related macular degeneration patients treated with anti-VEGF intravitreal injections. The scanning parameters are: 49 scan pattern, pattern size: 5.8 \(\times \) 5.8 mm, distance between B-scans: 121 \(\mu \)m, size X: 512 pixel, size Z: 496 pixel, the pixel size is 11.44\(\,\mu \)m and 3.87 \(\mu \)mm in X and Z direction, respectively.

Manual ILM and RPE layer segmentation was performed by ophthalmologists for 7 image sequences. This was considered as a ground truth for evaluating the boundary layer detection method on these 7 image volumes.

3.2 Results and Discussion

We implemented our proposed method in MATLAB, using the Image Processing Toolbox. We evaluated our retinal layer detection algorithm in two different ways. We compared the results of our algorithm against the manual delineations and we also compared the proposed method for automatic detection of subretinal fluid and cyst with some published methods from the literature.

Firstly, we consider the result of localization of major layers. We calculated the mean, maximum and standard deviation of boundary errors for every surface. The 7 curves shown in Fig. 6 depict the error histogram for those OCT volumes where manual annotation was available. Each curve aggregates the boundary errors in the 49 scans (slices) of a study. It can be observed that the largest error is between 1 and 4 pixels in most cases and Table 1 asserts to this statement. As presented in Table 1, the maximal distance between manually segmented and automatically detected layer boundary is 19 pixels (ca. 73.5 \(\mu \)m). This deflection comes from two sources, namely, the substantial jumping between B-scans and layer distortions due to the disease. Unfortunately, we could not exploit 3D information directly to segment the retina layers, because there are some anomalies among slices of the OCT volume, due to the image acquisition and registration process (within the device’s software).

Error histogram of 7 image sequences.

Full size image

In Sect. 2, we have presented two methods from the literature for the segmentation of cysts, as well as our proposed method, which is also suitable for delineating fluid areas. Unfortunately, expert annotation was not yet available to evaluate segmentation results, so we compare visually the outputs of the algorithms.

Figure 7 illustrates segmentation results by the algorithms in some slices. The method developed by Wieclawek detected fewer possible cyst regions, which may be due to the fact that the given threshold only keeps the actual light points. The disadvantage of this is, that important areas may be lost during processing. The other method from literature by Wilkins yields almost the same segmentation results, but in many cases it holds false objects, because the thresholds are not dynamically defined. In contrary, our method uses dynamic requirements based on a priori information.

We tested the re-implemented earlier published methods on images in which cysts and fluids may also appear. They can also detect these regions because these segments are also lighter object in the image. Our algorithm, however, as it can be seen in Fig. 7(c), can also distinguish these two types of structures from each other.

Illustration of detected cyst (red curve) and subretinal fluid (blue curve) regions by the described algorithms. Columns: (a) Wilkins et al., (b) Wieclawek, (c) Proposed method. Top three rows contain only cysts, while the bottom two scans also have fluid areas. (Color figure online)

Full size image

We presented a novel algorithm for the detection of subretinal fluid areas and cysts and we compared it with two methods from the literature for cyst localization. After having seen the results, the medical colleagues believe that digital image processing can help the quantitative assessment of the OCT features of AMD by providing automatic tools to detect abnormalities and to describe by objective metrics the current state and longitudinal changes during disease evolution and treatment. Using SD-OCT to follow up changes of subretinal fluid and cysts volume will become a useful tool in detecting subtle changes during the treatment process. Further studies are planned to evaluate these new tools in a cohort of AMD patients.