The algorithm uses the following steps for detecting vascular tortuosity: (1) detecting vessels on imaging, (2) eliciting the vessel frame with the help of a thinning method, (3) determining vessel crossovers and bifurcations, and (4) measurement of local and global tortuosity in fundal photographs.

The algorithm that we used for revealing vessels was based on the local Radon transform, previously described by Pourreza et al. [10]. In this algorithm, a green plate from the RGB image is first used, having the most prominent contrast with the background of the image. Then, the image is divided into several overlapping areas and converted to a Radon transform using a round mask on the image.

The local Radon transform algorithm is primarily used to determine the presence or absence of vessels; then, if a vessel is present under the coverage area, the Radon transform is used to determine the length of the vessel as a line. The vascular framework is obtained using the data. Fig. 1 shows the result of the algorithm, displaying the vessels. We extracted the vascular structure of each photo up to the third bifurcation of major vessels in at least two quadrants.

The next step to study vessel tortuosity is to extract the vessel meshwork from the images. To accomplish this, we employed the thinning algorithm used by Guo and Hall [11].

In this study, to evaluate vessel tortuosity, we first used some artificial pictures with sinusoidal waves resembling fundal images; then, after achieving acceptable results in our tortuosity analysis, the algorithm was used on true fundal images.

Areas of vessel crossovers and bifurcations were first diagnosed and excluded from the analysis, because it was not possible to evaluate vascular tortuosity in these regions. To achieve this goal, we used an algorithm that draws a circle on the vascular meshwork and counts the places where vessels are crossed, marking them as n (n = 3 and n = 4 were defined as bifurcations and crossovers, respectively in that region). These regions were excluded from the vascular skeleton extracted from fundal images. Fig. 2 shows the result of vascular meshwork extraction with regards to exclusion of bifurcation and crossover sites.

After these steps, vascular tortuosity was calculated by estimating the curvature. For a continuous curve [y = f(x)], the curvature, κ, was calculated using the following relation: , where y′ and y″ are first and second derivatives of y.

For a non-continuous curve, estimation of the curvature is more difficult. For this reason, the curvature of these curves is estimated through different algorithms. In this study, we used the algorithm proposed by Frette et al. [12]. In order to estimate the curvature at any point, a circle is drawn with the center being on that certain point with a radius of b. With the crossover of the circle and the curve, the area between the curve and the circle, the green area A, can be calculated as shown in Fig. 3.

Frette et al. showed that with drawing this circle, a good estimation of the curvature could be calculated by:

When b is constant, the curvature is proportional to A. To find a formula that is proportional to the curvature of a certain curve and also exhibits a curvature of 1 for a flat line, the following was proposed by Frette et al.: , where A_tot is the total area of the circle, A is the largest area where the curve and the circle cross over, and is the estimated curvature of the curve.

This formula is calculated for each point of the vascular meshwork. Finally, the following formula is used for calculating the total tortuosity of a distinct vessel: , where m is the number of points where curvature is calculated and k_i is the curvature at point i.

Therefore, tortuosity is always greater than 1 for a curve.

To estimate vessel diameter, we used an algorithm based on the distance transform presented by Maurer et al. [13]. In this algorithm, the smallest distance of each vessel point to non-vessel point in the fundal images is measured. For evaluating the algorithm, initially, artificial images with lines of known diameter were used. The artificial lines were chosen to have both different diameters and direction. The results of artificial line diameter measurement by the above algorithm had an accuracy of 98.5%. Next, true fundal images were then evaluated.

Two quantities, including mean vessel tortuosity and dilation, were extracted for each fundal image and a property vector was mapped for each fundal image.

In the next step, the software categorized these images based on their property vector into one of three groups: plus, non-plus, and pre-plus.

The software classified the images based on three classifiers: the support vector machine (SVM), the k-nearest neighbor, and the multilayer perceptron network (MLP), which obtained the highest accuracy rate in comparison with the consensus of three ROP experts.

In this software, the classifier matches the property vector of each image and categorizes it to either plus, nonplus, or pre-plus group. In other words, the purpose of the classifiers is to determine the relation between the property vector and the image classification based on the opinion of the three experts.

Initially, every classifier passes a learning phase. In the learning phase, the classifier makes a connection between the property vector of the fundal images and an obvious diagnosis from the three subgroups of ROP, allowing the software to classify images having the most contiguity with the three ROP expert opinions.

The k-nearest neighbor is a classic model of one-class classifiers. This classifier is most useful when just one class of data is present or when obtaining the other classes is not cost beneficial. The goal of this classifier is to separate one class from the others by describing it; for example, when just plus images are labeled without any recognition of non-plus or pre-plus labeling, this classifier can simply detect plus images from among all the images.

The SVM is a model of two-class classifiers. In the learning phase of the software, the boundary of the two classes is made by the test samples given to the software. The rest of the samples are then classified based on this boundary.

The MLP was the final classification method used in this study. Inspired by the human brain, the MLP is an artificial neuronal network consisting of different layers and zones. The function of the artificial neuronal network layers is to classify the data according to its layers. Normally, three layers of neurons are present in the network: (1) the input layer: this receives the raw data, which were the property vectors in our study; (2) the hidden layer; and (3) the output layer: this is the output of the classifier (fundal images classification: plus, pre-plus, and non-plus).

Selecting the appropriate architecture and learning algorithm are the two issues when working with an artificial network like MLP. Appropriate architecture includes selecting a suitable number of layers and neurons that results in precise classification of the data set.

For classification, some of the images are used for reducing the classifier and others for testing them. We used k-fold cross validation for selecting and testing the images. In this method, the data is divided into different segments (the number of segments is represented by k) and then classified in k stages. In each stage, one of the segments is used as the learning data while the rest are employed as testing data. Then, for each stage (k), the mean accuracy is measured.

For the first time, the Imperial College London recently used a semi-automated image analysis for the diagnosis and quantification of retinal vessels on RetCam images, which resulted in the introduction of an algorithm known as RISA.

This algorithm analyzes images based on geometric factors such as the first and second derivatives of the intensity of the image, maximum gradient, and principal curvature.

RISA finds the maxima in the scale space, and has an extra constraint on these maxima in that the gradient must be very low as a starter point for region growing. A multiple-pass region-growing procedure is used, which progressively segments the vessels using feature information from the 8 neighboring pixels [17].

RISA was first used for evaluating arterial tortuosity in ROP patients, after which arteriolar and venous curvature, diameter, and tortuosity were shown to increase in patients with plus disease [18,19,20].

Gelman et al. [19] studied arterial and venous tortuosity and dilation using RISA. They evaluated tortuosity using the integrated curvature (pixel / radius) and tortuosity index (unitless).

Integrated curvature constitutes all angles of the vessel, standardized by the vessel length. Tortuosity index is defined as the vessel arch divided by the chord of the vessel from the beginning to its end.

In comparison, the mean vascular diameter is measured by the total vessel area divided by the length of that vessel presented in pixels [19].

In 2007, Chiang et al. [16] showed that in comparison with 22 expert opinions, RISA had a sensitivity and specificity of 70% for analyzing arterial and venous tortuosity, and a sensitivity and specificity of 60% in analyzing the vascular diameter. When using a linear combination for vessel tortuosity and dilation, sensitivity and specificity increased to 94%.

ROPtool is another software program presented by Duke University that is used for the evaluation of vessels in ROP patients. In a comparison between this software and the opinion of six pediatric ophthalmologists, a sensitivity of up to 97% was reported for the categorization of ROP patients to plus, pre-plus, or non-plus disease for 185 images taken with a RetCam [21,22].

ROPtool is developed from a program that derives the tubular structures from three-dimensional images (like cerebral vascular images taken by magnetic resonance angiography) [23].

It is based on centerline extraction; when given a start point such as a mouse click, its associated intensity ridge can be found by mapping image intensity to height [24].

With the optic disc in center and a radius of 30°, Wallace [21] used ROPtool to derive a smooth curve for each vessel, and calculated tortuosity as the true length of the vessel divided by the length of the smooth curve. In the newer version of this software, the distance between the optic disc and the fovea is used as a normalizing criterion for images with different magnifications. Additionally, both tortuosity and dilation were used to evaluate patients with plus disease. A sensitivity and specificity of 91% and 86% for vessel tortuosity and 78% and 84% for vessel dilation respectively were reported by Kiely et al. [25] in 2010 using ROPtool.

Another software program, computer-aided image analysis of the retina (CAIAR), which is an extension of the RISA algorithm, was introduced with a more precise and sensitive image segmentation compared to RISA. CAIAR extracts vessels from images using a scale-space framework with the highest probability for model fitting. This software uses filters sensitive to ridge-like structures in four different scales and then matches the best fitting of the vessel with the software output in order to extract the vascular structure. After that, the vessel diameter, length, and direction are analyzed using a Gaussian profile. With a second derivative of the Gaussian filter, measurements are checked to see whether the greatest contrast is in the middle of the vessel [26].

Shah et al. [27] used RetCam images of infants to validate and compare the output of CAIAR with that of clinician graders.

In another study using the Nidek NM200D noncontact camera (with a smaller field-of-view and higher-resolution images) and a vessel-fragmentation program called Vasculomatic ala Nicola ver. 1.1 (IVAN; Department of Ophthalmology and Visual Science, University of Wisconsin-Madison, Madison, WI, USA), the risk factors for developing plus disease from pre-plus disease were evaluated using CAIAR. Of 30 eyes, 11 had progression to plus disease. The mean width and tortuosity values of the three widest or most tortuous vessels predicted which eyes would require treatment [27,28].

With the review of programs used for extracting and analyzing the vascular architecture with different methods and formulations used for tortuosity and dilation, and by comparing them against the opinions of ROP experts, differences in sensitivity and specificity were elucidated. A comparison of the widely-used programs in the ROP field was done by Wittenberg et al. [29] in 2012.

Standardization of similar programs for labeling images as plus or non-plus is always a challenge. In a study performed in 2007 for image-based diagnosis of plus disease, 22 experts had similar opinions in only 22% of the images (seven of 34 images). In comparison with others, each expert had a mean kappa of 0.19 to 0.66, and in comparison with a standard reference determined by the overall opinion of 22 experts, the sensitivity and specificity of diagnosing plus disease ranged from 0.31 to 1.0 and from 0.57 to 1, respectively for a fundal image [7].

A study by Wallace et al. [8] showed a 27% difference among three experts on ROP when classifying 181 images to three grades: plus, non-plus, and pre-plus. Other studies have confirmed this difference in the opinions of different ROP experts.

The kappa coefficient is usually used to show the agreement between experts; for example, the kappa coefficient for agreement among the opinions of pathologists on a diagnosis of residual esophageal carcinoma in tissue samples is 0.87 [8]. Because of the difference among experts, we used the consensus of three experts to educate and standardize our software.

In the current study, we separated and elucidated vessels based on the local Radon transform as a novel technique, excluding points of crossover and bifurcation. Therefore, the accuracy for vascular architecture detection was im proved theoretically. However, it would be useful to design a comparative study to confirm this theory.

Various methods have been presented for the evaluation of vascular tortuosity. These methods can be classified into four categories: (1) methods that are based on the arch and chord length, (2) methods that are based on the curvature, (3) methods that are based on the angular changes, and (4) methods that are based on transformations.

Previous versions of ROPtool evaluated tortuosity based on the arch and chord length of the vessels, which overestimates tortuosity, especially for vessels with a gradual change in the curvature, in comparison with the newer smooth vessel curvature method [30].

In order to evaluate the vessel tortuosity algorithm, we first used artificial sinus waves. Then, after the algorithm achieved acceptable results, real fundal images were evaluated and vessel tortuosity was calculated by a formulation based on the curvature estimation. The results were very similar to the data set analysis based on chord length method reported by Grisan et al. [31]. As there was no access to the image database of other studies using methods including the integrated curvature in RISA or the smooth vessel curvature in ROPtool, it would be useful to conduct further studies comparing different methods for analyzing the tortuosity index.

For estimating vascular diameter, a new algorithm with a higher accuracy proposed by Maurer et al. [13] was used. This novel method can estimate the vascular diameter with an acceptable accuracy, but further study is essential to compare it with other methods used for vascular diameter estimation, such as mean vascular diameter based on pixels.

Standardization of the images with different magnifications can be done based on the optic disc size or by the distance between the optic disc and the fovea, as in our study. Although standardization with the latter method does not have error induced by different optic nerve head sizes and can be more accurate, it has potential for error and bias in image analysis, which is a limitation of our study as in previous works.

The algorithm in this study automatically evaluated vessel tortuosity and dilation in images, both together and separately, for the diagnosis of plus disease. A great advantage is that this process is not operator-dependent, although some factors that are important for some experts for diagnosing plus disease including vascular congestion, vascular branching, or appearance of peripheral disease can have an impact on the accuracy rate. These factors were not evaluated in this study.

Previous studies have concluded that determining the type of vessels, either arterioles or venules, is not important because a significant correlation exists between treated and non-treated eyes based on the most tortuosed vessel, regardless of the vessel type [9,32]. Therefore, in this study, we did not separate arterioles from venules. Thus, the tortuosity of both arterioles and venules were included in the overall tortuosity measure, as a previous study by Wallace [21] showed that it is not always possible to accurately differentiate between them, and it is common for both arterioles and venules to become tortuous as plus disease develops.

Another factor that affected the accuracy in this study was that the images were analyzed with centration on the optic nerve; therefore, peripheral retinal changes which are important to the diagnosis of plusplus disease were not included.

Keck et al. [33] used the arch length and chord ratio for evaluating vessel tortuosity. In their study, 22 experts evaluated 34 fundal images. They did not use vessel diameter because of different image magnifications. They found that arterial tortuosity had more significant changes when moving from the optic disc towards the peripheral retina. Therefore, peripheral vascular changes can affect the categorization of the patients and should be included in future studies [33].

A study performed in 2010 showed that experts had a lower cut-off point for labeling a patient's fundal images as plusplus disease in comparison to standard images of plusplus disease. This issue could lower the sensitivity of the algorithm; therefore, unlike previous studies, we did not use standard plusplus images (like Cryo-ROP images) and used our data set photographs for plusplus disease categorization for both the retina experts and software instead [15].

The classifiers used in this algorithm can improve diagnosis theoretically. Although the accuracy rate with the MLP classifier was higher versus the two other classifiers, it was still lower than previously used software programs like ROPtool and RISA. Software design improvement with larger data sets and more expert opinions could enhance the accuracy of plus disease detection by this algorithm in future studies.

Because all the images that were analyzed in this study were well-focused and of good quality, it is not clear whether the software shows similar efficacy when lower quality images are used.

Another source of error was the different pressures used for the RetCam contact camera, which can affect tortuosity and vascular diameter. The mode and time of imaging, the imaging angle, and the use of mydriatics can also cause measurement bias [34].

In conclusion, the new automated algorithm used for the diagnosis and screening of patients with plus ROP in this pilot scheme had an acceptable accuracy. With more improvements, this algorithm can be widely used in the future, especially in centers without skilled ROP staff.