Automated quantification of choriocapillaris anatomical features in ultrahigh-speed optical coherence tomography angiograms

: In vivo visualization and quantiﬁcation of choriocapillaris vascular anatomy is a fundamental step in understanding the relation between choriocapillaris degradation and atrophic retinopathies, including geographic atrophy. We describe a process utilizing ultrahigh-speed swept-source optical coherence tomography and a custom-designed “local min-max normalized masking” algorithm to extract in vivo anatomical metrics of the choriocapillaris. We used a swept-source optical coherence tomography system with a 1.6 MHz A-scan rate to image healthy retinas. With the postprocessing algorithm, we reduced noise, optimized visibility of vasculature, and skeletonized the vasculature within the images. These skeletonizations were in 89% agreement with those made by skilled technicians and were, on average, completed in 18.6 s as compared to the 5.6 h technicians required. Anatomy within the processed images and skeletonizations was analyzed to identify average values (mean ± SD) of ﬂow void radius (9.8 ± 0.7 µm), ﬂow void area ( 749 ± 110 µm 2 ), vessel radius ( 5.0 ± 0.3 µm), branch-point to branch-point vessel length ( 26.8 ± 1.1 µm), and branches per branch-point ( 3.1 ± 0.1 ) . To exemplify the uses of this tool a retina with geographic atrophy was imaged and processed to reveal statistically signiﬁcant ( p < 0.05 ) increases in ﬂow void radii and decreases in vessel radii under atrophic lesions as compared to atrophy-free regions on the same retina. Our results demonstrate a new avenue for quantifying choriocapillaris anatomy and studying vasculature changes in atrophic retinopathies. processing should aim to quantify atrophic retinopathies comparatively rather than absolutely. This study demonstrates the combination of ultrahigh-speed OCTA imaging of the choriocapillaris with local min-max normalized masking, a novel algorithmic tool for processing and quantifying anatomy in en face images of the choriocapillaris. When angiograms of the choriocapillaris are processed with LMNM, a visually clear, easily quantiﬁed representation of the choriocapillaris results. To conﬁrm the reliability of LMNM, we compared choriocapillaris tracings produced by it to traces made by humans and found them to correspond well. Then, we demonstrated that LMNM processed images can be easily quantiﬁed to extract values for anatomical metrics of the choriocapillaris and that these metrics are signiﬁcantly diﬀerent under most lesions in images of a retina aﬀected by GA. With further reﬁnement the combination of ultrahigh-speed OCTA and LMNM image processing may yield improved visualization and quantiﬁcation of the morphological changes in the choriocapillaris associated with a variety of atrophic retinopathies.


Introduction
The choriocapillaris is a ∼10 µm thick, flattened, single-layer capillary plexus fed by the choroid which underlies the Bruch's membrane and facilitates metabolic-exchange with the outer retina [1][2][3]. Morphological changes within the choriocapillaris have been linked to degradation of the retinal pigment epithelium (RPE) and photoreceptors in a multitude of retinal diseases [4]. Age-related macular degeneration (AMD) is estimated to affect 6.5% of Americans over 40 years old and is responsible for over half of all vision loss amongst Caucasians Americans [5]. In both AMD and geographic atrophy (GA), a late-stage form of AMD, dropout of the choriocapillaris is observed locally to formation of lesions or drusen [4,[6][7][8][9][10]. In diabetic retinopathy and late-stage glaucoma, global capillary constriction and reduction of capillary density have been observed in the choriocapillaris [11][12][13]. Though choriocapillaris changes are associated with many retinopathies, their role in disease progression remains unclear. We believe that high resolution in vivo choriocapillaris imaging will permit detailed quantification of subtle disease-related morphological changes, potentially developing disease biomarkers and improving our understanding of the specific role of the choriocapillaris in disease progression.
Visualizing the choriocapillaris has historically proven challenging. Histological studies, while effective, do not permit functional or longitudinal studies [1,14]. Fundus photography and scanning light ophthalmoscopy cannot penetrate the RPE to visualize subretinal tissue like the choriocapillaris. The fenestrated endothelium of the choriocapillaris leaks most conventional dyes, making florescent dye angiography challenging [15,16]. By contrast, optical coherence tomography angiography (OCTA) is designed for noninvasive cross-sectional vascular imaging and is better suited to image the choriocapillaris [17][18][19]. However, numerous challenges exist. The RPE scatters and distorts the scanning beam passing through it [20]. This has been largely overcome using systems containing swept-source lasers (SS-OCT) with wavelengths centered around 1 µm to reduce scattering, which is stronger at shorter wavelengths [21,22]. Conventional SS-OCTA only provides a grainy visualization of the macular region of the choriocapillaris referred to as a confluent "meshwork" where the presence and absence of signal are respectively assumed to be capillaries and avascular regions dubbed "flow voids" [21][22][23][24]. This lateral resolution limitation exists in part because the OCTA's scanning laser has insufficient scanning speed to resolve the high flow-speed blood of the choriocapillaris, resulting in blood movement saturating the OCTA signal [25,26]. In the past several years, SS-OCTA systems utilizing rapidly (>400 kHz) scanning lasers, have overcome this hurdle, producing en face images of choriocapillaris with significantly improved vasculature visibility [26][27][28]. Concurrently, improved angiographic processing methods, as well as image refinement techniques such as volume averaging and deep learning, have been used to improve choriocapillaris visualization in conventional low-speed OCTA systems [29][30][31][32][33]. The combination of high-speed SS-OCTA with innovative post-processing is likely to improve image quality and facilitate better identification of disease-related changes.
In an earlier study, we developed an ultrahigh-speed SS-OCT system (UHS-OCT) with a 1.6 MHz scanning laser capable of reliably resolving the vasculature within even the highly dense macular choriocapillaris of healthy retinas [26]. In this study, we present a method for further improving the visibility of the choriocapillaris in these images and propose a set of five anatomical metrics that can be used to quantitatively summarize the images. These metrics are flow void radius, flow void area, vessel radius, branch-point to branch-point vessel length, and branches per branch-point. We validated this analytical method by implementing it on OCTA image analysis of four healthy retinas and one retina with GA. The data we collected suggest that changes in these metrics may be biomarkers for GA.

Subjects
Five human subjects were imaged under a protocol approved by the UC Davis Institutional Review Board. Written informed consent was obtained for all subjects prior to imaging. Four normal subjects and one patient diagnosed with GA were imaged in one eye with a custom UHS-OCT system. All eyes imaged were dilated with 1% Tropicamide 2.5% phenylephrine eye drops.

Image acquisition and angiographic processing
The system used to image subjects was a custom-built instrument described previously [26], whose core component was the Fourier-domain mode-locked (FDML) swept-source laser (FDML-1060-750-4B-APC, OptoRes GmbH, Munich, Germany). The central wavelength of the laser source was 1063 nm and sweeps across a bandwidth of 83 nm at a rate of 1.64 million times per second. This rate is significantly faster than most commercial OCT systems and allows for a rapid B-scan rate of 2 kHz [26]. Acquired OCT volumes consisted of 392 A-scans per horizontal B-scan, 10 sequential B-scans at each location, and 392 vertical positions. Final volumes were 392 × 392 × 512 voxels, and the lateral and axial resolutions in the eye were 9.1 µm and 7.2 µm, respectively. The lateral imaging area was 1.2 mm × 1.2 mm, approximately 4 • in a standard 60 D. Data processing to generate the OCT and OCTA volumes, automatically segment the layers, and generate the en face vascular images was performed with techniques demonstrated previously [26,34]. Subjects 1-4 were respectively imaged in 10, 15, 15, and 8 locations. Images were located between 3 • temporal to 9 • nasal and 3 • inferior to 3 • superior to the fovea with 3 degrees separation between each image center. The subject with GA was imaged in 4 overlapping areas starting at the fovea and moving nasally. These 4 volumes were manually stitched together and a ∼1.0 mm × 1.1 mm section with well-defined lesions was selected for analysis.

En face image processing and skeletonizing
To further resolve the capillary plexus, the OCT volumes were processed with a custom Python 2.7 script, as we show in Code 1 (Ref. [35]). Each flattened volume was axially sliced into a series of en face images of the retina. For each, the ∼10 µm layer containing the choriocapillaris was manually identified and averaged to produce a single 8-bit image (pixel value range [0,255] ADU) in which vessels were visible. High-frequency noise was removed from the image through thresholding (cutoff: 30 ADU) and repeated small-kernel gaussian blurring (repeats: 3, kernel size: 3 pixels, σ: 3 pixels). Then, this image was run through a custom-designed processing pipeline similar to a max-min power-weighted filter [36].
In this algorithm which we have dubbed "local min-max normalized masking" (LMNM), the image was broken into single-pixel-wide horizontal and vertical slices. For each slice, the locations of all minimums and maximums were graphically identified. To further reduce noise, maximums within 3 ADU of either of their adjacent minimums' values were excluded; and, for adjacent maximum pairs within 5 pixels of another, the smaller peak was removed. For each adjacent minimum-maximum-pair, pixels between the two were assigned either 255 ADU or 0 ADU depending on whether they were respectively above or below an 85% threshold of the difference between the maximum and minimum pixel value. The binarized horizontal slices were recombined into an image highlighting the vertical component of the vasculature. The same was done with the vertical slices, identifying horizontal vasculature. These two representations were combined using element-wise maximum filtering, producing a single image masking all vasculature. As the choriocapillaris is a singular connected network, it should be represented by a single unified mask. So, to remove masking errors, masked objects were categorized by size and those less than 800 pixels 2 in area were removed.
Next, the LMNM image was morphologically skeletonized. Polynomial-fitting branch-path extrapolation was used to complete partial branches over 6 pixels in length and prune branches under six pixels in length to repair insufficient and excessive areas of skeletonization respectively. The capillary skeletonization was converted into a partially connected graph datatype whose vertices were vessel branch points and edges were vessels. All thresholds detailed above were manually optimized to maximize qualitative agreement between the image and the graph representation of the choriocapillaris. Key steps in the processing pipeline are shown in Fig. 1.

Processing pipeline verification
To assess the reliability of the post-processing pipeline, the skeletonized vessel maps were compared to vessel traces created by three expert OCT technicians. Two members of the UC Davis Visual Field and OCT Reading Center and one member of our lab's clinical staff were selected because of their regular interaction with images of retinal vasculature including OCT angiograms. We developed custom software in Python 2.7 to record manual tracing. Each grader used it to trace 25% of each of four choriocapillaris angiograms from Subject 1. Resulting traces were each compared to vessel skeletonizations of the same choriocapillaris region automatically generated by the LMNM pipeline described above. Every point on the automated vessel skeletonization was characterized as a hit, miss, or false positive. Points on the automated vessel skeletonization existing within five pixels (to allow for the imprecision of manual point-and-click based tracing) of a pixel in the human-made vessel skeletonization, were marked as a hit if the line between the two matched pixels was contained entirely within the masked area of region's LMNM mask. All other points on the automated vessel skeletonization were defined as false positives. All pixels in the manual vessel skeletonization which were not matched to a hit on the automated vessel skeletonization were marked as a miss. By dividing the number of hit, miss, and false positive pixels by the sum of pixels in the human-generated skeletonization, hit, miss, and false positive ratios were determined to assess the reliability of the automated skeletonization and LMNM processing which produced it. The time it took the technicians and program to respectively perform the skeletonization was recorded.

Anatomical metrics generation
For each image, five anatomical metrics of the choriocapillaris were computed through a combined analysis of the LMNM vessel mask, vessel skeletonization, and vessel-connection graph. These metrics include flow void radius, flow void area, vessel radius, branch-point to branch-point vessel length, and branches per branch-point. All metrics were calculated locally and then averaged. Flow void radii were determined by inverting the binary vessel mask, converting the resultant to a Euclidian distance map, morphologically identifying each object in the image, and then recording the maximum pixel value of each object. Flow void area was determined by inverting the binary vessel mask, morphologically identifying each object in the image, and recording their pixel area. To determine vessel radius, the binary vessel mask was converted to a Euclidian distance mask and pixels aligned with the vessel skeletonization (which bisects the vessel mask) had their distance values recorded. Vessel length was determined by the length of edges in the vessel-connection graph. Branch number was determined by the number of edges connected to each node in the vessel-connection graph. All metrics except branch number were converted to real length with the known 3.1 µm/pixel scanning density of the system. Metric averages were calculated for each image. Before they were calculated for the GA subject, poorly segmented retinal vasculature was removed with manual masking, and a binary mask of GA lesions was generated by manually tracing areas of low average OCT signal in the RPE layer. These regions displayed increased brightness in the OCTA choroid image indicative of signal hypertransmission through the RPE, a staging standard for GA [37]. Averages of anatomical metrics were calculated separately for the region of the choriocapillaris under the mask and the region outside it.

Visual effectiveness of UHS-OCTA and post-processing
All 4 healthy subjects were imaged using both a Zeiss AngioPlex OCT system and our group's UHS-OCT. Imaging sessions took between 0.5 and 1 hour to complete. Imaging and LMNM post-processing of Subject 1 was representative of all four normal subjects and will therefore be the source of the images displayed herein. Imaging of Subject 1-4 at the macula with a Zeiss AngioPlex OCT system (Fig. 2(A)), at the layer of the choriocapillaris (Fig. 2B), show some variation in angiographic signal, but no distinct vessels, even when zoomed in to a 1.2 mm × 1.2 mm area (Fig. 2C). When the same 1.2 mm × 1.2 mm areas were imaged with our custom system (Fig. 2D), individual capillaries were visible. In a 0.4 mm × 0.4 mm subsection of these images, (Fig. 1(A)) individual vessels are visible but grainy. LMNM processing appears to identify local vessel cross sections as maxima and converts them into binary maps isolating horizontal ( Fig. 1(B)') and vertical ( Fig. 1(B)") elements of the vasculature. When both binary masks are combined, visually clear, minimal-noise representations of the vessels (Fig. 1(C)) are the result. The subsequent skeletonizations (Fig. 1(D)), when overlaid with the original flattened choriocapillaris images they derive from (Fig. 1(E)), appear to accurately trace vasculature.

Processing pipeline verification
While each 1.2 mm × 1.2 mm choriocapillaris angiogram and resultant vessel tracing (Fig. 3) appears visually like vasculature, we confirmed this qualitative evaluation by comparing the programmatic traces to human-created traces of the same regions. Three skilled technicians were asked to independently trace a total of sixteen 0.3 mm×0.3 mm sections of choriocapillaris images from 4 locations on Subject 1's retina. These same regions were processed and skeletonized via LMNM. It took the technicians an average of 20.9 ± 4.3 min to trace each section while it took the processing pipeline an average of 1.2 ± 0.24 s to complete the same task on a laptop with an i7700HQ CPU running at 3.1 GHz. The resulting traces were compared via signal detection theory [38], in which hits were defined as locations where both methods found vasculature, misses were defined as locations where only technicians found vasculature, and false positives were defined as locations where only the post-processing algorithm found vasculature. The average per-subsection hit, miss, false positive percentages were 89.3 ± 1.8%, 10.7 ± 1.8%, and 6.8 ± 4.1% respectively. When an equivalent analysis was performed to compare the trace of any two of the three technicians, the average per-subsection hit, miss, and false positive percentages were 86.7 ± 3.1%, 13.3 ± 3.0%, and 4.0 ± 1.3% respectively.

Anatomical metrics of healthy retinas
Anatomical metrics were extracted from the post-processed choriocapillaris images from the four healthy subjects. Flow void radius, flow void area, vessel radius, branch-point to branch-point vessel length, and branches per branch-point were quantified. After assigning values for each metric to each pixel in the image, they were visualized as maps for all healthy subjects. Maps of Subject 1's anatomical parameters, except branch number, are shown in Fig. 4. Average values for these metrics were collected for each image and per-subject averages were calculated (Table 1). There is no statistically significant (p>0.05) variance in any of these metrics between the subjects imaged. The value of all five metrics was compared across all locations imaged and there was no statistically significant (p>0.05) correlation between any parameter and retinal eccentricity.

Anatomical metrics of retina with geographic atrophy
The atrophic lesions of the GA subject were identified in projections of the RPE layer from the volumetric OCT intensity image and manually traced (Fig. 5(A)). Of the seven lesions traced, two were excluded because they were significantly obscured by shadows of retinal vasculature. The choriocapillaris en face angiogram showed poorer vasculature resolution than aniograms of healthy subjects but vessels can still be distinguished (Fig. 5(B)). Anatomical metrics were quantified separately for the choriocapillaris both beneath and distal to these lesions. The proportionally high OCTA signal intensity of the choroid under lesions supported their characterization (Fig. 5(C)). In Table 2 and 3, these metrics are reported for the lesion-free area of the image, the average of all analyzed lesions, and each of the 5 lesions analyzed. As compared to the lesion-free area, there were statistically significant (p=0.05) increases in flow void radius (Fig. 5(D)) and decreases in vessel radius (Fig. 5(E)) in four of the five analyzed lesions.

Discussion
Our custom, 1.6 MHz scanning-rate, 1 µm wavelength OCTA reveals vascular structure in the choriocapillaris seen in neither commercial OCTA systems or conventional, sub − 400 kHz scanning-rate, 800 nm wavelength research-grade systems [26][27][28]. The images collected from normal subjects for this study demonstrate that our system is capable of resolving the vessels of even the dense macular choriocapillaris [26].
Previous work has shown that automated statistical analysis methods can detect changes in choriocapillaris density even from images which are not well-resolved [6,[8][9][10][11][12]. However, we believe that the improved visibility offered by our high-speed OCTA system lends itself to the development of advanced metrics based on the visible structure. The benefit of such metrics is twofold. First, they may provide novel insight into the mechanism of choriocapillaris reduction in atrophic retinopathies, be it capillary thinning, withdrawal, organized morphological restructuring, or otherwise. Second, they may provide greater sensitivity to vasculature loss and therefore a more sensitive biomarker for disease progression. What differentiates our laboratory's system from other published ultrahigh-speed OCTs and likely provides this increased imaging precision is a significantly increased B-scan rate. When images of the choriocapillaris generated with this system are processed via local min-max normalized masking, vasculature is better differentiated from noise and the resulting network of vessels bears visual similarity with histology [39,40]. While visual appraisal of what is and is not a vessel is subjective, having an objective means of defining vasculature is advantageous for both research and clinical applications involving visualization of the choriocapillaris. This algorithm has potential applications in images of retinal vasculature as well, although its parameters would likely require tuning. To make LMNM processing available for future studies, it has been included within the supplementals of this report [35].
To reaffirm our subjective visual appraisal of LMNM post-processed images, signal reduction theory was used to compare choriocapillaris vessel skeletonizations generated programmatically to those made by experienced technicians. Correspondences between the automated traces and manual traces were at least as high as correspondence between manual traces from different technicians-90% in both cases. As should be expected, automated tracing was significantly faster and less labor intensive than manual tracing, taking 1/1000th the time it takes a human, bringing the time frame for tracing a full 1.2 mm × 1.2 mm choriocapillaris en face image from 5.6 hours to 18.6 seconds. At that speed, our approach is practical for studies with large patient populations. Although we did not explore reproducibility of LMNM tracing directly, we expect that it, like most automated data analysis, will be at least as reproducible as manual tracing.
Next, we performed a proof-of-concept that LMNM post-processing and skeletonization possess novel analytical applications. We used images generated by LMNM post-processing to assess choriocapillaris flow void radius, flow void area, vessel radius, vessel length and branch number of the four normal subjects. On average, flow void radii were 9.8 ± 0.7 µm. If flow voids were approximately circular, an average area of ∼300 µm 2 would be expected. An average flow void area of 749 ± 110 µm 2 was observed, indicating that flow voids are elliptical, which increases the ratio of their perimeters to areas. Increasing this ratio for avascular regions may aid oxygen diffusion, maximizing metabolic perfusion and oxygenation of the retina. Vessel radius was determined to be 5.0 ± 0.3 µm on average. While previous histological studies have reported vessel radii for the choriocapillaris of approximately 10 µm, those same studies record radii of choriocapillaris vascular lumens as 4 µm − 8 µm [2,14,41]. Since OCTA visualizes areas of blood flow rather than tissue, apparent vessels visualized by OCTA represent vascular lumens. Thus, the average "vessel radius" reported in this study is consistent with histology.
While there was no significant variance of the proposed metrics either between healthy subjects or across imaged eccentricities, there were interesting local correlations between metrics. In locations where a branch point's number of branches exceeds 3, vessel length is locally decreased, and vessel radius is locally increased as compared to adjacent areas (Fig. 4). These three changes are consistent with histological images of feeder vessels, where choroidal vessels feed into the choriocapillaris plexus. Such local variations could be a potential means of identifying feeder vessels, but further investigation is required.
To determine whether LMNM post-processing and skeletonization of our UHS-OCT angiograms can identify potential biomarkers for atrophic retinopathies, a retina with GA was studied. Since GA results in RPE dropout and hypertransmission of OCTA signal from the choroid [37], the seven GA lesions studied were selected because they presented these traits. As compared to average anatomical metrics for the four healthy retinas studied, the GA subject's lesion-free tissue was less resolved and less similar to histology. This part of the choriocapillaris had statistically equivalent flow void area, vessel radius, vessel length, and branch number, but had a reduced average flow void radius ( Table 2). The differences in appearance and flow void radius could be attributed to age-related morphological changes, disease-related morphological changes, or natural inter-retinal anatomical variation, but could also be a consequence of reduced image quality. To mitigate the potentially confounding influence of image quality, we will focus only on anatomical differences between sub-lesion and lesion-free regions of the GA subject's choriocapillaris. As subsequent studies use novel OCT systems to further resolve images of disease-affected vascular tissue, reliable values for anatomical metrics may be obtained using LMNM processing.
Five of the seven GA lesions were considered for analysis. Two lesions were excluded from further processing due to obstruction from retinal vasculature shadows. Under all the processed lesions, the choriocapillaris was visibly disrupted. Quantification revealed that four out of five of these regions had a minor but statistically significant (exceeding 95% confidence intervals) elevation in average flow void radii and reduction in vessel radii. The largest lesion, lesion b in Fig. 5, also showed statistically significant (p<0.05) elevations in average flow void area and vessel length. These observations could respectively be explained by partial vessel constriction in developing lesions and vascular withdrawal in late-stage lesions, though this must be verified by further study. These findings suggest both that local changes in flow void radius and vessel radius may be a component of the choriocapillaris thinning reported in GA [27], and that skeletonization of LMNM processed images is a viable means of analyzing potential biomarkers. Future studies may use this approach to develop clinically relevant methods to identify GA-affected choriocapillaris tissue.
There are limitations to this study. Firstly, our sample size is only four normal subjects and one GA subject. We present these images and associated quantification as a proof-of-concept, motivating further work in high-speed OCTA, LMNM processing, and vasculature skeletonization as a means of quantifying the anatomy of the choriocapillaris. Secondly, the image acquisition and post-processing possess several potential but minor sources of error. Angiograms produced by OCT are subject to a degree of optical blur proportional to their focused spot size. Previous work has shown that the dominant spatial frequency present in the choriocapillaris is well within the optical bandwidth of the system, given it's 9.1 µm lateral resolution [26]. Thus, the effects of blur may be slightly impacting spatial metrics such as vessel radius, but are not significantly impacting the vessel skeletonization or the topological metrics such as branch number calculated from it. Similarly, all masking-based image processing techniques, LMNM included, can potentially over or under sample signal based on user-defined parameters of the algorithm. The parameters we manually optimized produce capillary skeletonizations like those traced by technicians and produce anatomical metrics in agreement with histology. But, in future studies, information in the OCT intensity volume may be usable to define more objective post-processing parameters. Due to these potential errors, the absolute anatomical metric values presented within this report may be affected, but the comparative appraisal of these metrics, both beneath and distal to GA lesions, should be largely unaffected. Future studies seeking to use UHS-OCT combined with LMNM processing should aim to quantify atrophic retinopathies comparatively rather than absolutely.

Conclusions
This study demonstrates the combination of ultrahigh-speed OCTA imaging of the choriocapillaris with local min-max normalized masking, a novel algorithmic tool for processing and quantifying anatomy in en face images of the choriocapillaris. When angiograms of the choriocapillaris are processed with LMNM, a visually clear, easily quantified representation of the choriocapillaris results. To confirm the reliability of LMNM, we compared choriocapillaris tracings produced by it to traces made by humans and found them to correspond well. Then, we demonstrated that LMNM processed images can be easily quantified to extract values for anatomical metrics of the choriocapillaris and that these metrics are significantly different under most lesions in images of a retina affected by GA. With further refinement the combination of ultrahigh-speed OCTA and LMNM image processing may yield improved visualization and quantification of the morphological changes in the choriocapillaris associated with a variety of atrophic retinopathies.