December 2010
Volume 51, Issue 12
Free
Glaucoma  |   December 2010
An In Silico Model of Scanning Laser Tomography Image Series: An Alternative Benchmark for the Specificity of Progression Algorithms
Author Affiliations & Notes
  • Neil O'Leary
    From the Department of Optometry and Visual Science, City University London, London, United Kingdom; and
  • David P. Crabb
    From the Department of Optometry and Visual Science, City University London, London, United Kingdom; and
  • David F. Garway-Heath
    the National Institutes of Health Research (NIHR) Biomedical Research Centre for Ophthalmology, Moorfields Eye Hospital NHS (National Health Service) Foundation Trust and UCL Institute of Ophthalmology, London, United Kingdom.
  • Corresponding author: David P. Crabb, Department of Optometry and Visual Science, City University London, Northampton Square, London EC1V 0HB, UK; david.crabb.1@city.ac.uk
Investigative Ophthalmology & Visual Science December 2010, Vol.51, 6472-6482. doi:10.1167/iovs.10-5355
  • Views
  • PDF
  • Share
  • Tools
    • Alerts
      ×
      This feature is available to Subscribers Only
      Sign In or Create an Account ×
    • Get Citation

      Neil O'Leary, David P. Crabb, David F. Garway-Heath; An In Silico Model of Scanning Laser Tomography Image Series: An Alternative Benchmark for the Specificity of Progression Algorithms. Invest. Ophthalmol. Vis. Sci. 2010;51(12):6472-6482. doi: 10.1167/iovs.10-5355.

      Download citation file:


      © ARVO (1962-2015); The Authors (2016-present)

      ×
  • Supplements
Abstract

Purpose.: There is no gold-standard measurement of glaucomatous structural progression against which to validate software progression algorithms. A computer model was developed and validated to simulate stable series of Heidelberg Retina Tomograph II (HRT; Heidelberg Engineering, Heidelberg, Germany) images, with realistic topographic variability, suitable for benchmarking false-positive rates of progression algorithms.

Methods.: Three confocal image stacks were selected from each of five sets of HRT II scans, obtained within 6 weeks in 127 eyes of 66 patients. For each eye, a simulated series was propagated from one baseline confocal stack by adding fixational eye movements, photon-counting, and electronic measurement noise. Simulated confocal stacks were imported into the HRT software to generate topography images. Real and simulated image comparisons were quantified with the mean pixel height standard deviation (MPHSD), image cross-correlation (CC) of pixel-wise variability maps, and the rim area (RA) coefficient of variation (CV).

Results.: The mean difference (95% limits of agreement; LoA) in MPHSD between real and simulated images was 3.5 μm (−20.9 to 28.8 μm) within mean topographies and 2.0 μm (−5.4 to 9.3 μm) between mean topographies. The mean CC between real and simulated spatial variability maps was 0.58 within mean topographies and 0.54 between mean topographies. The mean difference (95% LoA) between real and simulated mean topography RA CV was −2.1% (−17.6% to +13.4%). Variability about anatomic features was well reproduced.

Conclusions.: Simulation realistically reproduces variability in real, stable images acquired over a short period. Stability in clinical datasets is uncertain, whereas in these modeled series, it is certain. This method provides benchmark datasets on which the specificity of progression algorithms can be tested.

Confocal scanning laser tomograph (CSLT) technology, typified by the commercially available Heidelberg Retina Tomograph II (HRT, Heidelberg Engineering GmbH, Heidelberg, Germany), gives three-dimensional images of the optic disc and peripapillary retina. This established imaging modality provides reproducible results 1 and is widely used in the assessment of structural damage in the glaucomatous optic disc. The principal goals of imaging with application to glaucoma are to assist the user in discriminating between normal and glaucomatous optic discs and to identify and track disease progression. 
Although the HRT II discriminates between normal discs and those with glaucoma reasonably well, 2,3 the real promise of CSLT technology is in detection of change in disc structure over time, offering the clinician another tool for glaucoma management. The repeatability of HRT measurements has been quantified and used to derive limits beyond which change cannot be accounted for by measurement variability. 4 Age-related change has been quantified to help differentiate age effects from disease progression. 5 Statistical techniques have been developed to detect progression based on population 6,7 or individual subject 8 10 variability. However, there is little evidence to suggest that one method of detecting change is better than another, and the clinician is left to wonder what method of detecting change is best to use. 
Most research assessing how well these techniques can identify and track glaucomatous progression with the HRT has focused on examining agreement between structural and functional measures of progression which In practice, involves assessing the diagnostic performance of HRT information in predicting visual field (VF) changes. 11 14 VF changes, when used as a reference standard in this way, can confound results by dissociating factors between structural and functional changes. 15 As a consequence, the relative proportions of associated and independent behavior are not well known, and the temporal sequence of structural and functional glaucomatous change cannot be well-defined. 
An alternative reference standard against which to assess the performance of CSLT imaging is change identified by glaucoma experts from optic disc or retinal nerve fiber layer (RNFL) photographs. Stereometric parameters have poor to moderate diagnostic performance as predictors for expert assessed changes in RNFL photographs. 16 Expert, but subjective, clinician assessment of ONH change using the HRT software progression reports has poor agreement with assessment of series of stereoscopic optic disc photographs. 17 Depending on the progression criteria applied and the expert observer, only moderate concordance between the HRT Topographic Change Analysis 9 (TCA) and observed change in monoscopic (81% 18 and 44%–71% 19 ) and stereoscopic (65% 20 ) optic disc photographs has been reported. It has also been shown that using stereophotograph-assessed change as the reference standard does not help to determine which HRT progression-detection algorithm, of a variety of methods, best identifies glaucomatous change. 21 It may be that HRT and stereophotograph assessment are identifying different aspects of structural progression. 
With stereoscopic disc photograph or visual fields as the reference standard for progression/stability, the apparent poor discrimination between stable and progressing eyes by HRT progression algorithms may be due to either poor performance of the algorithms or inaccuracies in the reference standard itself. Poor between-expert agreement in identifying progression based on stereoscopic photographs 22 24 illustrates the shortcomings of stereophotography as a reference standard. For this reason, the concept of a simulation mimicking stable serial images is a hugely useful one. Progression algorithms can be applied to the simulated stable series and criteria selected to obtain low false-positive rates and the false-positive rates of different progression algorithms can be compared with each other, knowing that the image series are truly stable. Realistic and sophisticated simulation is common in other areas of medical imaging—in particular, neuroimaging. 25,26 In these fields, it is widely considered the best way to validate automated postprocessing algorithms, such as assessing the performance of the alignment of nonrigid images and volumes over time. 27,28 However, with regard to the detection of glaucomatous progression, the practice of realistically simulating longitudinal image series has been neglected, with some exceptions, in scanning laser polarimetry (SLP) 29 and elementary work on the HRT (O'Leary N, et al. IOVS 2009;50:ARVO E-Abstract 2249), 8 and the predominate focus has been on simulation of VF series. 30 34  
In this study, we developed a virtual platform in which components of the optic disc topographic image formation process of the HRT are simulated computationally, taking into account features relevant to topographic variability, such as eye movements and other sources of noise during and between HRT examination. This goal is approached by using empiric and theoretical parameters of these noise sources to propagate HRT series from real baseline scans. The resulting simulated HRT images will be analogous to the real HRT images obtained in clinics with variability similar to that inherent in real, stable series. The key advantage of this approach is having exact knowledge that the morphology of the original structure is not changing, because all images are derived from the same baseline scan. This simulation can then provide a test bed for comparing false-positive rates of HRT change detection methods. 
Methods
Data
HRT II scans acquired for a study on test–retest measurement variability were used; the data are fully reported on elsewhere. 35 In summary, they consisted of the images of both eyes of 74 patients (43 ocular hypertensive [OHT], 31 primary open-angle glaucoma [POAG]) who underwent five HRT II examinations within 6 weeks. OHT was defined as intraocular pressure (IOP) more than 21 mm Hg on two or more occasions and a baseline VF test (Humphrey 24-2 full-threshold; Carl Zeiss Meditec, Inc., Dublin, CA) with an Advanced Glaucoma Intervention Study (AGIS) score of 0. POAG was defined as a consistent AGIS VF score greater than 0, and a pretreatment IOP greater than 21 mm Hg on two or more occasions. All HRT II data were processed in the latest HRT III software viewer version (ONH Viewer Module 3.1.2.15; Heidelberg Engineering). All subjects had previous experience with scanning laser tomography. This study adhered to the tenets of the Declaration of Helsinki and had local ethics committee approval and the subjects' informed consent. 
From 148 eyes of 74 patients, 127 eyes of 66 patients were selected on the basis of available raw data (five computable mean topographies and a baseline scan with original confocal image stacks). 
Simulation Development
The HRT II acquires a confocal image stack composed of optical sections taken through a sequence of focal planes 62.5 μm apart in a horizontal scanning angle of 17.5° and vertical scanning angle of 15° around the ONH surface. These are aligned to compensate for eye movements occurring during scan acquisition. Each optical section image is sampled as an array of 384 rastered horizontal line scans, each sampling 448 pixels. The intensity of light reflected in each optical section for a location on the ONH/retinal surface is used to calculate the surface height at that location (approximately the depth position of the optical section with the highest reflectance intensity at that location). The matrix of surface height values of the central (384 × 384) sampled locations forms the topography image with the information from the 64 horizontally peripheral pixels used only before topography formation. 
The simulation uses the unprocessed (unaligned) confocal image stack as the basis for simulation. For purposes of simplicity and clarity, a scan triplet is used to refer to the three single image stacks (referred to as single scans) obtained during each acquisition or examination. The patient is repositioned at the HRT between image acquisitions. It is from this scan triplet that three single topographies are derived; these are aligned and averaged to form the mean topography
One baseline single scan for each eye is chosen as the seed for the generation of the simulated series. Each optical section is aligned with the section above and below it by using a cross-correlation algorithm, to reduce eye and head movements already potentially present in the seed data. Within each seed stack, various sources of noise, present during image acquisition, that potentially result in variability between the generated topography (single and mean) images is added in the following order: (1) within-examination eye movement, (2) within-examination head movement, (3) between-examination head movement, and (4) device noise, composed of quantum effect noise (shot/photon counting) and electronic noise (Johnson-Nyquist/thermal). Simulated single topographies are formed from resultant image stacks using the HRT III software; the simulated mean topographies are calculated from the simulated single topographies. The input and sequence of the addition of these sources of noise are shown in Figure 1
Figure 1.
 
Schematic of simulation: HRT single and mean topography formation formed from a single scan volume of CSLO optical sections. Processes are represented by gray boxes with intermediate data states indicated by white boxes and initial and final data states by rounded boxes.
Figure 1.
 
Schematic of simulation: HRT single and mean topography formation formed from a single scan volume of CSLO optical sections. Processes are represented by gray boxes with intermediate data states indicated by white boxes and initial and final data states by rounded boxes.
Within-Examination Eye Movements.
In HRT image acquisition, the patient is asked to fixate on a target before and during the acquisition. Thus, the eye movement model consists of the three main components of normal fixational eye movements: ocular microtremor, ocular drift, and microsaccades. Ocular microtremor is assumed to have an amplitude range of 0 to 2 arcmin (median, 17.5 arcsec), 36 with a frequency power spectrum from 0 to 150 Hz. 37 Ocular drift is assumed to have linear speed in an amplitude range of 3 to 12 arcsec, with durations spanning the intervals between microsaccades. 38 Microsaccades are assumed to have a mean amplitude of 30 arcmin and a mean period of 30 ms. 39 Mean speeds of microsaccades are calculated from the main sequence relationship 40 and are parametrically related to the amplitudes of microsaccades. Each microsaccade speed is nonlinear and is approximated to that of a half-Gaussian function. The mean frequency of microsaccades is assumed to be 1.5 per second. 41  
Certain assumed rules are enforced in the eye movement model: drift and microtremor take place simultaneously; neither drift nor microtremor can occur during microsaccades; no subsequent ocular drift can be in the same direction (±90°) as the previous microsaccade. This approximates the directional frequencies of eye movements observed in previous studies. 39 No initial directional preference was chosen for eye movements, although results have shown that horizontal and vertical eye movements display different characteristics. 42  
From these parameters and rules, a retinal trace (x, y, t; the path of the retina in transverse axes over time) is derived for the position of the retina. The sampling period of each optical section is approximately 24 ms (frequency, 42 Hz), and within each section, the rastered horizontal line-sampling period is approximately 62 μs (frequency, 16 kHz). As the latter frequency is many orders of magnitude above maximum frequencies recorded for eye movements, the horizontal line scan time was considered a sufficient sampling frequency for our eye movement model (i.e., movement occurs between line scans and between optical sections). The transverse (x and y) spatial sampling for each HRT pixel is approximately 2 arcmin in an emmetropic eye. 
Within-Examination Head Movements.
These are represented by scaling (magnification) and the translational and rotational movements of the retina to correspond with movements of the head. There is little research on head movement parameters during positioning on chin and head rests, but studies have shown that a mean of approximately 2.5% of head motion remains uncorrected by translational compensatory eye movements. 43 As a consequence, only translations along the z-axis (in depth) of scanning are included in the head movement model. This model is implemented by applying cumulative, random, axial translation (with appropriate scaling) and rotation to each sequentially scanned optical section throughout the scan triplet. 
Between-Examination Movement.
The global uniform transformation and change in overall reflectivity of the confocal stack mimics variability resulting from placement of the head anew on the head and chin rests. A standard deviation (SD) of 10 pixels is chosen (by assumption) for the normal distribution from which lateral translations of the follow-up stack relative to baseline can be chosen as an estimate for the operator uncertainty in centering the ONH on scanning. Changes in the autofocus settings of the machine (correct axial placement of central CSLO section) are modeled by random small scaling and depth translations. Differences in initial head angle placement in each examination are mimicked by rotations of the scan triplet about each axis with relation to the baseline. 
For all eye and head movements, implementation of deriving new intensity values at each new simulated pixel involves the defining of new coordinates for every pixel in the CLSO section stack under a given transformation. Filling in of unknown data between known intensity values at fixed pixel coordinates requires interpolation. Bilinear or trilinear interpolation was chosen for this process. It applies products of linear functions to the data of the neighboring 4 or 8 pixels in two- or three-dimensional space to obtain intensity values of new coordinates that do not coincide exactly with the coordinates of any pixel in the seed stack. If the transformation involves only transverse components, then bilinear interpolation is performed using only the values of the new coordinate's transverse neighbors. If axial components are involved in the transformation, then trilinear interpolation is performed using the axial neighboring pixel information in addition. 
In the rare cases in which there is no valid measurement from the reflectance signal from the focused spot, these areas remain invalid under the new coordinate system. Bilinear/trilinear interpolation will render neighboring values that are close enough also invalid under the coordinate transform. The effect of movement that causes the resampled pixels to fall outside the boundaries of the original optical section, known as windowing, is an important technical point to consider. Null values are used and the topography formation and alignment algorithms will ignore information at these pixels. This effect occurs mostly to the edges of the optical section images. As previously stated, there are 64 extra horizontally peripheral pixels available for use in each optical section before topography formation. The simulation uses the information of these extra pixels to remove the effects of windowing from horizontal axis movement. 
As the dwell time of the scanning beam at each sampled point is <0.5 × 10−6 seconds, the effect of optical blur or smearing of the optical point spread function (PSF) due to eye movements is considered negligible for this model. Even at the maximum simulated eye speed, the maximum displacement during this dwell time would be <0.00025 of the transverse spatial sampling interval between pixels (one pixel width). One pixel width is approximately equal to a typical full width at half maximum of the PSF in the HRT scanning system, as referred to in the imaging system's lateral resolution specification. 44  
Device Noise.
The operating principles of the HRT are such that retinal reflected light is gathered point by point and focused onto a detector. A single spot detector—an avalanche photodiode (APD)—is used for sequentially measuring the incident reflectance intensities. The principal sources of device noise in HRT image acquisition are assumed to be intrinsic noise 45 in the form of quantum or photon-counting (shot) noise, 46 along with electronic noise, which is primarily thermal noise (Johnson/Nyquist) and APD noise. 47 Photon-counting noise can be modeled as a Poisson process, where the intensity measured at the detector (proportional to number of photons striking the detector) for a sampled point can be chosen from a Poisson distribution with a mean intensity of λ. In this simulation, λ is given by the intensity of the pixel of the seed optical section after eye and head movements have been added. Thermal and APD noise which are intrinsic noise from the properties of electronic circuitry in general and semiconductor diodes, respectively, can be combined and approximated to white noise with a given amplitude. Each pixel in each optical section has added spatially independent noise sampled from a normal distribution with a mean of 0 and with the amplitude of the noise determined by the SD of the distribution. In this simulation, the SD is estimated at 2% of the maximum 256 stored intensity gray-level values in an 8-bit analog-to-digital converter. After noise is introduced, each pixel intensity value is rounded to the nearest 8-bit integer, to represent the digitization of the signal in the HRT optical section. 
Noise Follow-up Time Dependence.
The difference between one scan triplet and another is simply modeled by the change in initial head and eye position corresponding to a repositioning of the head onto the head and chin rest and a resetting of the acquisition centration on the ONH. The within-scan triplet noise parameters do not change with time, nor do the differences between examinations. These assumptions are supported by the finding of no difference in variability between examinations, with a range of intervals too short to be age-related (1 hour to 4 weeks). 48  
Implementation.
This simulation model was implemented in commercial software (MatLab, 2007a; The MathWorks Inc. Natick, MA). It contains the full technical specifications and is freely available on request from the authors. As indicated in Figure 1, all simulated optical stack data are imported into the HRT III software for topography calculation and alignment. This processing is identical with that of real data. 
Simulation Testing
The appropriateness of the simulated images was assessed quantitatively by comparing the following measures between the real and simulated topographies for each eye. 
Local Within-Examination Variability.
This is a map of within-examination pixel height standard deviations (PHSDw). Each pixel of this map, PHSDw(x,y), is the sample SD at that pixel across the three aligned single topographies comprising a mean topography. 
Local Between-Examination Variability.
This is a map of between-examination pixel height standard deviations (PHSDb). Each pixel of this map, PHSDb(x,y), is the sample SD at that pixel across all mean topographies in the series. 
Global Within-Examination Variability.
This is the mean within-examination pixel height standard deviation (MPHSDw). It is used as an overall measure of within-examination quality and is provided in the HRT III software. It is the geometric mean of all pixels in the PHSDw map of a mean topography. 
Global Between-Examination Variability.
This is the mean of between-examination pixel height standard deviation (MPHSDb). It was used in the study as a global measure of the variability between mean topographies in a series. It is the geometric mean of all pixels in the PHSDb map for each series. 
Clinical Measure of Variability: Neuroretinal Rim Area (RA) Variability of Mean Topographies.
Neuroretinal RA gives a numerical value to a clinically recognizable and meaningful feature. The loss of RA tissue provides a surrogate measure of the loss of retinal ganglion cell axons typical of glaucomatous damage. 49 It is the least variable of the stereometric HRT parameters and has been well characterized. 4,50 RA is thus a useful indicator for detecting glaucomatous progression. 51,52 In this study, the Moorfields fixed, standard reference plane (RP) was used to calculate RA. 53  
The RA coefficient of variation (CV) is calculated as the ratio of the SD to the mean of the RA across mean topographies in a series. The RA CV is compared between real and corresponding simulated mean topography series, to examine the performance of the simulation in mimicking the variability of a clinically meaningful parameter. 
Single-value metrics of variability (MPHSDw, MPHSDb, and RA CV) were compared by plotting the paired values (or the mean when appropriate) of simulated and real series and examining the mean difference. Mean differences (bias) and the 95% LoA between real and simulated series were calculated for these metrics. Analysis of proportional bias and agreement was also performed. 
Correlation in the Spatial Distribution of Topographic Variability between Real and Simulated Series.
The maximum value of the normalized cross-correlation coefficient (NCC), a standard image-processing metric for similarity of two images, 54 was used to compare the local variability measures (PHSDw and PHSDb maps) of real and simulated topographies. The term maximum is used, as PHSDw and PHSDb maps can be shifted along the x and y coordinates (with an upper limit of 20 pixels) to find the maximum NCC. This approach is logical, as the real and simulated series will have different baselines, and, as such, PHSD maps must be approximately aligned to each other. This metric captures the extent of mutual spatial information in the real and simulated variability maps and thus how well the location of topographic variability is reproduced by the simulation. NCC ranges from a minimum possible value of −1, representing perfectly anticorrelated image data, to a value of +1, representing perfectly correlated image data, where a value of 0 represents no correlation whatsoever. 
NCC values from comparing real and simulated mean topography PHSDw variability maps were benchmarked against NCC values from comparisons of the real mean topography PHSDw variability maps with each other. The NCC for a real–simulated series comparison was the mean of the 25 possible pairs of comparisons. The NCC for the real–real series comparison was the mean of 10 possible pairs of comparisons. 
Results
Of the selected 127 eyes of 66 patients, 63 were left eyes and 64 were right eyes. The median patient age at baseline examination was 70 years (interquartile range, 63–75 years). 
Mean differences (bias) and 95% LoA between real and simulated data for MPHSDw, MPHSDb, RA CV, and NCCs for PHSDw maps are displayed in Table 1. The results of analysis of proportional effects on agreement between real and simulated series are also included in Table 1, with statistical significance of linear relationships between biases and average measurements and statistical significance of linear relationships between 95% LoA and average measurements. In addition, a log transformation was used on the mean values and before differences were calculated between real and simulated measures—an approach designed to ensure robustness against the effect of heteroscedasticity and potential outliers on any measured proportional effects. 55 These characteristics are clear in the Bland-Altman plots of Figures 2 and 5. Accordingly, log-transformed significance of proportionally dependent bias and 95% LoA are included in Table 1
Table 1.
 
Measures of Difference between the Simulated and Real Image Series for Global, Regional, and Local Measures of Variability
Table 1.
 
Measures of Difference between the Simulated and Real Image Series for Global, Regional, and Local Measures of Variability
Measures of Variability Measures of Difference
Uniform Significance of Proportional Dependence
Linear Log-Transformed
Bias 95% LoA Bias 95% LoA Bias 95% LoA
MPHSDw 3.5 μm −20.9 to 28.8 μm P = 0.10 P < 0.01 P = 0.61 P < 0.01
MPHSDb 2.0 μm −5.4 to 9.3 μm P < 0.01 P < 0.01 P = 0.35 P = 0.27
RA CV −2.1% −17.6 to −13.4% P = 0.09 P = 0.08 P = 0.40 P = 0.07
NCC PHSDw 0.052 0.039 to 0.065 P < 0.01 P = 0.25 P < 0.01 P = 0.87
The Bland-Altman plots in Figure 2 show the agreement between paired real and simulated eye series for MPHSDw and MPHSDb, respectively. Agreement between paired real and simulated topography series for RA CV values for global RA is shown in the Bland-Altman plot in Figure 3. Linear proportional biases and 95% LoA found to have statistical significance in either linear or log-transformed domains (P < 0.05) are displayed on all Bland-Altman plots. 
Figure 2.
 
Bland-Altman plots showing series-wise agreements between (a) average within-examination days (MPHSDw) for real and simulated mean topography data and (b) agreement of between-examination MPHSDb for real and simulated mean topography data. The mean difference (bias) of average MPHSDw was 3.5 μm (95% LoA, −20.9 to 28.8 μm). The mean difference of MPHSDb was 2.0 μm (95% LoA, −5.4 to 9.3 μm). Uniform 95% LoA illustrate only approximate limits of agreement, as the heteroscedasticity of the data is apparent. Dashed lines: statistically significant linear proportional bias for MPHSDb and the significant linear proportional increase of the widths of the 95% LoA for MPHSDw and MPHSDb.
Figure 2.
 
Bland-Altman plots showing series-wise agreements between (a) average within-examination days (MPHSDw) for real and simulated mean topography data and (b) agreement of between-examination MPHSDb for real and simulated mean topography data. The mean difference (bias) of average MPHSDw was 3.5 μm (95% LoA, −20.9 to 28.8 μm). The mean difference of MPHSDb was 2.0 μm (95% LoA, −5.4 to 9.3 μm). Uniform 95% LoA illustrate only approximate limits of agreement, as the heteroscedasticity of the data is apparent. Dashed lines: statistically significant linear proportional bias for MPHSDb and the significant linear proportional increase of the widths of the 95% LoA for MPHSDw and MPHSDb.
Figure 3.
 
Bland-Altman plots showing series-wise agreement between RA CV for real and simulated data. Mean difference between values for real and simulated data were −2.1% (95% LoA: −17.6% to 13.4%). No proportional effects were found on the bias or LoA.
Figure 3.
 
Bland-Altman plots showing series-wise agreement between RA CV for real and simulated data. Mean difference between values for real and simulated data were −2.1% (95% LoA: −17.6% to 13.4%). No proportional effects were found on the bias or LoA.
The local variability correlation between real and simulated data is shown in Figure 4. The average, series-wise, maximum, NCC of PHSDw distribution between simulated and real data is shown in Figure 4a, with a mean of 0.58 (SD 0.12). The distribution of series-wise, maximum NCC of PHSDb maps between simulated and real data is shown in Figure 4b, with a mean of 0.54 (SD 0.10). The distribution of average, series-wise, maximum NCC between real PHSDw maps and other real PHSDw maps in the same series is shown in Figure 4c, with a mean of 0.64 (SD 0.14). 
Figure 4.
 
Distributions of maximum NCCs of pixel SD maps between real and simulated data (a) within examination (PHSDw) averaged over all pair-wise comparisons and (b) between-examination PHSDb. (c) The series average NCCs of all pair-wise pixel SD maps between mean topographies in the same series. The extremes are interpreted as follows: −1, perfectly negatively correlated; 0, uncorrelated; 1, perfectly positively correlated.
Figure 4.
 
Distributions of maximum NCCs of pixel SD maps between real and simulated data (a) within examination (PHSDw) averaged over all pair-wise comparisons and (b) between-examination PHSDb. (c) The series average NCCs of all pair-wise pixel SD maps between mean topographies in the same series. The extremes are interpreted as follows: −1, perfectly negatively correlated; 0, uncorrelated; 1, perfectly positively correlated.
Agreement of average NCC of real–real PHSDw maps against the average NCC of real–simulated PHSDw maps is shown by the Bland-Altman plot in Figure 5. The mean difference was 0.052 and the 95% LoA were 0.039 and 0.065. 
Figure 5.
 
Benchmarking of PHSDw NCC values between real and simulated data: Bland-Altman plot showing agreement of average NCC for real–real PHSDw map comparisons and average NCC for real–simulated PHSDw map comparisons. Solid lines: uniform bias of 0.052 and 95% LoA of 0.039 to 0.065; dashed line: statistically significant linear proportional bias.
Figure 5.
 
Benchmarking of PHSDw NCC values between real and simulated data: Bland-Altman plot showing agreement of average NCC for real–real PHSDw map comparisons and average NCC for real–simulated PHSDw map comparisons. Solid lines: uniform bias of 0.052 and 95% LoA of 0.039 to 0.065; dashed line: statistically significant linear proportional bias.
Examples of mean reflectance images, mean topographies, and corresponding PHSDw maps for real and simulated data are shown in Figures 6 and 7. Examples of series average mean reflectance images, series average mean topographies, and corresponding PHSDb maps for real and simulated data are shown in Figures 8 and 9
Figure 6.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 20 μm). (d) Corresponding simulated mean topography (MPHSDw, 22 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.55.
Figure 6.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 20 μm). (d) Corresponding simulated mean topography (MPHSDw, 22 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.55.
Figure 7.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 17 μm). (d) Corresponding simulated mean topography (MPHSDw, 16 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.37.
Figure 7.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 17 μm). (d) Corresponding simulated mean topography (MPHSDw, 16 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.37.
Figure 8.
 
Qualitative display of between-examination, local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 39 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 38 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.51.
Figure 8.
 
Qualitative display of between-examination, local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 39 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 38 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.51.
Figure 9.
 
Qualitative display of between-examination local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 9 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 9 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.73.
Figure 9.
 
Qualitative display of between-examination local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 9 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 9 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.73.
Discussion
There is no gold-standard method of detecting structural changes in the glaucomatous eye. One solution may be afforded by simulating series of optic nerve head measurements with well-estimated measurement variability characteristics. This study proposes a blueprint for such an approach, using images from CSLO technology. This simulation is based on empirically derived parameters from previous studies and (HRT) values from published reports. The usefulness of a simulation is critically dependent on the degree to which it emulates real data. The simulation presented in this study reproduces well both the spatial characteristics and the magnitude of variability in short-time series of HRT data for both within-examination acquisition and between examinations. Image correlation of real and simulated spatial variability maps tend to be slightly less than the correlation of real spatial variability maps with each other, but the quantitative and qualitative differences are very small. The simulation also generates realistic values for RA, which is a clinically useful parameter for monitoring change. The contribution that simulation can make—independent of and before the collection of real data for scientific studies or clinical trials that are conducted to assess the validity of change detection measures—cannot be overstated. Given any number of baseline patient scans, it enables unlimited modification of the length of series, frequency of testing, and variability over time of the data. It also may reduce the potential waste in cost and time involved in data collection, given that study design may be optimized in advance by testing preliminary hypotheses on simulated test-bed data. 
The simulations were benchmarked against real short-time series. Clinically significant change is highly unlikely to have taken place during the short study period, and age-related effects would also be almost nonexistent; the only changes in the series over time would be real fluctuations in the tissue morphology and variability in the measuring device and settings. This study has used parameters from previous studies on the movement of the eye (and, as a result, the movement of the ONH) during fixation from previous studies, along with the reasonable estimates of operator variability in initial ONH centration and the variability inherent in the device hardware in an attempt to make these assumptions as realistic as possible. Errors in estimating these parameters may result in inaccuracies in the simulation. Other potential sources of variability, such as the ocular pulse and intraocular pressure changes, have not been introduced, as they have not been sufficiently well characterized. The ocular pulse has been shown to contribute to changes in vessel caliber 56,57 and changes in the backward and forward movements of the retinal surface in vessel-free areas. 58 Long- and short-term IOP changes can also cause changes in the vessel caliber 59 along with changes in optic disc cupping 60 and size. 61 The simulation does not account for these effects, but it may be possible to do so using appropriate parameters. Despite these potential shortcomings, the simulations appear to realistically reproduce observed topographic variability. 
Previous attempts to simulate series of ONH or retinal images are few and far between. One approach involved simulating retardance image series obtained from SLP baseline images. 29 SLP differs from CSLO, in that it measures relevant anatomy in a different manner and has different sources of variability. Another simulation methodology, specific to HRT data, has used small misalignments of identical finalized or postprocessed single topographies to generate longitudinal series. 8 Although this approach reproduces some of the global variability inherent in no-change patient HRT data, it reproduces very little of the spatial and feature dependent variability (O'Leary N, et al. IOVS 2009;50:ARVO E-Abstract 2249). It is also unsatisfactory from a theoretical point of view, as it is difficult to estimate the relationship of acquisition noise factors and the variability of postprocessed topographies. In contrast, all the three-dimensional information acquired in the HRT scan is used the methodology presented in this study. Further, a key component of the simulation described in the present study is that the simulated scans are reimported into the HRT III software. Thus, any misalignments added both within examination and between examinations are subject to partial correction by the HRT III alignment algorithms, consistent with the procedure for real patient data. 
Although within-examination single topographies are more closely related than between-session topographies, because of the effects of between-examination repositioning of the patient and realignment of the scanner to the eye/ONH, the simulation assumes no autocorrelation in time (knowing one value does not provide you with any information about another value) between examinations. Also, the within-examination noise parameters do not change with time nor do the differences between examinations. These assumptions are justified on the basis of study results showing no difference in variability between scans with intervals varying from 1 hour to 4 weeks. 48  
There were a few cases in which variability was poorly reproduced in simulated series, as can be seen in Figure 2a. Here, the largest anomalies were above the line of agreement showing that, in these cases, the average, global within-mean topography variability was higher in the simulation than in the real series. This finding may indicate that the optical section stack used to seed the simulation had more motion noise or lower image quality than those in the rest of the real series. Although the simulation model presented herein allows some alignment before noise addition, it is never possible to attain a motion- and noise-free image stack. Thus, the simulation depends on the optical section stack selected to seed each simulated series. Apart from random differences in a small sample of five mean topographies, there may be other potential causes for discrepancies in reproducing the variability of real data, for example, some differences in RA variability may be due to small differences in contour line placement between the real and simulated baseline topographies. 
Windowing, as stated in the methodology, may have small systematic effects on differences in the calculations of variability between real and simulated data due to the potential null values on optical section and topography image edges. However, windowing caused by axial shifts should not affect the topography formation process, as any missing data will be in the low-intensity tails of each axial intensity profile. Topographic height information is mostly derived from the central maximum of each axial intensity profile. Vertical axis shifts have the potential to most affect the topography formation. Despite this, between-examination, vertical translations of scan triplets greater than 20 pixels (±2 SD: P = 0.05) between scans resulting in the loss of more than 20 pixels at the top or bottom of a mean topography will have a small effect on the alignment algorithms of the HRT and the calculation of variability values. 
The simulation makes the simplifying assumption that eye movement and measurement noise parameters are the same for all patients. This assumption seemed sufficient in reproducing many of the characteristics of real data. Furthermore, it is important to note that in real patient data, patients with various degrees of lens opacity were scanned. This variation allowed the effects of a wide range of image quality to be modeled, with some eyes having high within- and between-examination variability. That this was also observed in the simulated data suggests that much of the observed range of image variability is a consequence of the range in the quality of ocular media and the structure of the ONH themselves. 
Future work could usefully consider the formation process in the CSLO pixels (Gruppetta S, et al. IOVS 2009;50:ARVO:E-Abstract 319), allowing the modeling of other processes that occur before the summation of light at the detector for each pixel by modifying the PSF during and between acquisition sessions. These include changes in pupil diameter, in media opacities, and in uncorrected lower order aberrations (defocus and astigmatism). This simulation provides a model that can be used to test progression algorithm specificity appropriate to the baseline image quality, but cannot mimic systematically changing image quality over time. The long-term changes in the optical quality of the eye over time (increased scatter and absorption) resulting in a reduced and noisier light signal reflected to the device detector should be addressed, to mimic this feature in some topographic series over time. 
The introduction in the simulations of topographic changes associated with focal and diffuse glaucomatous damage in the ONH should be investigated. Adding change to a series of replicate HRT scans, and then overlaying measurement noise from this simulation, enables assessment of the sensitivity of change detection algorithms (as the underlying location, temporal sequence, and magnitude of the superimposed change is known). 
A stable series of measurements can be generated by simply adding random noise. Such series may also be used to test the specificity of statistical change detection algorithms developed for detecting progression. However, in reality, there are many other sources of variability in optic nerve head scans, such as eye movements, head movements, and device noise, the effects of which may be falsely detected as true morphologic change by these techniques. Such variability is manifestly more complex than simple additive noise, and therefore a more realistic, stable series simulation is needed to more accurately test the specificity of these algorithms. Stability in any clinical series of data is uncertain, whereas in a modeled series it is certain. 
In conclusion, the methodology outlined in this work has definite use in producing benchmark data on which the false-positive rates of statistical methods for detecting HRT progression can be examined. This method could also improve the confidence in detection of change by these techniques in clinical studies. Simulating series of progressive change in HRT series awaits further investigation. Still, the results of this study are presented as a proof of principle toward a solution of the seemingly intractable problem of developing reliable methods of detecting disease progression in glaucoma. 
Footnotes
 Supported by the Special Trustees of Moorfields Eye Hospital and by an unrestricted grant from Heidelberg Engineering. DGH received a proportion of his funding from the Department of Health's NIHR Biomedical Research Centre at Moorfields Eye Hospital and the UCL Institute of Ophthalmology. The views expressed in this publication are those of the authors and not necessarily those of the UK Department of Health.
Footnotes
 Disclosure: N. O'Leary, None; D.P. Crabb, None; D.F. Garway-Heath, Heidelberg Engineering (F), Carl Zeiss Meditec (F, R), Optovue (F)
The authors thank Steve Gruppetta (City University London), Ed White (Moorfields Eye Hospital), and Michael Reutter (Heidelberg Engineering) for advice and invaluable technical support. 
References
Sihota R Gulati V Agarwal HC Saxena R Sharma A Pandey RM . Variables affecting test-retest variability of Heidelberg Retina Tomograph II stereometric parameters. J Glaucoma. 2002;11:321–328. [CrossRef] [PubMed]
Burgansky-Eliash Z Wollstein G Bilonick RA Ishikawa H Kagemann L Schuman JS . Glaucoma detection with the Heidelberg Retina Tomograph 3. Ophthalmology. 2007;114:466–471. [CrossRef] [PubMed]
Ferreras A Pablo LE Larrosa JM Polo V Pajarín AB Honrubia FM . Discriminating between normal and glaucoma-damaged eyes with the Heidelberg Retina Tomograph 3. Ophthalmology. 2008;115:775–781.e772. [CrossRef] [PubMed]
Tan JC Garway-Heath DF Hitchings RA . Variability across the optic nerve head in scanning laser tomography. Br J Ophthalmol. 2003;87:557–559. [CrossRef] [PubMed]
See JLS Nicolela MT Chauhan BC . Rates of neuroretinal rim and peripapillary atrophy area change: a comparative study of glaucoma patients and normal controls. Ophthalmology. 2009;116:840–847. [CrossRef] [PubMed]
Fayers T Strouthidis NG Garway-Heath DF . Monitoring glaucomatous progression using a novel Heidelberg Retina Tomograph event analysis. Ophthalmology. 2007;114:1973–1980. [CrossRef] [PubMed]
Strouthidis NG Scott A Peter NM Garway-Heath DF . Optic disc and visual field progression in ocular hypertensive subjects: detection rates, specificity, and agreement. Invest Ophthalmol Vis Sci. 2006;47:2904–2910. [CrossRef] [PubMed]
Patterson AJ Garway-Heath DF Strouthidis NG Crabb DP . A new statistical approach for quantifying change in series of retinal and optic nerve head topography images. Invest Ophthalmol Vis Sci. 2005;46:1659–1667. [CrossRef] [PubMed]
Chauhan BC Blanchard JW Hamilton DC LeBlanc RP . Technique for detecting serial topographic changes in the optic disc and peripapillary retina using scanning laser tomography. Invest Ophthalmol Vis Sci. 2000;41:775–782. [PubMed]
Balasubramanian M Bowd C Weinreb RN . Clinical evaluation of the proper orthogonal decomposition framework for detecting glaucomatous changes in human subjects. Invest Ophthalmol Vis Sci. 2010;51:264–271. [CrossRef] [PubMed]
Artes PH Chauhan BC . Longitudinal changes in the visual field and optic disc in glaucoma. Prog Retin Eye Res. 2005;24:333–354. [CrossRef] [PubMed]
Philippin H Unsoeld A Maier P Walter S Bach M Funk J . Ten-year results: detection of long-term progressive optic disc changes with confocal laser tomography. Graefes Arch Clin Exp Ophthalmol. 2006;244:460–464. [CrossRef] [PubMed]
Miglior S Casula M Guareschi M Marchetti I Iester M Orzalesi N . Clinical ability of Heidelberg retinal tomograph examination to detect glaucomatous visual field changes. Ophthalmology. 2001;108:1621–1627. [CrossRef] [PubMed]
Saarela V Falck A Airaksinen PJ Tuulonen A . The sensitivity and specificity of Heidelberg Retina Tomograph parameters to glaucomatous progression in disc photographs. Br J Ophthalmol. 2010;94:68–73. [CrossRef] [PubMed]
Hudson CJW Kim LS Hancock SA Cunliffe IA Wild JM . Some dissociating factors in the analysis of structural and functional progressive damage in open-angle glaucoma. Br J Ophthalmol. 2007;91:624–628. [CrossRef] [PubMed]
Saarela V Airaksinen JP . Heidelberg retina tomograph parameters of the optic disc in eyes with progressive retinal nerve fibre layer defects. Acta Ophthalmol. 2008;86:603–608. [CrossRef] [PubMed]
Vizzeri G Weinreb RN Martinez de la Casa JM . Clinicians agreement in establishing glaucomatous progression using the Heidelberg Retina Tomograph. Ophthalmology. 2009;116: 14–24. [CrossRef] [PubMed]
Chauhan BC McCormick TA Nicolela MT LeBlanc RP . Optic disc and visual field changes in a prospective longitudinal study of patients with glaucoma: comparison of scanning laser tomography with conventional perimetry and optic disc photography. Arch Ophthalmol. 2001;119:1492–1499. [CrossRef] [PubMed]
Chauhan BC Hutchison DM Artes PH . Optic disc progression in glaucoma: comparison of confocal scanning laser tomography to optic disc photographs in a prospective study. Invest Ophthalmol Vis Sci. 2009;50:1682–1691. [CrossRef] [PubMed]
Kourkoutas D Buys YM Flanagan JG Hatch WV Balian C Trope GE . Comparison of glaucoma progression evaluated with Heidelberg retina tomograph II versus optic nerve head stereophotographs. Can J Ophthalmol. 2007;42:82–88. [CrossRef] [PubMed]
O'Leary N Crabb DP Mansberger SL . Glaucomatous progression in series of stereoscopic photographs and Heidelberg Retina Tomograph images. Arch Ophthalmol. 2010;128:560–568. [CrossRef] [PubMed]
Coleman AL Sommer A Enger C Knopf HL Stamper RL Minckler DS . Interobserver and intraobserver variability in the detection of glaucomatous progression of the optic disc. J Glaucoma. 1996;5:384–389. [PubMed]
Azuara-Blanco A Katz LJ Spaeth GL Vernon SA Spencer F Lanzl IM . Clinical agreement among glaucoma experts in the detection of glaucomatous changes of the optic disk using simultaneous stereoscopic photographs. Am J Ophthalmol. 2003;136:949– 950. [CrossRef] [PubMed]
Jampel HD Friedman D Quigley H . Agreement among glaucoma specialists in assessing progressive disc changes from photographs in open-angle glaucoma patients. Am J Ophthalmol. 2009;147:39–44.e31. [CrossRef] [PubMed]
Collins DL Zijdenbos AP Kollokian V . Design and construction of a realistic digital brain phantom. IEEE Trans Med Imaging. 1998;17:463–468. [CrossRef] [PubMed]
Kwan RK Evans AC Pike GB . MRI simulation-based evaluation of image-processing and classification methods. IEEE Trans Med Imaging. 1999;18:1085–1097. [CrossRef] [PubMed]
Crum WR Hartkens T Hill DLG . Non-rigid image registration: theory and practice. Br J Radiol. 2004;77:S140–153. [CrossRef] [PubMed]
Holden M Schnabel JA Hill DLG . Quantification of small cerebral ventricular volume changes in treated growth hormone patients using nonrigid registration. IEEE Trans Med Imaging. 2002;21:1292–1301. [CrossRef] [PubMed]
Vermeer KA Vos FM Lo B . Modeling of scanning laser polarimetry images of the human retina for progression detection of glaucoma. IEEE Trans Med Imaging. 2006;25:517– 528. [CrossRef] [PubMed]
Gardiner SK Crabb DP . Frequency of testing for detecting visual field progression. Br J Ophthalmol. 2002;86:560–564. [CrossRef] [PubMed]
Spry PGD Bates AB Johnson CA Chauhan BC . Simulation of longitudinal threshold visual field data. Invest Ophthalmol Vis Sci. 2000;41:2192–2200. [PubMed]
Vesti E Johnson CA Chauhan BC . Comparison of different methods for detecting glaucomatous visual field progression. Invest Ophthalmol Vis Sci. 2003;44:3873–3879. [CrossRef] [PubMed]
Gardiner SK Crabb DP . Examination of different pointwise linear regression methods for determining visual field progression. Invest Ophthalmol Vis Sci. 2002;43:1400–1407. [PubMed]
Jansonius NM . Progression detection in glaucoma can be made more efficient by using a variable interval between successive visual field tests. Graefes Arch Clin Exp Ophthalmol. 2007;245:1647–1651. [CrossRef] [PubMed]
Strouthidis NG White ET Owen VMF Ho TA Hammond CJ Garway-Heath DF . Factors affecting the test-retest variability of Heidelberg retina tomograph and Heidelberg retina tomograph II measurements. Br J Ophthalmol. 2005;89:1427–1432. [CrossRef] [PubMed]
Ratliff F Riggs LA . Involuntary motions of the eye during monocular fixation. J Exp Psychol. 1950;40:687–701. [CrossRef] [PubMed]
Spauschus A Marsden J Halliday DM Rosenberg JR Brown P . The origin of ocular microtremor in man. Exp Brain Res. 1999;126:556–562. [CrossRef] [PubMed]
Riggs LA Armington JC Ratliff F . Motions of the retinal image during fixation. J Opt Soc Am. 1954;44:315–321. [CrossRef] [PubMed]
Schulz E . Binocular micromovements in normal persons. Graefes Arch Clin Exp Ophthalmol. 1984;222:95–100. [CrossRef] [PubMed]
Zuber BL Stark L . Microsaccades and the velocity-amplitude relationship for saccadic eye movements. Science. 1965;150:1459–1460. [CrossRef] [PubMed]
Engbert R Mergenthaler K . Microsaccades are triggered by low retinal image slip. Proc Natl Acad Sci. 2006;103:7192– 7197. [CrossRef] [PubMed]
Liang J-R Moshel S Zivotofsky AZ . Scaling of horizontal and vertical fixational eye movements. Phys Rev E Stat Nonlin Soft Matter Phys. 2005;71:031909.
Ferman L Collewijn H Jansen TC Van den Berg AV . Human gaze stability in the horizontal, vertical and torsional direction during voluntary head movements, evaluated with a three-dimensional scleral induction coil technique. Vision Res. 1987;27:811–828. [CrossRef] [PubMed]
Technical Specifications of the HRT III. Heidelberg Engineering; 2010: http://www.heidelbergengineering.com/international/products/ hrt-glaucoma/specifications/ . Accessed January 2010.
Pawley JB . Sources of noise in three-dimensional microscopical data sets. In: Buetow DE Cameron IT Padilla GM eds, three-dimensional confocal microscopy:volume investigation of biological specimens. New York: Academic Press; 1994.
Goodman JW . Fundamental Limits in Photoelectric Detection of Light. Statistical Optics. New York: John Wiley & Sons; 1985:465–527.
Vasilescu G . Physical Noise Sources. Electronic Noise and Interfering Signals—Principles and Applications. New York: Springer; 2005:45–68.
Chauhan BC MacDonald CA . Influence of time separation on variability estimates of topographic measurements with confocal scanning laser tomography. J Glaucoma. 1995;4:189–193. [CrossRef] [PubMed]
Yucel YH Gupta N Kalichman MW . Relationship of optic disc topography to optic nerve fiber number in glaucoma. Arch Ophthalmol. 1998;116:493–497. [CrossRef] [PubMed]
Owen VM Strouthidis NG Garway-Heath DF Crabb DP . Measurement variability in Heidelberg Retina Tomograph imaging of neuroretinal rim area. Invest Ophthalmol Vis Sci. 2006;47:5322–5330. [CrossRef] [PubMed]
Leung CK-S Cheung CY-l Lin D Pang CP Lam DSC Weinreb RN . Longitudinal variability of optic disc and retinal nerve fiber layer measurements. Invest Ophthalmol Vis Sci. 2008;49:4886–4892. [CrossRef] [PubMed]
Jampel HD Vitale S Ding Y . Test-retest variability in structural and functional parameters of glaucoma damage in the glaucoma imaging longitudinal study. J Glaucoma. 2006;15:152–157. [CrossRef] [PubMed]
Poli A Strouthidis NG Ho TA Garway-Heath DF . Analysis of HRT images: comparison of reference planes. Invest Ophthalmol Vis Sci. 2008;49:3970–3975. [CrossRef] [PubMed]
Eskicioglu AM Fisher PS . Image quality measures and their performance. IEEE Trans Commun. 1995;43:2959–2965. [CrossRef]
Bland JM Altman DG . Statistical methods for assessing agreement between two methods of clinical measurement. Lancet. 1986;1:307–310. [CrossRef] [PubMed]
Chen HC Patel V Wiek J Rassam SM Kohner EM . Vessel diameter changes during the cardiac cycle. Eye. 1994;8:97–103. [CrossRef] [PubMed]
Knudtson MD Klein BE Klein R . Variation associated with measurement of retinal vessel diameters at different points in the pulse cycle. Br J Ophthalmol. 2004;88:57–61. [CrossRef] [PubMed]
Chauhan BC McCormick TA . Effect of the cardiac cycle on topographic measurements using confocal scanning laser tomography. Graefes Arch Clin Exp Ophthalmol. 1995;233:568–572. [CrossRef] [PubMed]
Shin DH Tsai CS Parrow KA Kim C Wan JY Shi DX . Intraocular pressure-dependent retinal vascular change in adult chronic open-angle glaucoma patients. Ophthalmology. 1991;98:1087–1092. [CrossRef] [PubMed]
Parrow KA Shin DH Tsai CS Hong YJ Juzych MS Shi DX . Intraocular pressure-dependent dynamic changes of optic disc cupping in adult glaucoma patients. Ophthalmology. 1992;99: 36–40. [CrossRef] [PubMed]
Poostchi A Wong T Chan KC . Optic disc diameter increases during acute elevations of intraocular pressure. Invest Ophthalmol Vis Sci. 2010;51:2313–2316. [CrossRef] [PubMed]
Figure 1.
 
Schematic of simulation: HRT single and mean topography formation formed from a single scan volume of CSLO optical sections. Processes are represented by gray boxes with intermediate data states indicated by white boxes and initial and final data states by rounded boxes.
Figure 1.
 
Schematic of simulation: HRT single and mean topography formation formed from a single scan volume of CSLO optical sections. Processes are represented by gray boxes with intermediate data states indicated by white boxes and initial and final data states by rounded boxes.
Figure 2.
 
Bland-Altman plots showing series-wise agreements between (a) average within-examination days (MPHSDw) for real and simulated mean topography data and (b) agreement of between-examination MPHSDb for real and simulated mean topography data. The mean difference (bias) of average MPHSDw was 3.5 μm (95% LoA, −20.9 to 28.8 μm). The mean difference of MPHSDb was 2.0 μm (95% LoA, −5.4 to 9.3 μm). Uniform 95% LoA illustrate only approximate limits of agreement, as the heteroscedasticity of the data is apparent. Dashed lines: statistically significant linear proportional bias for MPHSDb and the significant linear proportional increase of the widths of the 95% LoA for MPHSDw and MPHSDb.
Figure 2.
 
Bland-Altman plots showing series-wise agreements between (a) average within-examination days (MPHSDw) for real and simulated mean topography data and (b) agreement of between-examination MPHSDb for real and simulated mean topography data. The mean difference (bias) of average MPHSDw was 3.5 μm (95% LoA, −20.9 to 28.8 μm). The mean difference of MPHSDb was 2.0 μm (95% LoA, −5.4 to 9.3 μm). Uniform 95% LoA illustrate only approximate limits of agreement, as the heteroscedasticity of the data is apparent. Dashed lines: statistically significant linear proportional bias for MPHSDb and the significant linear proportional increase of the widths of the 95% LoA for MPHSDw and MPHSDb.
Figure 3.
 
Bland-Altman plots showing series-wise agreement between RA CV for real and simulated data. Mean difference between values for real and simulated data were −2.1% (95% LoA: −17.6% to 13.4%). No proportional effects were found on the bias or LoA.
Figure 3.
 
Bland-Altman plots showing series-wise agreement between RA CV for real and simulated data. Mean difference between values for real and simulated data were −2.1% (95% LoA: −17.6% to 13.4%). No proportional effects were found on the bias or LoA.
Figure 4.
 
Distributions of maximum NCCs of pixel SD maps between real and simulated data (a) within examination (PHSDw) averaged over all pair-wise comparisons and (b) between-examination PHSDb. (c) The series average NCCs of all pair-wise pixel SD maps between mean topographies in the same series. The extremes are interpreted as follows: −1, perfectly negatively correlated; 0, uncorrelated; 1, perfectly positively correlated.
Figure 4.
 
Distributions of maximum NCCs of pixel SD maps between real and simulated data (a) within examination (PHSDw) averaged over all pair-wise comparisons and (b) between-examination PHSDb. (c) The series average NCCs of all pair-wise pixel SD maps between mean topographies in the same series. The extremes are interpreted as follows: −1, perfectly negatively correlated; 0, uncorrelated; 1, perfectly positively correlated.
Figure 5.
 
Benchmarking of PHSDw NCC values between real and simulated data: Bland-Altman plot showing agreement of average NCC for real–real PHSDw map comparisons and average NCC for real–simulated PHSDw map comparisons. Solid lines: uniform bias of 0.052 and 95% LoA of 0.039 to 0.065; dashed line: statistically significant linear proportional bias.
Figure 5.
 
Benchmarking of PHSDw NCC values between real and simulated data: Bland-Altman plot showing agreement of average NCC for real–real PHSDw map comparisons and average NCC for real–simulated PHSDw map comparisons. Solid lines: uniform bias of 0.052 and 95% LoA of 0.039 to 0.065; dashed line: statistically significant linear proportional bias.
Figure 6.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 20 μm). (d) Corresponding simulated mean topography (MPHSDw, 22 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.55.
Figure 6.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 20 μm). (d) Corresponding simulated mean topography (MPHSDw, 22 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.55.
Figure 7.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 17 μm). (d) Corresponding simulated mean topography (MPHSDw, 16 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.37.
Figure 7.
 
Qualitative display of within examination, local variability for real and simulated pairs. (a) Real mean reflectance image. (b) Corresponding simulated mean reflectance image. (c) Real mean topography (MPHSDw, 17 μm). (d) Corresponding simulated mean topography (MPHSDw, 16 μm). (e) Log of PHSDw maps of real mean topography. Darker areas represent areas of higher variability. (f) Log of PHSDw maps of corresponding simulated mean topography. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.37.
Figure 8.
 
Qualitative display of between-examination, local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 39 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 38 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.51.
Figure 8.
 
Qualitative display of between-examination, local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 39 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 38 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.51.
Figure 9.
 
Qualitative display of between-examination local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 9 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 9 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.73.
Figure 9.
 
Qualitative display of between-examination local variability for real and simulated pairs. (a) Real series-average mean reflectance image. (b) Corresponding simulated series-average mean reflectance image. (c) Real series-average mean topography (MPHSDb, 9 μm). (d) Corresponding simulated series-average mean topography (MPHSDb, 9 μm). (e) Log of PHSDb maps of real mean topography series. Darker areas represent areas of higher variability. (f) Log of PHSDb maps of corresponding simulated mean topography series. The maximum normalized cross-correlation of the two maps in (e) and (f) is 0.73.
Table 1.
 
Measures of Difference between the Simulated and Real Image Series for Global, Regional, and Local Measures of Variability
Table 1.
 
Measures of Difference between the Simulated and Real Image Series for Global, Regional, and Local Measures of Variability
Measures of Variability Measures of Difference
Uniform Significance of Proportional Dependence
Linear Log-Transformed
Bias 95% LoA Bias 95% LoA Bias 95% LoA
MPHSDw 3.5 μm −20.9 to 28.8 μm P = 0.10 P < 0.01 P = 0.61 P < 0.01
MPHSDb 2.0 μm −5.4 to 9.3 μm P < 0.01 P < 0.01 P = 0.35 P = 0.27
RA CV −2.1% −17.6 to −13.4% P = 0.09 P = 0.08 P = 0.40 P = 0.07
NCC PHSDw 0.052 0.039 to 0.065 P < 0.01 P = 0.25 P < 0.01 P = 0.87
×
×

This PDF is available to Subscribers Only

Sign in or purchase a subscription to access this content. ×

You must be signed into an individual account to use this feature.

×