Pattern Recognition Analysis Reveals Unique Contrast Sensitivity Isocontours Using Static Perimetry Thresholds Across the Visual Field

Abstract

Purpose: To determine the locus of test locations that exhibit statistically similar age-related decline in sensitivity to light increments and age-corrected contrast sensitivity isocontours (CSIs) across the central visual field (VF). We compared these CSIs with test point clusters used by the Glaucoma Hemifield Test (GHT).

Methods: Sixty healthy observers underwent testing on the Humphrey Field Analyzer 30-2 test grid using Goldmann (G) stimulus sizes I-V. Age-correction factors for GI-V were determined using linear regression analysis. Pattern recognition analysis was used to cluster test locations across the VF exhibiting equal age-related sensitivity decline (age-related CSIs), and points of equal age-corrected sensitivity (age-corrected CSIs) for GI-V.

Results: There was a small but significant test size–dependent sensitivity decline with age, with smaller stimuli declining more rapidly. Age-related decline in sensitivity was more rapid in the periphery. A greater number of unique age-related CSIs was revealed when using smaller stimuli, particularly in the mid-periphery. Cluster analysis of age-corrected sensitivity thresholds revealed unique CSIs for GI-V, with smaller stimuli having a greater number of unique clusters. Zones examined by the GHT consisted of test locations that did not necessarily belong to the same CSI, particularly in the periphery.

Conclusions: Cluster analysis reveals statistically significant groups of test locations within the 30-2 test grid exhibiting the same age-related decline. CSIs facilitate pooling of sensitivities to reduce the variability of individual test locations. These CSIs could guide future structure-function and alternate hemifield asymmetry analyses by comparing matched areas of similar sensitivity signatures.

Visual field (VF) assessment in clinical practice is performed using a number of techniques, including static (such as using standard automated perimetry [SAP]) and kinetic perimetry.^1,2 Both techniques use different stimulus characteristics and psychophysical methods to measure the sensitivity to light increments across the VF, and thus, their output measurements are conveyed differently. A recent study³ has shown that in healthy subjects, sensitivity to static and kinetic targets are the same when the psychophysical methods are equated. As such, characteristics of one type of perimetry could then be applied to the other. For example, kinetic perimetry produces asymmetric isopters of equal contrast sensitivity, which are test size–dependent and have unique age-related changes, reflecting differences in the shapes of the Hill of Vision (HoV).^4–8 

Although studies in static perimetry have shown similar asymmetries in the HoV and age-related changes in sensitivity,^9–11 no study to date has sought to group static perimetry thresholds to form contrast sensitivity isocontours (CSIs) across the VF to describe groups of test points of statistically similar sensitivity depicted in a similar way to kinetic perimetry isopters. This principle is supported by evidence that static perimetry sensitivity thresholds are not statistically independent, as implied by the discrete and separate test locations, but are instead significantly correlated: neighboring points may represent locations of similar sensitivity, and increasing distance between points may represent different “regions” of sensitivity (i.e., a different CSI).¹² 

Potential applications of SAP CSIs include grouping points to alter a common limiting factor of perimetric results, namely, noise and variance. Practically, CSIs can potentially group test locations to reduce variance inherent to perimetric testing by pooling a greater number of sample sensitivities,^10,13–17 guiding test point density selection,¹⁸ enhancement of structure-function mapping to improve clinical correlations again by reducing the variability and noise that may be different between structural and functional (i.e., SAP)-based techniques,¹⁹ and reconciliation of static sensitivities and kinetic isopters. Subjectivity and noise contributing to poor structure-function relationships could be overcome by highlighting and grouping areas exhibiting the same anatomic²⁰ and sensitivity signatures. An approach to linking structure and function has been outlined in the past, and, as such, in combination with recent work showing grouped ganglion cell density signatures, identification of SAP CSIs may be useful for improved correlations.²¹ Furthermore, clinically, asymmetry analyses, such as the Glaucoma Hemifield Test (GHT), cluster points into corresponding regions across the midline, according to the retinal nerve fiber layer (RNFL) anatomy.²² CSIs can identify clusters of points that can then be compared with those used by the GHT to determine if the GHT clusters consist of locations exhibiting the same sensitivity signature. 

In the present study, we describe a novel method for clustering static perimetry sensitivity thresholds across the VF. As studies in kinetic perimetry have shown differences in age-related changes in isopters found using different stimulus sizes,^4,8 we first measured the static perimetry sensitivity thresholds of a large cohort of observers of varying ages to determine the effect of age on the thresholds measured using Goldmann sizes I to V (GI–V). We applied two methods of clustering, pattern recognition and hierarchical cluster analysis, to determine the number of unique rates of age-related decline in sensitivity across the VF, with clustered locations representing points with the same age-related change (age-related CSIs). Using these results, we age-corrected the sensitivities of all observers into a 50-year-old equivalent observer to determine the unique age-corrected CSIs for GI–V using the clustering methods, hypothesizing that these CSIs reflect the size-specific differences of the HoV. We also tested the hypothesis that size-specific correction factors produce similar age-corrected sensitivity thresholds as using a uniform correction factor (GIII) across all sizes as performed by previous studies,^9,21,23–25 which would be consistent with the absence of an age effect on photopic spatial summation characteristics in the VF.^9,26 Finally, we examined the concordance of age-corrected CSIs with test points identified by the RNFL distributions (i.e., the GHT clusters). 

Methods

Observers

Sixty healthy observers (mean age: 42.5 years, SD: 16.3 years; 29 males, 31 females) underwent VF testing on the 30-2 test grid of the Humphrey Field Analyzer (HFA) using GI–V in full-threshold mode (Table 1). Results from these observers have been reported, in part, in previous studies, in which the inclusion and exclusion criteria can be found.^9,24,25 Ethics approval was given by the relevant University of New South Wales Ethics committee. The observers gave written informed consent before data collection, and the research was conducted in accordance with the tenets of the Declaration of Helsinki. 

Table 1

 

Characteristics of Study Participants

 

Apparatus and Procedures

The HFA was used to measure sensitivity thresholds at the 75 points (including the fovea, and excluding the two points near the physiological blind spot) of the 30-2 test grid, as per previously described methods.^9,24,25 We report the output sensitivities of the HFA in dB (despite dB being a measure of attenuation specific to the instrument), as we wished to determine the CSIs using a clinically established technique, which is reflective of clinically obtained data.^24,25 

Determination of Age-related Changes in Sensitivity

Output sensitivity thresholds (in dB) were plotted as a function of age (years), to obtain the slope value signifying the age-related change in sensitivity at each location. A slope that deviated from 0 at a P < 0.05 level was considered significant. 

Cluster Analysis: Pattern Recognition

Cluster analysis was performed using pattern recognition, carried out using commercially available satellite imaging software (PCI Geomatica version 10; PCI Geomatics, Richmond Hill, Ontario, Canada). Principles of clustering analysis using pattern recognition software have been described extensively elsewhere in studies of retinal amino acid labeling.^27–29 More recently, its use has been described in structural measurements of the eye²⁰ and the VFs of both healthy subjects (Kalloniatis M, et al. 2016 Imaging and Perimetry Society Abstract, Udine, Italy) and patients with age-related macular degeneration (Choi A, et al. IOVS 2016;57:ARVO E-Abstract 6104) for clustering ganglion cell layer thickness and sensitivity thresholds, respectively. 

In short, it is an iterative procedure that begins with an arbitrary initial cluster boundary. The second and third steps of the procedure involve classification of each pixel to the closest cluster and then calculation of new cluster mean boundaries. Steps two and three are repeated, with assignment of new points and boundaries through the stack of images (N-dimensional space) until clusters do not significantly change with subsequent iterations (Supplementary Fig. S1). The algorithms attempt to minimize the errors between each pixel value and its assigned cluster center, mathematically the process of sums of squares distances, which is equivalent to minimizing the mean squared error. In the present study, the grayscale pixel value represented different levels of sensitivity. 

Images were generated using Adobe Creative Cloud Photoshop Version 6 (Adobe Systems Pty Ltd, Victoria, Australia). Output sensitivity values (in dB) at each test location within the 30-2 test grid were converted into pixel values ranging from 0 (lower dB) to 255 (higher dB) (Fig. 1). Each test location contains pixels with identical pixel values reflective of the average sensitivity value (dB value). Although the pixel location is reflected in the final result, the classification is based on clustering pixel values in N-dimensional space, irrespective of the location from which the pixel was derived. 

Figure 1

 

Method of pattern recognition analysis in the present study. (A) HFA sensitivity thresholds (in dB) were converted into grayscale pixel values (PVs) to obtain a grayscale map. Each square represented one test location within the 30-2. PVs ranged from 0 (darkest, lower dB) to 255 (lightest, higher dB). The two blind spot locations (crossed out) and the background has a PV of 255. (B) Pattern recognition analysis carried out using PCI Geomatica 10. The separabilities between classes were determined using the D_T statistic: the class means (e.g., m₁ and m₂) and distributions were compared across clusters to obtain a D_T value, which provides an approximation of the rate of correct classification. (C) Theme maps were generated once classes had a separability of 1.86 or more (>96% probability of correct classification), and were color coded for clarity.

The imaging process involves qualitative visual inspection in red-green-blue triplets to ascertain areas of correlation to subsequently guide the quantitative assessment.^27–29 To support this process, better separation of the pixel value range was achieved by applying a multiplication factor of 5.25 to average dB sensitivity values (Fig. 1A). We also tested different multiplication factors (2 and 4) to determine if integer scaling had an effect on the final pseudocolor map. The lower multiplication factors resulted in similar patterns to the 5.25 factor but required the use of a greater number of possible input classes in the Iterative Self-Organizing Data Analysis Technique Algorithm (ISODATA) algorithm before converging to a similar final classification result. Overall, fewer than 4% of all points exhibited a mismatch between the different multiplication factors, and therefore the choice of factor did not contribute to the results (i.e., this was within our defined criterion of 96% chance of correct classification for different theme classes). We therefore continue to report results using the 5.25 scaling. 

Unsupervised classification with ISODATA was used for cluster analysis. ISODATA is a method of K-means (a migrating means method) that additionally allows for a different number of clusters to be found, whereas the number of clusters is known a priori in K-means, particularly when using large data sets. In ISODATA clustering, there is automatic splitting of high-variance classes and merging of highly overlapping classes (i.e., that class numbers may be reduced or increased within the range input by the user) (see Supplementary Fig. S1 for additional details).³⁰ Unlike K-means, in ISODATA, clusters are merged if the distance between the cluster center (m) is less than a user-defined threshold, or if the number of pixels within the cluster is less than that limit.³⁰ The starting values used in the unsupervised classification method may result in very different classifications, although the mean squared errors may be similarly small. The desired number of classes was set at a high number to prevent artificial constriction of clusters, but other levels of desired classes were also examined. 

The signature separability of these classes was analyzed using the Transformed Divergence (D_T) statistic.³¹ The D_T statistic assumes that the underlying univariate theme class signal distributions are largely normal. There have been empirically defined probabilities of error by D_T value, spanning from D_T = 0 (approximately 13% probability of correct classification) to D_T = 2.0, which represents 100% probability of correct classification.³¹ Classes that had a D_T<1.86 were merged, until all signatures had a minimum D_T>1.86, which represents a >96% correct classification. After the final pseudocolor image was obtained, the colors of the classes were changed for clarity. 

N-stacks for Cluster Analysis

Cluster analysis was performed to address two questions. First, are there location-specific differences in the age-related decline in sensitivity (i.e., age-related CSIs) across the VF that differ across GI–V? Second, can age-corrected CSIs be determined using static perimetry sensitivity thresholds? Correspondingly, there were two types of N-stacks generated for cluster analysis: first, threshold data were separated by age groups into five (N = 5) decade-spanning brackets: 20 to 29 years (n = 19, mean 24.1 years), 30 to 39 years (n = 10, mean 33.9 years), 40 to 49 years (n = 10, mean 44.7 years), 50 to 59 years (n = 10, mean 55.5 years), and 60+ years (n = 10, mean 66.1 years) and by size (GI–V) for analysis (age-related CSIs); and second, following age-correction (see below), subjects were then randomly assigned to generate N groups representing pooled 50-year-old equivalent data to determine average age-corrected CSIs across each stimulus size. Because the pooled 50-year-old equivalent data were not limited by the decade-spanning age brackets, we tested different combinations of N-stacks to examine for differences in age-corrected CSIs. 

Age Correction of Contrast Sensitivity Thresholds

Previous studies have described methods for age-correcting sensitivity thresholds to facilitate pooling of threshold data of an equivalently aged patient.^9,21,23–25 In the present study, we used two methods of age correction: using a size-specific factor obtained from patients within the present cohort, such as that performed by Wall et al.,^14,32 and alternatively using a single correction factor (GIII) from a previously published study.¹⁰ Previous studies have applied a single-sized correction factor for thresholds obtained using other stimulus sizes,^9,21,23–25 as spatial summation characteristics are similar for different ages.^26,33–35 

Cluster Analysis: Hierarchical Cluster Analysis

Previous studies of the eye have used hierarchical cluster analysis to determine patterns of progression of disease.^36,37 We also used this statistical method to allow comparison with the results of pattern recognition analysis procedure described above. Hierarchical cluster analysis was performed using SPSS Statistics Version 22.0 (IBM Corporation, New York, NY, USA). The agglomerative technique was used to group thresholds into subgroups defined by locations (case label), beginning with each location being its own cluster, and continued until similar clusters were merged together. Within-groups linkage method was used. The first branch of the resultant dendrogram used to establish the initial clusters. 

Conventional clustering methods, such as K-means and hierarchical cluster analysis, present a number of possible clusters, k. There are a number of methods for determining the optimum number of clusters.^38–40 Due to the small number of expected clusters in the present study, k was determined by the separation of each cluster center (mean sensitivity) calculated by using d':  

\(\def\upalpha{\unicode[Times]{x3B1}}\)\(\def\upbeta{\unicode[Times]{x3B2}}\)\(\def\upgamma{\unicode[Times]{x3B3}}\)\(\def\updelta{\unicode[Times]{x3B4}}\)\(\def\upvarepsilon{\unicode[Times]{x3B5}}\)\(\def\upzeta{\unicode[Times]{x3B6}}\)\(\def\upeta{\unicode[Times]{x3B7}}\)\(\def\uptheta{\unicode[Times]{x3B8}}\)\(\def\upiota{\unicode[Times]{x3B9}}\)\(\def\upkappa{\unicode[Times]{x3BA}}\)\(\def\uplambda{\unicode[Times]{x3BB}}\)\(\def\upmu{\unicode[Times]{x3BC}}\)\(\def\upnu{\unicode[Times]{x3BD}}\)\(\def\upxi{\unicode[Times]{x3BE}}\)\(\def\upomicron{\unicode[Times]{x3BF}}\)\(\def\uppi{\unicode[Times]{x3C0}}\)\(\def\uprho{\unicode[Times]{x3C1}}\)\(\def\upsigma{\unicode[Times]{x3C3}}\)\(\def\uptau{\unicode[Times]{x3C4}}\)\(\def\upupsilon{\unicode[Times]{x3C5}}\)\(\def\upphi{\unicode[Times]{x3C6}}\)\(\def\upchi{\unicode[Times]{x3C7}}\)\(\def\uppsy{\unicode[Times]{x3C8}}\)\(\def\upomega{\unicode[Times]{x3C9}}\)\(\def\bialpha{\boldsymbol{\alpha}}\)\(\def\bibeta{\boldsymbol{\beta}}\)\(\def\bigamma{\boldsymbol{\gamma}}\)\(\def\bidelta{\boldsymbol{\delta}}\)\(\def\bivarepsilon{\boldsymbol{\varepsilon}}\)\(\def\bizeta{\boldsymbol{\zeta}}\)\(\def\bieta{\boldsymbol{\eta}}\)\(\def\bitheta{\boldsymbol{\theta}}\)\(\def\biiota{\boldsymbol{\iota}}\)\(\def\bikappa{\boldsymbol{\kappa}}\)\(\def\bilambda{\boldsymbol{\lambda}}\)\(\def\bimu{\boldsymbol{\mu}}\)\(\def\binu{\boldsymbol{\nu}}\)\(\def\bixi{\boldsymbol{\xi}}\)\(\def\biomicron{\boldsymbol{\micron}}\)\(\def\bipi{\boldsymbol{\pi}}\)\(\def\birho{\boldsymbol{\rho}}\)\(\def\bisigma{\boldsymbol{\sigma}}\)\(\def\bitau{\boldsymbol{\tau}}\)\(\def\biupsilon{\boldsymbol{\upsilon}}\)\(\def\biphi{\boldsymbol{\phi}}\)\(\def\bichi{\boldsymbol{\chi}}\)\(\def\bipsy{\boldsymbol{\psy}}\)\(\def\biomega{\boldsymbol{\omega}}\)\(\def\bupalpha{\bf{\alpha}}\)\(\def\bupbeta{\bf{\beta}}\)\(\def\bupgamma{\bf{\gamma}}\)\(\def\bupdelta{\bf{\delta}}\)\(\def\bupvarepsilon{\bf{\varepsilon}}\)\(\def\bupzeta{\bf{\zeta}}\)\(\def\bupeta{\bf{\eta}}\)\(\def\buptheta{\bf{\theta}}\)\(\def\bupiota{\bf{\iota}}\)\(\def\bupkappa{\bf{\kappa}}\)\(\def\buplambda{\bf{\lambda}}\)\(\def\bupmu{\bf{\mu}}\)\(\def\bupnu{\bf{\nu}}\)\(\def\bupxi{\bf{\xi}}\)\(\def\bupomicron{\bf{\micron}}\)\(\def\buppi{\bf{\pi}}\)\(\def\buprho{\bf{\rho}}\)\(\def\bupsigma{\bf{\sigma}}\)\(\def\buptau{\bf{\tau}}\)\(\def\bupupsilon{\bf{\upsilon}}\)\(\def\bupphi{\bf{\phi}}\)\(\def\bupchi{\bf{\chi}}\)\(\def\buppsy{\bf{\psy}}\)\(\def\bupomega{\bf{\omega}}\)\(\def\bGamma{\bf{\Gamma}}\)\(\def\bDelta{\bf{\Delta}}\)\(\def\bTheta{\bf{\Theta}}\)\(\def\bLambda{\bf{\Lambda}}\)\(\def\bXi{\bf{\Xi}}\)\(\def\bPi{\bf{\Pi}}\)\(\def\bSigma{\bf{\Sigma}}\)\(\def\bPhi{\bf{\Phi}}\)\(\def\bPsi{\bf{\Psi}}\)\(\def\bOmega{\bf{\Omega}}\)\begin{equation}{\rm{d}}^{\prime} {\rm{\ }} = {\rm{\ }}\left| {{\rm{\ }}\left( {{x_1}-{x_2}} \right){\rm{\ }}} \right|{\rm{\ }}/{\rm{\ }}(0.{\rm{5\ }} \times {\rm{\ }}{\left( {{\sigma _{\rm{1}}}^{\rm{2}} + {\rm{\ }}{\sigma _{\rm{2}}}^{\rm{2}}} \right)^{0.{\rm{5}}}}{\rm {,}}\end{equation}

where x₁ and x₂, and σ₁ and σ₂ represent the first and second cluster's mean and SD, respectively. Pairs of clusters where d' < 1 were combined systematically, beginning with the lowest d' value, until pairs of d' were all >1. Due to the relatively lower variability of sensitivity thresholds found using larger stimuli^14,32 (e.g., GIV and GV) compared with smaller stimuli^9,41 (e.g., GI and GII), we used an additional criterion for separable groups: the mean sensitivity of clusters needed to be at least 1 dB different, which, in practical terms, is the minimum threshold difference reported by the HFA instrument printout. Clusters with mean sensitivity that were within 1 dB were systematically merged, beginning with the smallest difference, until all clusters were at least 1 dB different.  

Statistical Analysis

Statistical analysis was conducted using GraphPad Prism Version 6 (La Jolla, CA, USA) and SPSS Statistics Version 22.0 (IBM Corporation). Outlier thresholds were identified and excluded using the ROUT Method set at Q = 10% (GraphPad Prism 6).⁴² A D'Agostino and Pearson omnibus normality test (α = 0.05) was performed on the normal cohort for each location, and this did not yield significant results, showing that the sensitivity data were normally distributed at all locations within the 30-2 test grid. 

After removal of outliers,^9,23–25 the normally distributed thresholds obtained from the present study contrasted with the results of Heijl et al.,¹⁰ who reported skewed data at each 30-2 test location. This study used a similar number of participants to that of Heijl et al.,¹⁰ who used n = 88 with two visits and n = 74 with three visits (whereby only second- and third-visit results were used). The reported differences in the present study may be due to the different inclusion criteria used here: we recruited patients with significant prior VF testing experience (to reduce the learning effect), and most patients had all data collection performed on the same day (to reduce intervisit variability).⁴³ This ensured that the obtained threshold values in the present study reliably reflected visual function, as errors due to subject variability were minimized. Importantly, the removal of outlier thresholds led to normally distributed data, allowing the application of the pattern recognition analysis, as class separation statistics is based on normally distributed data sets. 

Data were analyzed using descriptive statistics, paired t-tests, 1-way ANOVA, and 2-way repeated measures ANOVA, as sensitivity data were normally distributed. Post hoc analyses (Tukey's multiple comparisons with Dunn's corrections at α = 0.05) were performed when significant effects were found on ANOVAs. In addition to reporting statistically significant differences (P < 0.05), we also describe clinically relevant differences in sensitivity as equal to or greater than 1 dB (the minimum difference reported by the HFA printout). 

Results

Age-related Changes in the HoV

Figure 2 shows the mean sensitivity for each age group for GI–V along the horizontal meridian, with significant effects of age (average F_4,40 = 192.5, all sizes: P < 0.0001), and of eccentricity (average F_9,10 = 260.1, all sizes: P < 0.0001), similar to previous studies.^9,23 Qualitatively, there appeared to be two main “groups” of sensitivity change through age: the 20 to 49 age group, and the 50+ age group. A similar tendency was evident for the vertical meridian. Post hoc tests showed no difference within the 20- to 49-year age bracket, and no difference within the 50+ age bracket along both horizontal (average P = 0.5166) and vertical (average P = 0.6104) meridians for GI–V, suggesting two “groups” between which there is a marked decline in sensitivity: 20 to 49 years and 50+ years. 

Figure 2

 

Age-related changes across the horizontal HoV for GI–V, plotting sensitivity as a function of eccentricity for each age group. A negative eccentricity value indicates a nasal test location, whereas a positive eccentricity indicates a temporal location. Error bars indicate 1 SD.

Linear Regression Analysis of Sensitivity as a Function of Age

Per-decade decrease in sensitivity values (in dB) for each stimulus size was determined at each location within the 30-2 test grid (Fig. 3). These slopes were significantly different to 0 at all locations (P < 0.05), except at three locations in the superior hemifield for GI (average P = 0.069). The average R² value for the fit of these data was 0.27. There was a size-dependent effect on R² values (H₅ = 83.62, P < 0.0001), with GI having the lowest (0.20) compared with all other sizes. This was likely due to the greater variability of sensitivity thresholds measured using smaller stimulus sizes compared with larger stimuli.^9,14,32,41 

Figure 3

 

Sensitivity slope values (±SD) showing the change in dB per decade of the linear regression analysis for GI–V at each location within the 30-2 test pattern. The fovea result is offset in the upper-left corner for clarity. Three locations at which the slope value was not significantly different from 0 (average: P = 0.069) are highlighted in yellow; all others were significant (P < 0.05).

Slope values were smaller when measured using larger stimulus sizes (GIV and GV) compared with smaller sizes (GI–GIII), reflecting a slower age-related decline. One-way ANOVA revealed a significant difference in slope values between sizes (H₅ = 103.1, P < 0.0001), with post hoc analysis showing significant differences (P < 0.05) between all pairings except for GI/GIII (P > 0.9999), GII/GIII (P = 0.0953), and GIV/GV (P = 0.9844). These age-correction factors thus allowed conversion of sensitivity thresholds from all subjects into an equivalently aged observer. 

The 50-Year-Old Equivalent Observer

Both GII (P = 0.0138) and GIV (P < 0.0001) correction factors, but not GI, GIII, or GV (P value range: 0.1677–0.6766), produced statistically significant differences in sensitivity to those provided by Heijl et al.¹⁰; however, the magnitude of difference was not clinically relevant (i.e., <1 dB, mean [SD] differences: GII, 0.04 ± 0.13 dB; and GIV, −0.08 ± 0.13 dB). Hence, for subsequent analysis, we continued to report age-corrected sensitivity thresholds using correction factors provided by Heijl et al.¹⁰ 

We also compared the sensitivities of a subgroup of 10 subjects with a mean age of 50 years (range: 45–54 years) with the rest of the cohort following linear age correction for each size. Although there were statistically significant differences in sensitivity between the two groups (range of P values: 0.0295 to <0.0001), these were not clinically relevant (range: 0.17–0.48 dB), so we proceeded to use pooled, age-corrected values for subsequent cluster analysis. 

Cluster Analysis: Age-related CSIs

Pattern recognition analysis revealed 8, 9, 8, 7, and 7 resultant theme classes (age-related CSIs) representing locations that exhibited the same sensitivity change with age for GI–V respectively (Fig. 4, column 1; and Table 2). If excluding the “Extra” classes with only a small number of grouped locations (i.e., those at the very superior edge of the 30-2 test grid), then there were 7, 7, 7, 6, and 6 classes for GI–V, respectively. Using GI to GIII, the greater number of age-related CSIs represents more discrete groups of points with the same age-related sensitivity decline, particularly in the mid-peripheral region, in comparison with GIV to V. Within the central 10 to 20°, some CSIs did not form complete “rings,” particularly for GI to III, due to the relative sparseness of sampling by the 30-2 test grid (6° spacing). More peripheral CSIs had a greater rate of sensitivity change compared with more central locations (Fig. 4, column 2). 

Figure 4

 

Pseudocolor theme maps of locations exhibiting the same age-related sensitivity threshold decline for GI–V within the 30-2 test grid following pattern recognition analysis (fovea located in the center). Each color represents a different theme class according to their approximate spatial location. Column 1: Optimal theme maps (D_T values in Table 2). Column 2: Linear regression analyses (sensitivity decline in dB/years) for each theme class. Error bars indicate 1 SD. Columns 3–5: Results of restricting the number of clusters in pattern recognition analysis using ISODATA. Fewer classes resulted in gradual merging of the central and peripheral classes.

Table 2

 

D_T for Pairwise Comparisons Between Clusters Across GI–V Found Using Pattern Recognition Analysis Across Different Age Groups (as Per Fig. 4)

 

Restricting the number of classes allowed by the pattern recognition analysis algorithm resulted in fewer CSIs, but the concentric configuration was relatively well preserved, as adjacent regions were forced into groups sharing similar sensitivities (Fig. 4, columns 3–5). Merging classes, particularly within the central 10 to 20°, formed more “complete rings” by overcoming the wide test location separation. Hierarchical cluster analysis revealed a similar number of clusters and a small number of mismatches (Supplementary Fig. S2; Supplementary Table S1) when compared with the optimal number of classes (Fig. 4, column 1). 

Using the age-related CSIs, test locations were grouped to determine the age-related sensitivity change across different theme classes (Table 3). Similar to ungrouped analysis, there was a tendency for smaller stimuli to have a faster age-related rate of decline, and peripheral locations also changed more rapidly compared with central regions. 

Table 3

 

Linear Regression Slopes Following Grouping of Test Locations Based On Pattern Recognition Analysis Theme Maps in Figure 4. Fovea to Mid-Peripheral, and 2nd Outer to Extra Points Results Have Been Divided for Clarity

 

Cluster Analysis: Age-corrected CSIs

Following age-correction to a 50-year-old equivalent patient, we applied cluster analysis to determine the age-corrected CSIs for GI–V within the 30-2 test grid (Fig. 5). For the pooled analysis, we hypothesized that division of the total cohort (n = 60) into more groups (greater N-stacks) may result in a greater number of resultant classes, due to the increased variance with fewer members per group (i.e., more separation). However, we found no significant differences in the patterns or number of clusters with different levels of N (N = 3, 4, 5, 6, 7, 8, 9, 10). Thus, we divided the cohort into six equal groups of 10 subjects each for subsequent analyses. 

Similar to the results from the age analysis, pattern recognition revealed more age-corrected CSIs when examining smaller stimulus sizes compared with larger stimuli (Fig. 5, Table 4). Again, the greater number of theme classes found using smaller stimuli (GI–III), particularly in the central 10–20°, reflects the more pronounced change in the HoV (i.e., “steeper”) shown in Figure 2, in comparison with larger stimuli (GIV–V) which have less sensitivity change (i.e., a “flatter” HoV) (Fig. 5, column 2). The similarity between the theme classes found between the age analysis (Fig. 4) and the CSIs for the 50-year-old equivalent observer (Fig. 5) reaffirms that the distinct classes found were driven by differences in sensitivities across the VF, rather than just purely age-related effects. 

Figure 5

 

Pseudocolor theme maps of locations exhibiting the same sensitivity signature (i.e., age-corrected CSIs when corrected to a 50-year-old equivalent patient and pooled) for GI–V within the 30-2 test grid (fovea located in the center), color coded as per Figure 4. Column 1: Optimal theme maps (D_T values in Table 4). Column 2: Average sensitivity (dB) for each of theme classes. Error bars indicate 1 SD. Columns 3–5: Results of restricting the number of clusters, as per Figure 4.

Table 4

 

D_T for Pairwise Comparisons Between Clusters Across GI–V Found Using Pattern Recognition Analysis for a 50-Year-Old Equivalent Patient

 

Restricting the number of classes allowed by the pattern recognition analysis algorithm resulted in fewer CSIs, particularly in the central 10 to 20°, representing merging of regions within the HoV that have smaller differences in sensitivity (Fig. 5, columns 3–5). However, the concentric shape was still relatively well preserved with greater restriction of classes. Again, hierarchical clustering showed similar patterns (Supplementary Fig. S3; Supplementary Table S2) to the optimal classes found in Figure 5, column 1. 

Superimposition of the GHT Zones

One current clinically applied index for examining glaucomatous VF defects is the GHT, which groups test locations within the 24-2 test grid into superior and inferior mirrored zones that follow the RNFL distribution (i.e., an anatomical basis). The five zones of the GHT were superimposed on the optimal age-corrected CSI maps for GIII, and also the two other extreme stimulus sizes (GI and GV) (Fig. 6, column 1) to examine whether or not the anatomical basis of the GHT coincided with groups of points within the same age-related CSI (i.e., a contrast sensitivity basis of grouping). Corresponding zones of the GHT mirrored across the horizontal midline consisted of points from multiple theme classes, particularly in zones 3 to 5. Thus, these zones consist of more than one distribution of sensitivities, each with a unique mean and variance value, and are different to the anatomical grouping method. Restricting the number of permitted classes (as per Fig. 5) reduced the amount of discordance in theme classes between corresponding zones, but zones 4 and 5 still exhibited differences in their sensitivity signature (Fig. 6, columns 2–5). 

Figure 6

 

The five zones of the GHT superimposed on 50 year-old equivalent age-corrected CSI theme maps of GI, GIII, and GV for optimal classes (column 1), and following restriction of classes (columns 2–5). Symmetrical zones of the GHT are shown in the bottom right for clarity.

Discussion

Age-related Changes Across the VF

In the present study, we applied a pattern recognition analysis technique in a novel manner to VF data, showing CSIs within the 30-2 test grid to address whether or not there are location- and size-specific differences in the age-related decline in sensitivity. Pointwise analysis showing a nonuniform decline in sensitivity with age across the VF was consistent with previous studies,^10,44 and we have additionally provided, for the first time, stimulus size-specific regression slopes (Fig. 3). However, pointwise linear regression analysis of both VF¹⁰ and anatomic data^45–47 is confounded by the generally poorly fitting linear regressions, reflective of variability between individual points and individual patients. 

Grouping points into age-related CSIs potentially reduces this variability. These age-related CSIs (Fig. 4) predict a downward shift and narrowing of the profile of the HoV (Table 3). In kinetic perimetry in which stimulus contrast and stimulus size are modulated, this manifests as a greater inward constriction of isopters in the periphery compared with central locations.^4,5,8 Thus, the pointwise analysis and grouped age-related CSIs agree with the size- and location-specific sensitivity changes predicted by kinetic perimetry. Location-specific age-related sensitivity declines are consistent with studies suggesting nonuniform loss of retinal detector elements, whereby the peripheral retina has a faster rate of age-related retinal ganglion cell (RGC) loss compared with the macula.^48,49 

Anatomic studies have further suggested that a biphasic, rather than linear fit, is a better representation of age-related RGC changes.^48,49 More recently, Yoshioka et al.²⁰ also showed that a biphasic fit described age-related ganglion cell layer thickness decline better, compared with a linear regression fit. This would be consistent with studies showing stable sensitivity thresholds up until a critical age, with a more rapid decline thereafter.^44,50–53 However, nonlinear regression analysis in the present study did not show significant improvement in the quality of the fit over the linear regression fit (F-test of residuals: average F_2,4 = 1.487, P > 0.05 at all locations) despite grouping test locations, which may be due to a number of previously discussed reasons.^50–52 Specifically, the limited sampling of the 60+ years age bracket reduced the number of data points with which to fit the second slope of the nonlinear regression. Further studies with a greater age range of subjects would be informative. 

Interestingly, we showed that age-related sensitivity decline was greater when using small stimuli compared with large stimuli, yet previous studies have suggested no effect of age on spatial summation characteristics, which would be expected with size-dependent effects.^9,26,33,34 As the two-line spatial summation function relies on the accuracy of a small number of thresholds, variability, such as due to small stimuli and age, can lead to errors in the fit,^54,55 and these can mask subtle age-related changes in spatial summation characteristics that themselves may have significant variance when estimated using a two-line fit method (e.g., Fig. 3 of Khuu and Kalloniatis⁹ and Supplementary Fig. S1 of Phu et al.²⁵). Although the effect of test size on sensitivity change is generally small, the effect may be greater when age differences are large. 

Age-corrected CSIs Across the VF

Determination of age-related changes facilitated examination of age-corrected CSIs. Differences in test procedure means that concordance between static and kinetic perimetry may be difficult to establish.⁷ However, in combination with the work of Phu et al.,³ the present results suggest that kinetic perimetry performed using stimulus size and contrast levels corresponding to the age-corrected CSIs may result in overlap between the isopter and theme classes (Fig. 7). The kinetic perimetry isopters represent the boundary of the “isocontrast area” of equal sensitivity, which may be analogous to the concept of the “fusional area” (Panum's area) in binocular vision.⁵⁶ The age-corrected CSIs found in the present study should theoretically fall within the isocontrast area, as shown in the schematic in Figure 7. Discordance with kinetic perimetry isopters, such as in the “broken” age-corrected CSIs in the peripheral VF, may be improved if the 30-2 test grid test point density and extent are increased. Thus, these age-corrected CSIs could provide guidance for future studies investigating the concordance of static and kinetic perimetry predicted by Phu et al.³ 

Figure 7

 

Superimposition of kinetic perimetry isopters (from a single representative healthy subject) found using Goldmann size I and size II (equivalent to the HFA 24 dB [log Weber contrast 0.103] level of stimulus intensity) and a two-way Method of Limits (peripheral to central “inner” isopter and central to peripheral “outer” isopter) as per Phu et al.³ on the 30-2 spatial map. The colored test locations in the 30-2 test grid reflect age-specific CSIs found in Figure 5, and show approximate concordance with the isocontrast area bounded by the inner and outer kinetic isopters.

Preservation of concentric rings despite restricting the number of theme classes supports the initial hypothesis that the optimum theme maps are biologically meaningful. RGC density and VF sensitivity display concentric contours, representing the anatomical topography of the RGCs and the HoV, respectively.^21,57 We hypothesize that using theme classes derived from RGC density and theme classes derived from VF sensitivity at similar retinal spatial locations may represent biologically meaningful areas and hence may provide improved structure-function correlations. Previous pointwise correlations between structure and function may potentially be limited by variability or noise that differs between measurement devices (e.g., histology and optical coherence tomography for structural measurements, and SAP for functional measurements), and grouping test locations, in reducing variability, may provide better correlations. Grouping test locations may also assist in improving correlations due to the differences in test grid space between functional (e.g., the 6° spacing on the 30-2 VF test grid) and structural measurements (e.g., the approximately 3° width of each test point within the posterior pole analysis grid on the Spectralis OCT; Heidelberg Engineering, Heidelberg, Germany). Further work is required to compare pointwise and grouped structure-function correlations. 

Age-corrected CSIs Applied to the GHT

Cluster analysis also can test the choice of zones for examining VF asymmetries, such as the GHT.²² The GHT relies on the normative distribution of underlying sensitivity thresholds at each location within the zone (i.e., each test point is considered separately). Although some of these points exhibit the same sensitivity signature (and hence distribution), the present results suggest differences between corresponding superior and inferior test zones in terms of the overall theme classes examined. The grouping of different theme classes within the same zone may account for the asymmetries of the underlying normative distribution within the zones, as described by Åsman and Heijl.²² Such a method of grouping would not be restricted to mirrored zones within the superior and inferior hemifields, but would still be biologically meaningful as they represent locations with the same sensitivity rather than points with different sensitivity signatures that could potentially affect the normative distribution. Similar to grouping to perform regression analyses as a function of age, a greater number of pooled locations for asymmetry analyses may reduce variance and could possibly increase the sensitivity of defect detection (see Kalloniatis M, et al., 2017 American Academy of Optometry Conference Paper Session E354A for preliminary results). Future work could investigate the significance of these underlying differences between the zones examined by the GHT and other available cluster analyses, and those reflected in the pattern recognition theme maps. 

In a similar vein, these CSIs may be used in diseases, either through clustering of sensitivities to reduce variability or provide biologically meaningful reference locations, or by comparison of the theme map. The CSIs in the present study represent sensitivities from a cohort of healthy subjects: in patients with disease, reductions in sensitivity across the VF may result in patterns within theme maps that deviate from that of healthy subjects, hence providing another means by which pathologic changes in visual function may be identified. 

Pattern Recognition Analysis Versus Hierarchical Cluster Analysis

The few mismatches between pattern recognition and hierarchical cluster analysis may reflect slight differences in separation criteria (D_T versus d'); however, pattern recognition analysis affords two main advantages. First, it is able to manage large data sets (greater N-dimensions) and complex patterns of points, which would become increasingly unwieldy when manually examining cluster centers using hierarchical clustering. Second, it provides the degree of separability between clusters with an accompanying P value for correct classification without setting additional criteria. The differences in variability between sensitivity thresholds obtained with GI–V required the additional criterion of 1-dB separation. Although this reflects CSIs revealed by conventional static perimetry, the shapes of the CSIs are likely to change when using different measurement units and steps. This requires further investigation with a greater number of stimulus size and contrast combinations. 

Conclusions

In conclusion, cluster analysis revealed test locations within the 30-2 test grid with the same age-related and age-corrected signatures (CSIs). These theme maps could provide a means with which to average sensitivity thresholds across locations with the same sensitivity signature for grouped analyses, such as normative comparison, asymmetry analyses, static and kinetic perimetry concordance, and structure-function studies.^58,59 

Acknowledgments

The authors thank Cornelia Zangerl for technical assistance. 

Supported by a PhD scholarship provided by Guide Dogs NSW/ACT and an Australian Postgraduate PhD scholarship (JP), and the National Health and Medical Research Council of Australia (NHMRC 1033224). Guide Dogs NSW/ACT are partners in the NHMRC grant. Also supported by National Institutes of Health EY015128, EY02576, and EY014800; Vision Core; and an unrestricted grant from Research to Prevent Blindness to the Moran Eye Center. 

Disclosure: J. Phu, P; S.K. Khuu, P; L. Nivison-Smith, None; B. Zangerl, None; A.Y.J. Choi, None; B.W. Jones, None; R.L. Pfeiffer, None; R.E. Marc, None; M. Kalloniatis, P