1.

Introduction

Optical coherence tomography angiography (OCT-A) is an emerging imaging modality with which the retinal circulation can be visualized by computing the decorrelation signal on a pixel-by-pixel basis. Variants of OCT-A methods have been described in recent review articles.^1–4 We have recently validated the speckle variance approach to OCT-A against standard invasive techniques, such as fluorescein angiography (FA),^5,6 in which only the superficial capillaries can be distinguished due to excessive choroidal fluorescence⁷ and ex vivo histological analyses.^6,8–10

The noninvasive, in vivo visualization of the retinal microvasculature using OCT-A can be instrumental in studying the onset and development of retinal vascular diseases. For example, OCT-A has enabled the visualization of the deep plexus layer and furthered the understanding of diseases, such as paracentral acute middle maculopathy^11–13 and diabetic retinopathy.^14,15 Quantitative measurements, such as capillary density, can be used to stratify the risk of disease progression, visual loss, and also for monitoring the course of disease.^16,17 Due to projection artifact and poor contrast, it is often difficult to trace individual vessels in this layer when only one en face image is visualized. An additional challenge to this end is the small dimension and pulsatile flow of the retinal capillaries, making them less consistently visible and difficult to distinguish from the speckle noise relative to larger vessels. This limits the detection sensitivity for changes in the retinal microvascular circulation due to diseases, aging, or treatment. Methods for reliable visualization of the microvasculature in the OCT-A images are required for studies conducting longitudinal and cross-sectional quantitative analysis.

Our work focuses on improving the visible definition of retinal microvasculature in OCT-A by motion correction, registration, and averaging of sequentially acquired images. Detailed, high-quality OCT-A images are needed for clinical studies, such as comparisons of OCT-A with histology and fundus photography FA, and studying the shunting of vessels in a focal area, such as the inner ring of vessels in the foveal avascular zone (FAZ) or in glaucomatous focal defects.¹⁸

Serially acquiring and averaging multiple OCT-A images can be an effective solution for confirming the presence or absence of capillaries as the discontinuous appearance of the capillary vessels is beyond improvement simply by just applying image filtering.^3,19 A crucial step in the serial acquisition approach is the registration of multiple OCT-A images, the difficulty of which is compounded by the fact that an OCT-A image is acquired over multiple seconds and, thus, particularly susceptible to motion artifacts. The registration of sequential B-scans²⁰ can aid in attenuating small motion artifacts, but not the larger motion artifacts associated with imaging subjects with pathologies. OCT-A with eye-tracking has been implemented in commercial retinal imaging systems, although this increases hardware cost and complexity on which the sensitivity and reliability of motion detection also depend. Previous works on posthoc motion artifact removal by en face summed volume projection have been reported in the literature, see for example Refs. 3, 21, and 22. In the work by Hendargo et al.,²¹ two to three sets of orthogonal ($x$- and $y$-fast) volumes were acquired and divided into motion-free strips. The visualization and contrast of the vessels were improved by multiresolution Gabor filtering, and the strips were registered one-by-one, first globally by $x$- and $y$-translation that maximized the correlation in the overlapping region and locally by B-spline free-form deformation in the overlapping region. Zang et al.²² did not acquire orthogonal data sets, but instead serially acquired two OCT-A volumes in the same scan orientation that were divided into parallel motion-free strips. The strips were first registered by $x$- and $y$-translation and rotation that minimized the squared difference of “large vessels,” which were defined in the paper as pixels with decorrelation value greater than 1.3 times the mean value. This was followed by B-spline free-form deformation on “small vessels,” defined as pixels with decorrelation value less than 1.3 times and greater than 0.6 times the mean value. Both groups presented mosaicking of OCT-A images into wide-field views, which has been reported in other works as well.^19,23

In this work, we present averaging of up to 10 serially acquired OCT-A images with parallel strip-wise microsaccadic noise removal and localized nonrigid registration. Unlike the previous two methods,^21,22 which concentrated motion artifact removal and wide-field imaging, our purpose was to improve the contrast and signal to background of the capillaries in focal regions. The details of our methodology are presented below. In brief, the serially acquired OCT-A images were divided into microsaccade-free strips. The target strips were first aligned to a template image by $x$- and $y$-translation based on maximum cross-correlation, followed by affine registration using scale-invariant feature transform (SIFT),²⁴ a feature extraction method robust to scaling, orientation changes, illumination changes, and affine distortions.

The image warping and local distortion due to slower eye movements are less obvious and more difficult to model than the strong stripe artifacts from microsaccadic motion. Instead of free-form deformation,^21,22 our approach optimized the intensity value at each pixel location as the average of the values from each overlapping strip determined by translation and rotation of a windowed region in each strip. Thus pixel-wise correspondence across multiple OCT-A images was found by local neighborhood matching.

The remainder of this report is organized as follows. Section 2 describes our processing algorithm for OCT-A image averaging in detail, as well as quantitative metrics to evaluate the improvements to the image quality. The algorithm was tested on OCT-A images of six healthy volunteers, with the vessel visibility improvement qualitatively demonstrated in all, superficial, and deep plexus layers. Quantitatively, we evaluated the algorithm performance by contrast-to-noise ratio (CNR) and signal-to-noise ratio (SNR) with the background speckle noise information from the FAZ. We end with a discussion of the image quality improvements using our method of averaging serially acquired OCT-A images.

2.

Methods

All subject recruitment and imaging took place at the Eye Care Centre of Vancouver General Hospital. The project protocol was approved by the Research Ethics Boards at the University of British Columbia, Simon Fraser University, and Vancouver General Hospital, and performed in accordance with the tenets of the Declaration of Helsinki. Written informed consent was obtained from all subjects.

2.3.

Angiogram Registration

The algorithm overview is shown in Fig. 1. The 10 serially acquired en face images of all retinal layers were divided into microsaccade-free strips, which were then registered to a template image, first using rigid registration for the coarse alignment, followed by nonrigid registration for finer features. Transforms applied to the en face image of the full retinal thickness were then applied to both the superficial and deep layer angiograms.

Fig. 1

Overview of the strip-based registration algorithm for multiple serially acquired OCT-A images. Representative images are used to demonstrate the algorithm in Fig. 2.

2.3.1.

Microsaccade-free strip generation

For each eye, a microsaccade-free image from the 10 en face images was chosen as the template image. In the case that all images contained microsaccadic motion artifacts, a template was generated by stitching together microsaccade-free strips using the registration methods discussed below.

After the template image was chosen/generated, the remaining images were divided into strips between positions in the image corresponding to where the patient fixation was lost, which appeared as vertical white stripes in the en face image, as shown in Fig. 2. Strips less than 40 pixels wide often contained large drift artifact and were therefore discarded. If multiple microsaccade-free images existed per eye, the first was selected as the template, and the rest were divided into three equal-sized strips for registration. Each strip was zero padded to match the size of the template image and coarsely aligned to the template by $x$- and $y$-translation using maximum cross-correlation.²⁹

Fig. 2

Demonstration of the image stripping, coarse translation, affine registration, and nonrigid registration steps of the proposed algorithm. The template image (green) and registered strip (magenta) are shown as composite images where white regions indicate where the two images have the same intensities. The areas under the red, orange, and yellow boxes are further explored in Fig. 3.

2.3.2.

Strip-based affine registration

SIFT keypoints were automatically extracted from both the template image and each strip to be registered.²⁴ Briefly, keypoints are the locations of local scale-space extrema in the difference-of-Gaussian function convolved with the image. Further refinement to the keypoints can be made by assigning each keypoint an orientation to achieve invariance to image rotation. Finally, a local image descriptor is assigned to each keypoint using the image location, scale, and orientation as found above. Readers are encouraged to refer to Ref. 24 for a more detailed description of the SIFT algorithm.

As the SIFT feature descriptor is invariant to uniform scaling and orientation, it is ideal for identifying matching keypoints in noisy or speckled images, such as OCT-A angiograms. The calculation of Euclidean distances in MATLAB^® is computationally expensive, and therefore matching keypoints between the template and strip were identified as the closest corresponding keypoints by a small angle approximation to the Euclidean distance. Keypoints were considered matching if the ratio of the vector angles from the nearest to the second nearest match was less than a threshold value of 0.75. As the image had been coarsely aligned in the previous step, a second check was included to ensure that the matched keypoints were no more than 40 pixels distant in the $x$- or $y$-direction.

All strips that had a minimum of four matched keypoints were then transformed using an affine transform estimated using the matching keypoints as inputs to the estimate geometric transform function in MATLAB^®. This function iteratively compares an affine transformation using three randomly selected keypoints, where the transformation with the smaller distance metric calculated using the M-estimator sample consensus algorithm is used as the transformation matrix for the next comparison. The maximum number of random trials for finding the inliers was set to 5000 for improved robustness.

2.3.3.

Strip-based nonrigid registration

The vertical white lines in the target image in Fig. 2 mark the image discontinuities due to microsaccades accounted for by the strip-based affine registration; however, localized mismatch still remains in the aligned images after this affine registration step. The next step in our algorithm is to compensate for the smoother tremor and drift motions represented by image warping and distortion, by using nonrigid registration. Prior to nonrigid registration, a $2\times 2$ averaging filter was applied to both the template and the aligned strip to smooth any fine speckle that may affect the nonrigid registration. The template and aligned strips were both then zero padded by 15 pixels. For each pixel, in the strip, the normalized cross-correlation³¹ was calculated, defined as

Eq. (3)

$$x{\mathrm{corr}}_{\text{norm}}(s,t)=\frac{\sum _{x,y}[f(x,y)-\overline{f_{s,t}}][m(x-s,y-t)-\overline{m}]}{\sqrt{\sum _{x,y}{[f(x,y)-\overline{f_{s,t}}]}^2\sum _{x,y}{[m(x-s,y-t)-\overline{m}]}^2}},$$

where $f(x,y)$ is the $29\times 29\text{pixel}$ matrix field centered on the $(x,y)$ pixel of the template image, $\overline{f_{s,t}}$ is the mean of the image in the region under the mask, $m$ is the $15\times 15\text{pixel}$ mask matrix centered on the pixel of the strip, and $\overline{m}$ is the mean of the mask. This was also done for $-15\mathrm{deg}$, $-10\mathrm{deg}$, $-5\mathrm{deg}$, 5 deg, 10 deg, and 15 deg rotated field matrices. The pixel located at the index of the maximum normalized cross-correlation was then used as the registered pixel for the strip. Figure 2 shows a pictorial schematic of the registration pipeline described in this section. A smaller field of view demonstration of the coarse, affine, and nonrigid registration steps is shown in Fig. 3. The stack of registered strips could then either be combined by taking the mean or median to generate a higher quality image.

Fig. 3

Comparison of three different strips registered to the same template image using (a) coarse translation, (b) affine registration, and (c) nonrigid registration. The template image (green) and registered strip (magenta) are shown as composite images where white regions indicate where the two images have the same intensities.

3.

Results

A total of 10 eyes from six healthy volunteers (4 males and 2 females) aged $36.8\pm 9.3$ years were acquired according to the imaging protocol. A comparison of the template image and the final averaged OCT-A images for all retinal layers, as well as the superficial and deep vascular layers is shown in Figs. 4 and 5. In the template images, the vessels near the FAZ are relatively clear; however, it becomes harder to differentiate the vessels further toward the periphery. In contrast, the vessels in the averaged images are clearly seen throughout. Improvement in vessel visibility is particularly marked in the deep layer, where the OCT signal strength is weaker. Qualitatively the median images appear sharper than the mean images as the median averaging acts as a speckle reducing filter. However, the mean images appear more smooth than the corresponding median image.

Fig. 4

Template image, mean, and median averaged images (all retinal layers, superficial, and deep plexus) for subject 3 oculus dexter, a healthy male subject, 29 years of age.

Fig. 5

Template image, mean, and median images (all retinal layers, superficial, and deep plexus) for subject 6 oculus sinister, a healthy male subject, 56 years of age.

For quantitative comparisons of the template and final averaged images, the average CNR and SNR of the images were calculated. The average CNR of the angiograms with all retinal layers increased from $0.52\pm 0.22\mathrm{dB}$ using the template images to $0.77\pm 0.25\mathrm{dB}$ with the mean images and $0.75\pm 0.24\mathrm{dB}$ with the median images. Additionally, the average SNR of the angiograms with all retinal layers increased from $19.58\pm 4.04\mathrm{dB}$ using the template images to $25.05\pm 4.73\mathrm{dB}$ with the mean images and $25.02\pm 4.89\mathrm{dB}$ with the median images. The mean improvement of both the CNR and SNR was statistically significant ($p<0.01$) using a paired $t$-test.

To evaluate the change in perceptual quality per strip, the SSIM was calculated on each incremental averaged image of the template and registered strips. Although 10 volumes were acquired per eye, the number of microsaccades and strips less than 40 pixels in the corresponding en face images was different for all eyes, and therefore the number of strips used to generate the averaged images was not necessarily equal. The mean number of strips per eye was $21\pm 7$ strips. As seen in Fig. 6, the SSIM values show a rapid increase as the first few strips are registered and applied to the template image, and then the rate of improvement slows with additional registered strips. This trend was observed in both the mean and median averaged images.

Fig. 6

SSIM values for incremental averaged images for all eyes for the all retinal layers, superficial layers, and deep layers.

4.

Discussion

The major findings in this paper are as follows: averaging multiple registered sequentially acquired OCT-A images (1) qualitatively enhances the visualization of the retinal microvasculature networks, (2) increases the SNR and CNR of the angiograms, and (3) increases the perceptual visual quality when using SSIM as a metric.

After averaging multiple en face images, the vessels of the deeper capillary plexus are more readily identified, making quantification more reliable and thereby facilitating investigation of its role in the pathophysiology of retinal vascular disease. Although minimal projection artifact can still be seen in Fig. 4 corresponding to the larger superficial vessels, the overall qualitative condition of the en face images is improved.

The SNR and CNR both increased significantly by averaging the individual strips. Although both the SNR and CNR of the mean images are larger than those of the median images, there is no significant difference between the mean and median, and therefore no recommendation of an averaging method can be made based on these metrics.

The SSIM is a full reference metric where the final averaged image was taken to be the perfect quality reference image. Note that each of the 10 volumetric sets for each subject was divided differently, depending on the motion in a particular acquisition, and therefore a different number of strips were used for each averaged reference image. SSIM is a relative metric, and an SSIM value of 1 indicates that the image has reached the same quality as the reference image for that data set. The SSIM values cannot be compared across different images, for example, an SSIM value of 1 for different images does not indicate that they are all of the same quality. As shown in Fig. 6, the SSIM increases with each additional registered strip that is averaged to the template image. As expected when averaging images, the first few strips applied to the template affected the SSIM the most whereas the later strips provided only modest improvement to the SSIM. The deep plexus showed the greatest increase overall. By increasing the visibility of individual vessels, this technique has the potential to improve automated segmentation results thereby improving our ability to quantify capillary density in normal and diseased states.

Although we have demonstrated the ability to enhance the visualization of the retinal plexuses through averaging multiple sequentially acquired OCT-A images, we acknowledge several limitations of this work. This study assessed only relatively young subjects with clear ocular media and good fixation ability. The presence of media opacities in older subjects may limit the amount of capillary information that can be attained from images. Although the algorithm attenuates nonmicrosaccade motion in the registered strips, the template may contain distortions and image warping, which is not accounted for here.