Molecular Vision 2018; 24:767-777
<http://www.molvis.org/molvis/v24/767>
Received 31 July 2018 | Accepted 30 November 2018 | Published 01 December 2018
Molecular Vision 2018; 24:767-777
<http://www.molvis.org/molvis/v24/767>
Received 31 July 2018 |
Accepted 30 November 2018 |
Published 01 December 2018
Michael A. Simmons,1 Alexander V. Cheng,2 Silke Becker,2 Richard D. Gerkin,3 M. Elizabeth Hartnett2
1Department of Ophthalmology, University of Texas Southwestern, Medical Center, Dallas, TX; 2John A Moran Eye Center, University of Utah, Salt Lake City, UT; 3Department of Internal Medicine, University of Arizona College of Medicine - Phoenix, Phoenix, AZ
Correspondence to: M. Elizabeth Hartnett, John A Moran Eye Center University of Utah, Salt Lake City, UT; Phone: (801) 587-7173;
FAX: (801) 581-3357; email:
Purpose: The aim of this study was to create an algorithm to automate, accelerate, and standardize the process of avascular area segmentation in images from a rat oxygen-induced retinopathy (OIR) model.
Methods: Within 6 h of birth, full-term pups born to Sprague Dawley rat dams that had undergone partial bilateral uterine artery ligation at embryonic day 19.5 were placed into a controlled oxygen environment (Oxycycler, BioSpherix, Parish, NY) at 50% oxygen for 48 h, followed by cycling between 10% and 50% oxygen every 24 h until day 15. The pups were then moved into room air until day 18.5. Ten lectin-stained retinal flat mounts were imaged in montage fashion at 10x magnification. Three masked human reviewers measured two parameters, total retinal area and peripheral avascular area, for each image using the ImageJ freehand selection tool. The outputs of each read were measured as number of pixels. The gold standard value for each image was the mean of the three human reads. Interrater agreement for the measurement of total retinal area, avascular area, and percent avascular area was calculated using type A intraclass correlation coefficients (ICCs) with a two-way random effects model. Automated avascular area identification (A3ID) is a method written in ImageJ Macro that is intended for use in the Fiji (Fiji is Just ImageJ) image processing platform. The input for A3ID is a rat retinal image, and the output is the avascular area (in pixels). A3ID utilizes a random forest classifier with a connected-components algorithm and post-processing filters for size and shape. A separate algorithm calculates the total retinal area. We compared the output of both algorithms to gold standard measurements by calculating ICCs, performing linear regression, and determining the Dice coefficients for both algorithms. We also constructed a Bland–Altman plot for A3ID output.
Results: The ICC for percent peripheral avascular/total area between human readers was 0.995 (CI: 0.974–0.999), with p<0.001. The ICC between A3ID and the gold standard was calculated for three image parameters—avascular area: 0.974 (CI: 0.899–0.993), with p<0.001; total retinal area: 0.465 (CI: 0.0–0.851), with p=0.001; and the percent peripheral avascular/total area: 0.94 (CI: 0.326–0.989), with p<0.001. In the linear regression analysis, the slope for prediction of the gold standard percent peripheral avascular/total area from A3ID was 0.98, with R2=0.975. A3ID and the total retinal area algorithm achieve an average Dice coefficient of 0.891 and 0.952, respectively. The Bland–Altman analysis revealed a trend for computer underestimation of the peripheral avascular area in images with low peripheral avascular area and overestimation of peripheral avascular area in images with large peripheral avascular areas.
Conclusions: A3ID reliably predicts peripheral avascular area based on rat OIR retinal images. When the peripheral avascular area is particularly high or low, hand segmentation of images may be superior.
Retinopathy of prematurity (ROP) is a potentially blinding disorder that involves delayed and disordered retinal vascular development with subsequent abnormal vasoproliferation into the vitreous, ultimately leading to retinal detachment. It has been associated with high supplemental oxygen at birth and fluctuations in oxygenation delivered to the retina of preterm human infants [1]. It has been described by a two-stage hypothesis in which normal retinal vascular growth stalls or newly developed vessels are injured (Phase 1), followed by abnormal intravitreal vasoproliferation (Phase 2) [2].
Animal models play an essential role in ROP research particularly in studying cell-environment interactions and signaling events surrounding retinal vascular development and abnormal vasoproliferation that would otherwise be unsafe, unethical and logistically impractical in preterm human infants. Several animal models exist. All models involve exposing newborn, full-term animals to supplemental oxygen under various protocols to inflict types of oxygen-induced retinopathy (OIR), each type of which shares features with human ROP [2,3].
The rat 50/10 model is one common model of OIR. It involves exposing newborn rat pups to supplemental oxygen levels that fluctuate between 50% and 10% inspired oxygen every 24 h [4,5]. One key similarity between the rat OIR model and human ROP is the creation of a central retinal vascular zone with an avascular periphery. In contrast, other models (including the mouse OIR model) compromise already-developed physiologic vascularity (vaso-obliteration), which manifests as a central avascular area (see Figure 1) and may mimic what occurs after high oxygen use in some premature infants [6]. The ratio of peripheral avascular to total retinal area is a clinically relevant parameter [7] and a standard endpoint measured in OIR rat models [8,9].
For many years, the gold standard method for calculating peripheral avascular retinal area in rat OIR models has been hand-traced segmentation of the vascular perimeter [10] in digital images of flat-mounted, enucleated eyes. In this “manual” approach, researchers trace the total retinal area and total peripheral avascular area using a freehand selection tool in image analysis software, such as ImageJ, and report these areas in units of pixels. This is a time-consuming process that introduces inter-rater variability and complicates comparisons across experiments and between laboratories.
Recently, two groups have published algorithms to automate analysis of images from the mouse OIR model [11,12]. Both groups demonstrated that the performance of their algorithms was equivalent to human expert tracings. Their work represents a significant step forward toward the automation and standardization of retinal analysis in OIR research. However, because of the differences between mouse and rat retinal vascularization in OIR models, Xiao et al. note that their algorithm does not support retinal avascular area segmentation in rat OIR models [11]. Automating the identification of peripheral avascular area in rat retinas represents a distinct computational challenge with no readily available solution.
In this article, we present a novel, machine learning-based, and fully automated function, automated avascular area identifier (A3ID), to detect the retinal peripheral avascular area in retinal flat mounts from a modified OIR rat model. We also present a gold standard data set of retinal images with human annotations to validate the algorithm.
All animal experiments involved in the creation of the gold standard image set for this study were performed in accordance with the Guidelines and Regulations for the Care and Use of Laboratory Animals of the University of Utah in Salt Lake City and the Association for Research in Vision and Ophthalmology Statement for the Use of Animals in Ophthalmic and Vision Research. Experimental protocols were evaluated and approved by the Institutional Animal Care and Use Committee of the University of Utah. Animals had access to food and water ad libitum and were kept on a 12 h:12 h light-dark cycle.
Retinal images in this data set were obtained from rat pups subjected to an intrauterine growth retardation adaptation of the previously described rat OIR model [4,13-15]. In this adaptation, pregnant dams underwent bilateral partial surgical ligation of the uterine arteries on embryonic day 19.5 to reduce fetal blood flow. Pups were born naturally about 2 d later. Within 6 h of birth, the pups and their dams were exposed to 50% O2 for 48 h, followed by cycling between 10% and 50% O2 every 24 h until postnatal day 15. The animals were subsequently kept in room air and euthanized on postnatal day 18.5, after which retinal flat mounts were prepared [15]. Flat mount processing involved enucleation, fixation in 4% paraformaldehyde for 1 h, and dissection of the retina after removal of the cornea, lens, and vitreous. Retinas were separated from the RPE, choroid, and sclera. Alexa-Fluor 568-labeled isolectin GS-IB4 from Griffonia simplicifolia (2.5 μg/ml, Invitrogen, Grand Island, NY) was used to stain the retinal vasculature. Retinas were flattened onto slides, and an inverted fluorescence microscope (Olympus IX81, Olympus Corp., Tokyo, Japan) was used to capture a series of images at 10x magnification. Individual images were joined together using Metamorph 7.0 software (Molecular Devices Inc., Sunnyvale, CA).
Ten complete retinal montage images were included in the gold standard image set. All images were from the supporting data from a previous publication [15]. Three separate reviewers traced the full retinal outline and the vascular perimeter in each image using the freehand selection tool in ImageJ. The peripheral avascular retinal area was then calculated by subtracting the vascular retinal area from the total retinal area. One reviewer analyzed the entire set two additional times in a blinded fashion to permit calculation of intra-rater reliability. The output for all reads was the retina and peripheral avascular area in pixels. The final gold standard measurement was derived from the average of the first read for each expert for each image. For one set of reads, the actual tracings were also saved as region of interest (ROI) files to allow for direct visualization of the measured areas.
A3ID is a method for identifying the peripheral avascular area in images of rat retinas in OIR models written in ImageJ Macro with 400 lines of code. It is designed to be used in the Fiji (Fiji Is Just ImageJ) image processing platform developed by ImageJ [16]. As a method, A3ID can also be called in any other integrated development environment or even from a command line interface. We provide a separate algorithm for determining the total retinal area. Both algorithms and user documentation are readily available upon asking the corresponding author.
The input for A3ID and our total retinal area algorithm is a single, closely cropped montage image of a lectin-stained retinal flat mount from the rat OIR model (Figure 2A). The algorithm first converts images to an 8-bit format. Each image is then processed with the contrast-limited adaptive histogram equalization algorithm to normalize local contrast. Because many laboratory computers lack sufficient processing power to run a complete image, we also include a step in which the algorithm divides each image into four equal quarters that are then run separately.
The first step in identifying the peripheral avascular area was to segment the images into three categories—vasculature, peripheral avascular retina, and background. To accomplish this, we used four sub-images (i.e., small selections of the gold standard data set) to train our classifier. Each selection contained representative sections of each category. The Trainable Weka Segmentation Tool in Fiji [17] was used to learn the differences between each category. This tool enables a user to apply several segmentation algorithms from the Weka data mining tool to a given image set using a variety of image features. The user must first manually identify portions of the image representative of the category that they wish segmented by tracing it on the image with a free-hand selection tool (Figure 3). The total number of attributes, or manual tracings, made across all three classes, was 43. These attributes included 2958 instances, or manually classified pixels.
After developing the library of sub-images and creating the tracings and instances for training, we sequentially identified different combinations of features and classifiers available on Weka, performing ad hoc analysis of their output to identify a combination that achieved an acceptable balance of both accuracy and speed. The final training feature set consisted of Gaussian blur, Sobel filter, difference of Gaussians, variance, and derivatives and used the fast random forest classifier. Figure 2C shows a section of a fully segmented image.
The segmentation step described above identifies all the avascular areas of the retina (green, Figure 2C). However, it does not permit differentiation between the avascular spaces that lie between blood vessels and the peripheral avascular retina. To segment the peripheral avascular retina, we used the IJBlob Shape Filter plugin in Fiji [18], which is a connected-components algorithm. This algorithm searches a given image for components with similar color that are connected and excludes these components, producing an image similar to the one in Figure 2D.
We employed two final quality-assurance steps to differentiate the peripheral avascular spaces from those within the vascular retina. First, we introduced an area threshold to identify the larger avascular areas (which were more likely to be peripheral avascular areas) and omit the smaller avascular spaces (which were more likely to be found within the vascular retina). Second, we used the Shape Filter plugin again to remove the more circular areas because we observed that circular spaces are more likely to be spaces within the vascular retina, whereas peripheral avascular spaces generally have a rectangular shape. The remaining area after these two steps is the final peripheral avascular area (Figure 2E).
We determined the total retinal area with a separate algorithm involving two steps. First, we selected a pixel intensity threshold to exclude the non-retinal background from the brighter retinal pixels. The threshold we used was a pixel intensity of three. Second, we applied an area filter to exclude areas with fewer than 5×10^6 pixels. The result was a single “silhouette” of the retina from which we calculated the retinal area in pixels.
All statistical analyses for this project were performed using IBM SPSS Statistics for Mac, Version 24.0 (IBM Corp., Armonk, NY). We calculated type A intraclass correlation coefficients (ICCs) with a two-way random effects model. We calculated the precision (positive-predictive value), recall (sensitivity), and accuracy by defining a true positive as any pixel for which the algorithm’s classification agreed with manual segmentation. False positives, false negatives, and true negatives were determined in like fashion. We report the accuracy using a Dice similarity coefficient (also called F-score or F1 measure), which is the harmonic mean of the precision and recall.
We evaluated the inter-rater and intra-rater agreement in the gold standard data set by calculating ICCs between separate graders (inter-rater correlation) and repeat measurements by the same grader (intra-rater correlation) for total retina area, total avascular area, and percent avascular area. Table 1 displays the results of this analysis.
We evaluated the output of A3ID and the total retinal area algorithm in three ways. First, we calculated type A ICCs similar to our evaluation of the gold standard data set and performed linear regression analysis. The gold standard value in this comparison was the average of the first three expert measurements. The results of this analysis are displayed in Table 2.
It is possible that area measurements may be highly correlated but still inaccurate if an algorithm consistently misclassifies a proportionate amount of each image. To evaluate this possibility, we analyzed both algorithms’ classifications of each pixel in the gold standard data set and calculated precision and recall. We report the results of this evaluation along with the overall Dice similarity coefficients in Figure 4.
As a final evaluation step, we performed a post hoc direct visual comparison between algorithm output and manual tracings. Representative images of the output of A3ID are shown in Figure 5 (peripheral avascular area) and Figure 6 (retinal area). To help us analyze trends in algorithm output, we constructed a Bland–Altman plot for A3ID, which depicts differences between expert measurements and algorithm output across the range of total peripheral avascular area present in the gold standard dataset (Figure 7).
A3ID is a machine learning-based and fully automated tool for segmentation of retinal images from rat OIR models. It replaces the time-consuming process of manually segmenting retinal images. When combined with our algorithm for total retinal area, A3ID enables automated evaluation of the percent avascular area for a given image. The use of a single automated tool such as A3ID will benefit the OIR research community by increasing efficiency and standardizing measurements between individuals and across laboratories.
A3ID achieves an excellent average intraclass correlation of 0.974 when compared with human measurements from a gold standard data set. Likewise, the combination of A3ID and our algorithm for total retinal area achieves an intraclass correlation of 0.940. The segmentations themselves appear more intricate than manual segmentations. Furthermore, A3ID’s average Dice coefficient of 0.891 for peripheral avascular area compares favorably with the median Dice coefficient of 0.870 reported by Xiao et al. in their work on mouse retina avascular area identification [11]. Our algorithm for total retinal area maximizes recall by including a “silhouette” of the total retina. This includes the ciliary body and frequently incorporates some artifacts as well. Nevertheless, even with these sources of error it achieves a Dice coefficient of 0.952.
One limitation affecting the performance of A3ID is its variable handling of the ciliary body, which may be present in some images and not in others. A3ID inconsistently includes and excludes the ciliary body from its calculations of total peripheral avascular area. Erroneous inclusion of the ciliary body may lead to false increases in either the total retinal area or the total peripheral avascular area or both. Another source of error is that the protocol for preparing the retinal flat mounts, which involves cutting the spherical retina into a flat cloverleaf, creates artifact avascular spaces. These artifacts falsely elevate the total avascular and peripheral avascular area measurements. In the future, A3ID could be improved by including a step to directly segment and remove the ciliary body or leaflet-associated artifacts when present.
Small, isolated peripheral avascular spaces constitute an additional source of error. These small “island” segments of peripheral avascular area are isolated from the main peripheral avascular spaces by tears in leaflets, debris within the avascular area, and blood vessels growing very close to the edge of the retina. Although these islands of peripheral avascular retina should be included in the measured total peripheral avascular area, they are occasionally incorrectly excluded by A3IDs post-processing filters. Another source of error in the performance of A3ID stems from our choice to crop each image into four leaflets during processing. This step was necessary to improve the speed of the algorithm, but it introduces error because the contour of the retina is disrupted, which occasionally produces artifacts at the new boundaries between leaflets.
One final and important limitation is that A3ID was developed using images from a modification of the rat OIR model involving intrauterine growth restriction (IUGR) via ligation of the uterine arteries in the dam. Although the patterns of vascular abnormalities are consistent between IUGR images and other OIR images, we have not evaluated whether this may limit the performance or generalizability of the algorithm to other rat OIR models. Likewise, A3ID is intended only for rat OIR models and is not validated as a tool for peripheral avascular area detection in other animal models. A3ID does not reliably detect retinal avascular area in mouse OIR models (data not shown).
In this manuscript, we describe the use of A3ID in determining the peripheral avascular retinal area; however, our classifier initially segments all avascular spaces, including the small perivascular spaces that lie between the retinal vessels. The ratio of these perivascular spaces to the vascular retinal area is important and may reflect compromised physiologic vascularity, such as from high oxygen use or following treatment with anti-angiogenic agents, including intravitreal antibodies against vascular endothelial growth factor [19]. Although we include post-processing filters in A3ID to exclude the perivascular spaces, one future direction would be to rebuild the algorithm to quantify these spaces in addition to the peripheral avascular retinal area.
A3ID is a fully automated algorithm developed as a method for ImageJ that segments the total retinal area, total avascular area, total peripheral avascular area, and percent peripheral avascular area using an unmodified image of a retinal flat mount as its only input. A3ID reliably predicts peripheral avascular area based on rat OIR retinal images. At extremes of area measurements, hand segmentation of images may be superior. A3ID reliably predicts peripheral avascular area from rat OIR retinal images. When the peripheral avascular area is particularly high or low, hand-segmentation of images may be superior.A3ID and the gold standard data set used in its validation are assets available to those researching OIR upon asking the corresponding author.
We are grateful to Katy Brown and Randy Brown for their work on annotating retina images, Eric Kunz for his support in the laboratory, and to Zhiyong Lu, PhD for his review of the manuscript. Some data in this article were previously presented at the Association for Research in Vision and Ophthalmology Annual Meeting 2018. Images used in this study were analyzed and reported in Becker et al., cited below [15]. None of the authors have any commercial conflicts of interest, relevant to this work, to disclose. This work was supported by an unrestricted grant from Research to Prevent Blindness to the Department of Ophthalmology and Visual Sciences at the University of Utah, and by the following grants: R01EY017011, R01 EY015130, T35EY026511 (to MEH). Imaging support was provided by the National Eye Institute Vision Core grant EY014800.