Quantifying Variability in Longitudinal Peripapillary RNFL and Choroidal Layer Thickness Using Surface Based Registration of OCT Images.

PURPOSE: To assess within-subject variability of retinal nerve fiber layer (RNFL) and choroidal layer thickness in longitudinal repeat optical coherence tomography (OCT) images with point-to-point measurement comparison made using nonrigid surface registration.
METHODS: Nine repeat peripapillary OCT images were acquired over 3 weeks from 12 eyes of 6 young, healthy subjects using a 1060-nm prototype swept-source device. The RNFL, choroid and the Bruch's membrane opening (BMO) were segmented, and point-wise layer thicknesses and BMO dimensions were measured. For each eye, the layer surfaces of eight follow-up images were registered to those of the baseline image, first by rigid alignment using blood vessel projections and axial height and tilt correction, followed by nonrigid registration of currents-based diffeomorphisms algorithms. This mapped all follow-up measurements point-wise to the common baseline coordinate system, allowing for point-wise statistical analysis. Measurement variability was evaluated point-wise for layer thicknesses and BMO dimensions by time-standard deviation (tSD).
RESULTS: The intraclass correlation coefficients (ICCs) of BMO area and eccentricity were 0.993 and 0.972, respectively. Time-mean and tSD were computed point-wise for RNFL and choroidal thickness and color-mapped on the baseline surfaces. tSD was less than two coherence lengths of the system 2ℓ = 12 μm at most vertices. High RNFL thickness variability corresponded to the locations of retinal vessels, and choroidal thickness varied more than RNFL thickness.
CONCLUSIONS: Our registration-based end-to-end pipeline produced point-wise correspondence among time-series retinal and choroidal surfaces with high measurement repeatability (low variability). Blood vessels were found to be the main sources contributing to the normal variability of the RNFL thickness measure. The computational pipeline with a measurement of normal variability can be used in future longitudinal studies to identify changes that are above the threshold of normal point-wise variability and track localized changes in retinal layers in high spatial resolution. TRANSLATIONAL RELEVANCE: Using the registration-based approach presented in this study, longitudinal changes in retinal and choroidal layers can be detected with higher sensitivity and spatial precision.

Introduction

Optical coherence tomography (OCT) has been widely adopted in ophthalmology for noninvasive, in vivo imaging that provides cross-sectional and volumetric visualization of the inner structures of the eye. Measurement of various shape parameters in OCT images allows for quantitative analysis of the anatomy and investigation of structure—function relationships. Longitudinal studies of OCT images can be especially illuminating because it allows for a direct connection between changes in the anatomy and disease progression, aging, or treatment.

A challenge involved in analyzing multiple longitudinal OCT data is ensuring anatomical correspondence at different time points such that the measurements are indeed comparable. Given a parameter defined by a metric (thickness, area, depth) and a target structure (retinal nerve fiber layer [RNFL], Bruch's membrane opening [BMO], optic cup), ambiguity can arise from either the definitions or the method of measurement. Parameters such as BMO area, which are based on a relatively well-defined landmark, are less ambiguous, whereas parameters defined by image orientation and distance (for example, RNFL thickness in the inferior region) are inherently limited in precision and more variable by image acquisition and processing. Retinal layer thickness is often measured in a limited number of sample points or as an average in a predefined area or a sector. Such approaches aim to mitigate the effect of spatial mismatch between acquisitions and increase the measurement repeatability; however, these also reduce the measurement resolution and scope, and small, localized changes that could be clinically important indicators of disease or treatment may be undetected.

In parallel to the advances in OCT image acquisition and very high resolution longitudinal, within-subject imaging, there is a need for algorithms that can provide point-wise correspondence across such images and reveal the finer patterns of change over time that may be lost in the sectoral-averaging or limited sampling such as in circular peripapillary scans. Image registration techniques have been used in OCT for averaging consecutive scans[1-4] and motion artefact removal.[5-8] Rigid registration of longitudinal OCT images in 2D by manually segmented blood vessels[9] and in 3D using scale invariant feature transform (SIFT) features[10] has been explored. For nonrigid registration, rigid registration by blood vessels was followed by graph-cut based alignment of individual A-scans,[11] and in another work, global affine transformation by the fovea and the inner and outer retinal boundaries was performed and A-scans were smoothly deformed to refine retinal layer boundary matching.[12] The methods above do not provide one-to-one correspondence across surfaces now have the point-wise comparison of metrics such as RNFL or choroidal thickness demonstrated across multiple longitudinal measurements.

We have previously presented a computational framework for establishing point-to-point correspondence between retinal surfaces from multiple OCT images through nonrigid surface registration.[13,14] In Gibson et al.,[13] intensity-based 3D rigid and nonrigid registrations were performed on optic cup surfaces followed by spherical mapping and spherical demons algorithm. Lee et al.[14] extended upon this method by introducing current-based nonrigid registration of the surfaces prior to the demons registration. Our approach presents a generalized and versatile framework that can be used to track arbitrary measurements on retinal surfaces in high spatial resolution. These point-to-point correspondence techniques overcome the noise issue in high resolution tracking using naïve alignment, and are exquisitely suited to detecting localized changes via optimal mapping across multiple OCT data sets. The resulting measurements are, by their very design, also sensitive to errors that could be introduced in the measurement pipeline, such as segmentation or alignment errors. The measurements would also reflect each individual's day-to-day physiological variability, including the diurnal fluctuations in blood pressure and effects of sleep and diet. The magnitude and spatiality of these effects at each point-wise location on retina surfaces and whether some subjects are more variable than others are yet to be explored.

This paper presents a framework and results of longitudinal analysis of point-wise peripapillary measurements in healthy subjects using the nonrigid registration techniques presented in our previous work.[14] Starting with the working hypothesis that there is no variability in measurements such as RNFL/choroidal thickness made from longitudinal repeated scans, we demonstrate experiments to assess this normal-range short-term variability of the registration-based peripapillary measurements. The experiments were designed to focus on healthy subjects imaged repeatedly over a short span of 3 weeks to limit the observed changes to normative noise, excluding substantial and significant changes that may occur over a longer period of time. The study design enables selective quantification of variability from physiological factors as well as processing artefacts, such as image misregistration, from different stages of the processing pipeline, summarized in Figure 1. It is important to quantify the “noise floor” of the end-to-end system in order to establish the reference threshold for distinguishing and detecting changes due to disease progression or effectiveness of a treatment from measurement noise. The physiological and artefactual factors affecting the measurement variability at different stages of the pipeline are discussed with a view to quantifying these sources and potentially minimizing them in future studies.

Figure 1

Processing pipeline for registration-based repeatability analysis of optic nerve head morphometrics in OCT images. The images are acquired over the baseline time point t0 and follow-up time points t1–tT. The images are processed through motion correction and smoothing, and segmentation of retinal layers and BMO. Retinal surfaces are generated and cropped, and measurements are made in each surface in the original coordinates. All follow-up surfaces are registered to the baseline surface, yielding new mappings of the follow-up measurements in the baseline (template) coordinates. The measurements from t0 to tT are now in the common baseline coordinates and can be analyzed point-wise.

Methods

Participants and Image Acquisition

The study followed the tenets of the Declaration of Helsinki, and informed consents were obtained from the participants. Ethics review for the study was approved by Simon Fraser University (SFU) and University of British Columbia (UBC).

Twelve eyes from six healthy female participants (age: 29.5 ± 3.4) were included in the study. All OCT imaging was performed at the Eye Care Centre at Vancouver General Hospital over a period of 2 months. Each participant underwent a 24-2 or 30-2 visual field test and obtained a mean deviation (MD) greater than −1. The demographic information for each subject is included in Table 1. Each eye was imaged nine times over a period of 3 weeks, with the minimum of 22 hours and maximum of 6 days and 2 hours between any two consecutive scans. A total of 108 such OCT image volumes were included in the study.

Table 1

Subject Demographics

The images were acquired using a custom swept-source (SS) OCT system with a 1060-nm wavelength source designed and built by the Biomedical Optics Research Group at SFU.[15,16] Each image consisted of 400 B-scans, with 400 A-scans per B-scan and 1024 pixels per A-scan. The axial voxel resolution was 2.7 μm, the axial coherence length was 6 μm, and the lateral resolution ranged from 11.9 μm to 14.5 μm depending on the eye's axial length. The image field of view was centered at the optic nerve head. The A-scan rate of 100 kHz resulted in approximately 1.6 seconds of acquisition time per volume. In each imaging session, multiple volumes were acquired for each eye, which were later examined in B-scans, orthogonal scans, and enface projections for any severe acquisition artifacts. The image with the best image quality in terms of the signal-to-noise ratio and structural visibility, motion artifacts, and the relative central position of the optic nerve head in the image field of view was selected for each eye for further processing.

Image Processing, Segmentation, and Measurement

The images were corrected for axial motion using cross-correlation between adjacent frames. Three-dimensional bounded-variation smoothing was applied to reduce the effect of speckles and enhance the visibility of anatomical structures in the image. Next, the anterior and posterior surfaces of the RNFL and choroid were automatically segmented using a 3D graph-cut based algorithm.[17,18] BMO was manually segmented on 80 2D radial slices extracted from the image volume. For BMO measurements, an ellipse was fitted to the segmented BMO points using principal component analysis and least square fitting. All following references to BMO refer to this best-fit BMO ellipse. RNFL and choroid surfaces were cropped at 0.25 mm from the BMO ellipse since the layer boundary becomes ambiguous near the BMO. The surfaces were additionally cropped to annuli with the outer boundary at 1 mm from the BMO. This was performed to account for the difference in the BMO location relative the image frame across the images and the different lateral resolutions across the eyes depending on the axial lengths.

After segmentation, BMO area, BMO eccentricity, RNFL thickness, and choroid thickness were measured for each image. The BMO area and eccentricity were obtained from the fitted BMO ellipse. The RNFL and choroid thickness was measured at each point of the posterior surface as the closest 3D Euclidean distance to the anterior surface.

Registration

Prior to nonrigid registration, all corresponding surfaces of the same eye were rigidly aligned in three steps. The translation, scaling, and rotation in the lateral plane were corrected by intensity-based affine registration of the enface sum-voxel image of each follow-up to the baseline image. The axial positions of the surfaces were normalized by subtracting the mean axial coordinate in each surface. Finally, any remaining tilt difference between the baseline and follow-up surfaces was compensated by rotating the latter at the surface centroid, with the inner and cross products of the normal vectors of the two surfaces as the rotation angle and axis, respectively.

The nonrigid registration between the baseline and each of the preregistered follow-up surfaces was performed by current-based deformation followed by spherical demons registration.[14] First, the surfaces were represented by currents of geometric measure theory, which allows one to measure the closeness of two surfaces via a norm in a Hilbert space.[19] The optimal transformation of the target (follow-up) to the template (baseline) was then found as a smooth diffeomorphism, bringing the two surfaces into close proximity. Vertex-wise correspondence was achieved by spherical demons registration,[20] which maps all follow-up surfaces to the common coordinates of the baseline surface, allowing for vertex-wise temporal statistics of the layer thickness in the baseline coordinates.

Figure 2 visualizes an example of the RNFL registration. The first row shows the nine segmented and cropped RNFL posterior surfaces and thickness mapping acquired over 3 weeks, from the baseline time point t0 and subsequent follow-up time points t1 to t8. The follow-up RNFL surfaces were registered to the baseline RNFL, establishing point-to-point correspondence between the baseline RNFL and each follow-up RNFL. The second row of Figure 2 shows RNFL thickness of t1 to t8 mapped onto the registered RNFLs. The vertex-wise RNFL thickness difference map between the baseline RNFL and each follow-up is shown in the third row and in the last row with thresholding by the axial coherence length of the OCT system. Figure 3 shows the registration for the choroid of the same eye. As in Figure 2, the first row shows the nine segmented and cropped posterior choroidal surfaces and thickness maps acquired over 3 weeks, from the baseline time point t0 and subsequent follow-up time points t1 to t8. The follow-up choroidal surfaces were registered to the baseline choroid, establishing point-to-point correspondence between the baseline choroid and each follow-up choroid. The second row shows choroidal thickness of t1 to t8 mapped onto the registered choroids. The vertex-wise choroidal thickness difference map between the baseline choroid and each follow-up is shown in the third row and in the last row with thresholding by the axial coherence length of the OCT system. The figures show the registration step largely preserves the original surface topology and thickness maps without introducing artefacts (Row 2) while establishing the point-wise correspondence between the baseline and follow-up surfaces, resulting in spatially detailed time-difference maps (Row 3, 4).

Figure 2

Variability analysis of RNFL thickness over 3 weeks. Row 1: RNFL thickness maps on original RNFL surfaces acquired at nine time points in 3 weeks (t0: baseline, t1 – t8: follow-up). Row 2: Follow-up RNFLs were registered to the baseline, establishing point-to-point correspondence between each follow-up RNFL and the baseline RNFL serving as the common template. RNFL thickness is shown remapped onto the registered RNFLs. Row 3: Vertex-wise RNFL thickness difference between the baseline and each follow-up RNFL. Row 4: Difference maps in Row 3 are shown thresholded by the axial coherence length of the system ℓ = 6 μm. The registration step enables point-wise comparison of the follow-up data to the baseline while preserving the original surface topology and thickness maps.

Figure 3

Variability analysis of choroidal thickness over 3 weeks. Row 1: Choroidal thickness maps on original choroidal surfaces acquired at nine time points in 3 weeks (t0: baseline, t1 – t8: follow-up). Row 2: Follow-up choroids were registered to the baseline, establishing point-to-point correspondence between each follow-up choroid and the baseline choroid serving as the common template. Choroidal thickness is shown remapped onto the registered choroids. Row 3: Vertex-wise choroidal thickness difference between the baseline and each follow-up choroid. Row 4: Difference maps in Row 3 are shown thresholded by the axial coherence length of the system ℓ = 6 μm. The registration step enables point-wise comparison of the follow-up data to the baseline while preserving the original surface topology and thickness maps.

Statistical Analysis

Variability of retinal layer thickness measurement was measured by each vertex, eye, and for all subjects. For N eyes, each imaged T times, and L layers measured for thickness, the number of vertices of the baseline surface of layer l= 1…L and eye n = 1…N is given by . All registered follow-up surfaces consist of the same number of vertices by the point-to-point correspondence achieved by the registration process above. The number of vertices M ranged between 40,000 to 50,000, depending on the lateral voxel resolution of the image and size of the BMO, with the lateral voxel resolution of 400 × 400.

Descriptive time-statistics of mean and standard deviation were computed vertex-wise. For example, the time-mean (tMean) and time-standard deviation (tSD) of RNFL thickness of the right eye (OD) of Subject 1 at vertex m were given by where is the RNFL thickness of Subject 1OD at the vertex m at the time point t. For each eye and layer, the vertex-wise time-statistics tMean and tSD were visualized using surface color maps. Distribution of tSD over all vertices was analyzed by histograms, cumulative distributive function plots, 95th percentile, the percentage of tSD less than two coherence lengths of the imaging system (2ℓ = 12 μm), median, and mean of the tSD given by as the average of the tSD over all vertices in the eye n and layer l.

Finally, for comparison of our results with previous studies on repeatability of retinal layer thickness measurements, a single-value average thickness for the layer l at the time-point t of the eye n was computed as and the one-way random effect analysis of variance (ANOVA) intraclass correlation coefficients (ICCs) for the layer thickness was computed with at all time-points t = 1…T for all eyes n = 1…N.

The longitudinal variability of the BMO area and eccentricity measurement was assessed in each eye by the time-mean, minimum, and maximum, standard deviation, and coefficient of variation over all time points. BMO area and eccentricity were also assessed for repeatability by ICCs based on one-way random effect ANOVA. In addition to session-wise repeatability, BMO area, and peripapillary mean RNFL and choroidal thickness were checked for effects of the acquisition time of the day (diurnal change) and the time from baseline (longitudinal change over 3 weeks) by general regression models.

Results

BMO Dimension Repeatability

Tables 2 and 3 summarize the BMO area and eccentricity measurements over 3 weeks. In both tables, Subject 2 stands out with larger BMO areas and eccentricities, and greater standard deviations and coefficients of variation (COV). The ICC of the BMO area was 0.993, with 0.985 to 0.998 as the 95% confidence interval (CI). The ICC of the BMO eccentricity was 0.972, with 0.939 to 0.991 as the 95% CI. Overall, the BMO area and eccentricity measurements were stable.

Table 2

BMO Area Repeat Measure Statistics

Table 3

BMO Eccentricity Repeat Measure Statistics

Figure 4 plots BMO area by the imaging time of the day. In Figure 4a, data from all eyes are plotted. Each curve represents a single eye, and the points on the curve represent the measurements at different time points. Left and right eyes of the same subject can be distinguished as two curves with the same time points. In Figure 4b, the measurements were normalized by subtracting the average value from each eye, and plotted by subject. Intereye similarity can be observed although no consistent cross-sectional diurnal pattern is present. In agreement with the values presented in Table 1, S2 OD and OS show much larger BMO area measurements and fluctuation of the values throughout the day. Figure 5 plots BMO area by the number of days from baseline. In Figure 5a, data from all eyes are plotted. As in Figure 4a, each curve represents a single eye, and the left and right eyes of the same subject are shown as two curves with the same time points. In Figure 5b, the normalized BMO area measurements are plotted by subject. In the general regression analysis, there was no significant change in the BMO area over 3 weeks, and the time of the day did not affect BMO area.

Figure 4

Diurnal pattern of BMO area. (a) Each curve represents a single eye, and the points represent the measurements over 3 weeks acquired on different days and times of the day. Fellow eyes can be distinguished as two curves with the same time points. Subject 2 (light and dark green curves) shows significantly larger BMO area and measurement variability than other subjects. No consistent diurnal pattern is observed across the eyes. (b) BMO area measurements were normalized by subtracting the average value from each eye and plotted for each subject. Red curve represents the right eye (OD) and blue curve represents the left eye (OS). As in (a), no consistent diurnal trend is visible.

Figure 5

BMO area measurements over 3 weeks. (a) Each curve represents a single eye, and the points represent the measurements acquired on different days and times of the day. Fellow eyes can be distinguished as two curves with the same time points. Subject 2 (light and dark green curves) shows significantly larger BMO area and measurement variability than other subjects. No consistent longitudinal pattern is observed across the eyes. (b) BMO area measurements were normalized by subtracting the average value from each eye and plotted for each subject. Red curve represents the right eye (OD) and blue curve represents the left eye (OS). As in (a), no consistent longitudinal trend is visible.

RNFL Thickness Point-Wise Variability

The registration process yields point-wise correspondence between the baseline surface and the all follow-up surfaces, enabling point-wise statistics. In Figure 6, the left panel shows the 3-week vertex-wise time-average RNFL thickness over nine time points, color-mapped on the baseline RNFL surfaces for each eye. Even after averaging nine repeat images, the characteristic pattern and intereye similarities are well preserved, suggesting a high degree of spatial matching among the registered surfaces. The right panel in Figure 6 shows the 3-week tSD of RNFL thickness color-mapped on the baseline surfaces.

Figure 6

Left: Voxel-wise time-average of RNFL thickness from nine repeat measurements over 3 weeks, color-mapped on the baseline surfaces. The surfaces were cropped to annuli between 0.25 mm and 1 mm from the BMO. The averaged plots show the individual RNFL thickness patterns with high OD-OS similarity. Right: Vertex-wise tSD of RNFL thickness over the same measurements, color-mapped on the baseline surface as an indicator of measurement variability. S3 and S4 show higher variability than other eyes, patterned radially from BMO.

Distribution of the vertex-wise RNFL thickness tSD was visualized in aggregated histograms in Figure 7 and cumulative distribution function (CDF) plots in Figure 8. In both figures, the red vertical line marks two coherence lengths.

Figure 7

Aggregated histograms of vertex-wise RNFL thickness tSD over 3 weeks. The red vertical line marks two coherence lengths of the imaging device, 2ℓ = 12 μm. For all histograms, the majority of the vertices have tSD < 2ℓ.

Figure 8

Aggregated CDF plots of the vertex-wise tSD of RNFL thickness over 3 weeks. The red vertical line indicates two coherence lengths of the imaging device, 2ℓ = 12 μm. The blue horizontal line marks the 95th percentile. In four out of six subjects (8 out 12 eyes), the 95th percentile tSD is less than 2ℓ. For the rest, the 95th percentile tSD is between 2 and 4ℓ.

Figure 7 shows the majority of vertices with RNFL thickness tSD well below 2ℓ for all eyes. For Subjects 1, 2, 5, and 6, histograms peaked at less than 5 μm. For Subjects 3 and 4, histograms peaked between 5 and 10 μm. In Figure 8, the blue horizontal line marks the 95th percentile. For Subjects 1, 2, 5, and 6, more than 95% of the vertices have tSD less than 2ℓ. For Subjects 3 and 4, between 75% and 90% of the vertices have tSD less than 2ℓ.

The vertex-wise RNFL thickness tSD was further summarized over all vertices in Table 4 by the mean, median, and 95th percentile tSD values, and the percentage of vertices with tSD < 2ℓ. The mean and median tSD are less than 10 μm for all eyes, and tSD values at the 95th percentile are less than 12 μm for Subjects 1, 2, 5, and 6, and between 16 and 21 μm for Subjects 3 and 4. The percentage of vertices with tSD less than 2ℓ was over 96% for Subjects 1, 2, 5, and 6, and between 79% and 91% for Subjects 3 and 4.

Table 4

Point-Wise RNFL Thickness tSD

The ICCs of mean RNFL thickness was 0.997 with the 95% CI 0.993 to 0.999. In general regression analysis, mean RNFL thickness did not significantly change with the acquisition time of the day or the time from the baseline.

Choroidal Thickness Point-Wise Variability

In Figure 9, the left panel shows the 3-week vertex-wise time-mean choroidal thickness over nine time points, color-mapped on the baseline choroidal surfaces for each eye. The irregular thickness pattern is due to the choroidal blood vessels delineated at the posterior boundary of the layer. As in the time-average RNFL thickness maps, the shape features and intereye similarities are preserved after averaging of nine repeat measurements, suggesting a high degree of spatial matching among the registered surfaces. The right panel in Figure 9 shows the 3-week vertex-wise tSD of choroidal thickness color-mapped on the baseline surfaces.

Figure 9

Left: Voxel-wise time average of choroidal thickness from nine repeat measurements over 3 weeks, color-mapped on the baseline surfaces. The surfaces were cropped to annuli between 0.25 mm and 1 mm from the BMO. The averaged plots show the individual choroidal thickness patterns with high OD-OS similarity. Right: Vertex-wise tSD of choroidal thickness over the repeat measurements as an indicator of measurement variability. The distinct patterns reflect the segmentation of the posterior choroidal boundary as the outline of the posterior-most choroidal vessels, and indicate vascular variation as a source of variability in choroidal thickness measurement.

The distribution of the vertex-wise choroidal thickness tSD was visualized in aggregated histograms in Figure 10 and CDF plots in Figure 11. In both figures, the red vertical line marks two coherence lengths. Figure 10 shows majority of vertices with choroidal thickness tSD well below 2ℓ, and histogram peaks less than 10 μm. Compared to the histogram plots for RNFL thickness tSD in Figure 7, Figure 10 shows choroidal thickness tSD distributions are wider and more similar across the subjects, with all histogram peaks below 10 μm. In Figure 11, the blue horizontal line marks the 95th percentile. Generally larger tSD in choroidal thickness than in RNFL thickness is also apparent here, with all subjects except for Subject 6 showing less than 95% of vertices with tSD less than 2ℓ.

Figure 10

Aggregated histograms of vertex-wise choroidal thickness tSD over 3 weeks. The red vertical line marks two coherence lengths of the imaging device, 2ℓ = 12 μm. For all histograms, the majority of the vertices have tSD < 2ℓ.

Figure 11

Aggregated CDF plots of vertex-wise tSD of choroidal thickness over 3 weeks. The red vertical line indicates two coherence lengths of the imaging device, 2ℓ = 12 μm. The blue horizontal line marks the 95th percentile. Except for one subject (S6), all eyes show the 95th percentile tSD greater than 2ℓ but less than 4ℓ.

The vertex-wise choroidal thickness tSD was further summarized over all vertices in Table 5 by the mean, median, and 95th percentile tSD values, and the percentage of vertices with tSD < 2ℓ. The mean and median tSD are less than 2ℓ for all eyes. The tSD value at the 95% ranges between 13 and 16 μm except for Subject 3 OD (21.5 μm) and Subject 5 OS (23.9 μm). The percentage of vertices with tSD < 2ℓ ranges between 88% and 93%, except for Subject 3 OD (77.9 %) and Subject 5 OD (74.1 %).

Table 5

Point-Wise Choroidal Thickness tSD

The ICCs of mean choroidal thickness was 0.999 with the 95% CI 0.997 to 1.000. In general regression analysis, mean choroidal thickness did not significantly change with the acquisition time of the day or the time from the baseline.

The registration mean tSD for both RNFL and choroidal thickness before and after the registration steps was summarized over all vertices in Table 6. The tSD is lower in all postregistration cases.

Table 6

Comparison of RNFL and Choroidal Thickness Mean tSD Pre- and Postregistration

Discussion

This paper presents quantitative analysis of the variability of peripapillary measurements from longitudinal OCT images taken from healthy eyes over 3 weeks using nonrigid retinal surface registration based on our earlier published methods.[13,14] Twelve normal eyes of six subjects were imaged nine times over 3 weeks, and BMO dimensions, RNFL thickness, and choroidal thickness were measured. It should be noted that the subjects used in this analysis are all young and healthy. As such, the variability results may differ in the population that generally requires serial OCT images. The nonrigid surface registration used in this analysis establishes point-wise correspondence between the baseline and follow-up retinal surfaces, allowing for spatially detailed statistics. The results showed the range of measurement variability in the data due to physiological variations and errors in the processing pipeline, which will be useful to account for in our future longitudinal studies on tracking localized retinal layer changes over time. This framework is general, and shows that prior to longitudinal studies, assessing the “noise-floor” or the normal variability of any acquisition and measurement pipeline is an important task.

BMO, delineated in radial OCT frames as the termination points of the Bruch's membrane, has been shown as the more anatomically consistent landmark structure than the conventional optic disc measurements in stereo disc photography or other 2D projection images of the ONH.[21,22] In this study, the mean BMO area and mean BMO eccentricity ranged from 1.24 mm2 to 5.65 mm2 and 1.07 to 1.26, respectively. Subject 2 had significantly larger BMO area than other subjects, which may be due to the subject's high myopia and associated long axial length.[23,24] Subject 2 also had the largest BMO area variability. It is difficult to postulate the cause for this as Subject 2 measurement were extreme outliers and the rest of the data did not show a clear correlation between mean BMO area and BMO area variability (P = 192). A larger cohort spanning a range of BMO sizes including high myopes would be required for further analysis. In a previous repeatability study using a commercial Cirrus HD-OCT,[25] the optic disc area and eccentricity test–retest ICCs were 0.992 and 0.991, respectively. Our results were comparable within the 95% CI with ICCs of 0.993 for BMO area and 0.972 for BMO eccentricity.

Previous works on intervisit repeatability of mean RNFL thickness measurements reported ICC and COVs of 0.96% and 5.1%,[26] 0.97% and 2.7%,[27] and 0.95% and 1.13%.[28] With the RNFL thickness averaged at each time point, the ICC of our measurements was 0.997 with 0.993 – 0.999 as the 95% CI, and the COV was 1.58%, suggesting that the repeatability was comparable to that of the commercial machines. In the point-wise statistics, we found the mean RNFL thickness standard deviation was less than 10 μm for all eyes, with the majority of the 95th percentile tSD of less than three coherence lengths.

The averaged RNFL thickness maps (Fig. 6, left panel) show the distinct hourglass patterns. In the corresponding RNFL thickness tSD color maps (Fig. 6, right panel), vessel-like radial patterns were observed in some eyes. The same RNFL thickness tSD color maps were plotted over the baseline OCT enface images in Figure 12. The figure confirms the areas with high measurement variability overlapping with the locations of the retinal blood vessels (RVs). Several recent studies investigated the effect of RVs in RNFL measurements: RVs contributed to RNFL thickness measurements[29,30] with the local thickness maxima correlating with the locations of major RVs,[31] and RV distribution was shown to contribute to interindividual variability of circumpapillary RNFL thickness measurements among healthy volunteers.[32,30] However, in another work,[33] adjustment for major blood vessels did not decrease interindividual variability of circumpapillary RNFL thickness, and the location of major blood vessels did not affect the false-positive rate in RNFL color codes in Cirrus OCT analysis.[34] The spatial correspondence of the retinal vessels with nerve fibers also demonstrates the blood vessel development influenced by axonal density.[35-37]

Figure 12

RNFL thickness tSD color maps overlaid with baseline enface OCT images. The figure shows the regions with high RNFL thickness measurement variability correspond to the locations of RVs.

Figure 12 shows that in our data RVs contributed to intraindividual measurement variability. Pereira et al.[30] reported superior and inferior regions as the regions with the highest impact of retinal vessel density on circumpapillary RNFL thickness, which agrees with our result. However, previous studies on reproducibility of RNFL thickness measurements using Stratus OCT[38] and Cirrus OCT[27] consistently found the nasal sector to have the highest measurement variability, suggesting sources of variations other than large RVs in the superior and inferior regions. In a recent study, Radcliffe et al.[39] reported positional shift of RVs associated with functionally progressive glaucoma, neuroretinal rim loss, and disc hemorrhage. The registration pipeline in this study may be utilized to detect such localized changes by mapping RV positions directly on the registered surfaces.

The choroidal thickness measurements exhibited overall greater variability than the RNFL thickness measurements, as shown in Figures 10 and 11 and Table 4. This is expected considering the choroid is a vascular layer, and the vessels in RNFL were shown to contribute to the RNFL thickness measurement variability as discussed above. The choroidal vasculature is structurally more complex than the larger retinal vessels, and the intensity fall-off at the imaging depth, limitation in the lateral resolution, shadowing, and light scattering by the retina make it difficult to clearly visualize the choroidal vasculature in order to compare the layer enface image and measurement variability as for RNFL in Figure 12. The pattern in the tSD maps in Figure 9, similar to the outline of large blood vessels at the posterior border of the choroid, suggests that the variability is affected by the vasculature. However, in the same figure, there is no visible connection between the choroidal thickness (left panel) and measurement variability (right panel).

Previous studies on choroidal thickness measurement repeatability were mainly concerned with the difficulty in delineation of the posterior boundary of the choroid by manual raters.[40-42] The SS system with 1060-nm wavelength source used in this study allows for better imaging of deeper structures, and the automated segmentation algorithm was evaluated in our previous studies for choroid segmentation.[43,44] A recent paper on a prototype SS-OCT system from Topcon with 1050-nm wavelength source and automated layer segmentation and thickness measurement software[45] reported ICCs of 0.99 with the 95% CI 0.98 to 0.99 for average optic disc choroidal thickness in three consecutive measurements. The ICCs for averaged choroidal thickness in our data was 0.99 with the 95% CI 0.997 to 1.000. We note that the specification of the posterior boundary of the choroid can vary depending on the rater or system in different studies. In the Topcon prototype, the posterior boundary of choroid was delineated as a smooth curve encompassing the blood vessels. In our measurements, the posterior boundary was delineated as the outline of the most posterior vessels. In the point-wise statistics, the mean choroidal thickness SD was less than 11 μm, and the majority of eyes showed 95th percentile SD less than three coherence lengths with two exceptions, S3 OD and S5 OD.

Tan et al.[46] observed a significant diurnal variation in manually delineated foveal choroidal thickness with a mean decrease of 33.7 μm. In the general regression analysis on our sample, the acquisition time of the day did not significantly affect the measurement. Possible reasons for this difference include that our choroidal thickness was averaged over the peripapillary region instead of being measured at a single foveal location, and the study was not specifically designed for observing diurnal patterns (the subjects were imaged at arbitrary times on different days). A more deliberate data acquisition followed by the registration process may provide a detailed visualization and statistics of diurnal variation in the choroid. A comparative study of the choroidal diurnal patterns between the foveal and peripapillary regions can be also conducted.

The mean tSD of RNFL thickness was correlated to that of choroidal thickness in mixed model analysis (P = 0.029). This correlation was driven by S3 OD and S4 OD/OS, which showed much larger vessel-oriented RNFL tSD than other eyes (Fig. 6, right panel). Excluding these three eyes, the correlation was not significant. Also, S5 OD had both the largest mean choroidal thickness tSD and smallest mean RNFL thickness tSD. This suggests that any factor that affects both layers, or the image volume as a whole, such as uncorrected global misalignment, rotation, or motion artefact, was not a significant part of the thickness variability. The result also implies that the correlation between the vessel-induced fluctuations in RNFL thickness and those in choroidal thickness may not be simple and direct, and in comparing RNFL and choroid, for example in diurnal patterns, it would be beneficial to separate RNFL to vascular and nonvascular regions for more accurate analysis.

Each stage in the processing pipeline (Fig. 1) may add synthetic variability to the inherent biological variations. In the image acquisition stage, lateral motion and aberrations in the images can be potential sources of measurement variability. Additionally, the image quality in terms of signal-to-noise ratio can also have an effect on each subsequent stage in the pipeline. We also observed that “edge-clipping” in few surfaces contributed to their measurement variability. Figure 13 shows the nine RNFL surfaces and thickness mapping for Subject 2 OS (first row), registration of the eight follow-up surfaces to the baseline surface (second row), and the difference maps between each follow-up and the baseline RNFL (third and fourth rows). In the first row, we note that at t3, t5, and t8, the temporal (left) side of the surfaces are clipped, being out of the image frame. The surface registration matches these incomplete surfaces to the original surface, and in the third row there are artefact RNFL thickness differences along the left edge of the surfaces at these time points. The effect of this can be seen in the RNFL thickness SD map of S2 OS in Figure 6, and in the choroidal thickness SD map of the same eye in Figure 9. Such artefacts can be avoided by better centering the image in reference to the BMO, removing the clipped surfaces from analysis, and using a larger field of view for those with larger BMO like Subject 2 with high myopia.

Figure 13

RNFL surface registration in Subject 2 OS. Row 1: RNFL thickness maps on the original RNFL surfaces acquired at nine time points in 3 weeks. At t3, t5, and t8, the RNFL surfaces are clipped at the temporal (left) side. Row 2: Follow-up RNFLs were registered to the baseline. The clipped surfaces are artificially matched to the unclipped baseline surface. Row 3: Vertex-wise RNFL thickness difference between the baseline and each follow-up RNFL. Row 4: Difference maps in Row 3 are shown thresholded by the axial coherence length of the system ℓ = 6 μm. In both Rows 3 and 4, the clipped regions show large difference with the baseline RNFL.

Segmentation errors are additional sources of measurement variability. Human error in the manual segmentation may have affected the BMO measurements. To lessen the effect of outliers, an ellipse was fitted to the segmented BMO points, which may have introduced error of its own. The automated segmentation technique in this paper was shown to have high agreement with a trained manual rater.[43] All automated segmentation results in this study were reviewed and corrected by a trained manual rater. This mitigates the effect of large machine-errors in the segmentation, but it may introduce additional human errors. Alonso-Caneiro et al.[47] demonstrated variations in axial retinal layer thickness measurements due to tilt or curvature in OCT images. In the current paper, the thickness was measured as the closest distance between the posterior and anterior surfaces of the layer, which, as Alonso-Caneiro and the coauthors note, may not have a unique mapping between the surfaces in certain cases. Although the closest distance method is likely more robust to the errors experienced by the axial measurement method due to image tilt or curvature, and the issue of unique mapping is likely negligible where the retinal layer is smooth and planar, it would be worthwhile to investigate if alternative methods of thickness measurement, such as Laplacian thickness, would affect the measurement variability.

The registration step occurred in rigid (affine transformation) and nonrigid (diffeomorphisms) stages. The results of the intensity-based rigid registration were qualitatively examined before being passed on to the nonrigid registration module. Errors in the rigid registration occurred mainly due to intensity variations and motion artefacts in the images. These were corrected by manually specifying registration landmarks. The nonrigid registration method was validated in our previous work[14] using clinical optic nerve head OCT images and artificial phantom data, which showed errors less than the coherence length at all vertices.

Conclusion

The study demonstrated spatial matching achieved by a registration-based pipeline among repeated longitudinal RNFL and choroidal surfaces, and assessed the measurement variability in normal subjects. Starting with the working hypothesis that there is no variability in measurements such as RNFL/choroidal thickness made from longitudinal repeated scans, the data collected showed that the thickness measurements were overall highly repeatable considering the image coherence length given that the measurements were collected and compared point-wise and not spatially averaged. However, larger variability due to blood vessels and subject-wise variability were observed, which showed that sources of variability do exist, and are spatially and individual dependent. In the future longitudinal studies, it may be beneficial to obtain a set of repeat measures at the onset of study for each subject, and establish individual “baselines” of measurement variability relative to which longitudinal changes due to disease or treatment progression can be ascertained.