1. Introduction

Accommodation is a dynamic process. It is well known that the power of the lens changes rapidly and continuously [1–3]. Variations in the steady-state accommodative response are typically less than 0.25 D in magnitude, and are commonly known as microfluctuations [4]. The microfluctuations of accommodation play an important part in the variability of the optical quality of the eye and have therefore attracted much study. Campbell et al. [1] first proposed a description of the main features of the commonly recorded accommodation signal: a low frequency component (< 0.5 Hz), which corresponds to the drift in the accommodation response of the subjects, and a peak at higher frequency, usually observed in the 1–2 Hz band. This description of the frequency composition was also adopted in later studies [3, 5–9]. An area of continued debate is the possible roles that microfluctuations play in the function of accommodation, and the question of whether they are involved in the accommodative control system [3]. The microfluctuations have been shown to be detectable by a normal observer [10], yet their exact role in the accommodative system is not fully understood. It has been debated whether the microfluctuations could be intrinsically related to the accommodative control system with characteristics tuned according to the viewing conditions, or whether they simply represent uncontrolled interference. Another possibility is that the microfluctuations could have characteristics that are independent of the control system, but still provide cues that assist control [3].

The high frequency components are very likely to be a mechanical property of the accommodation system, as they exhibit correlation with the cardiopulmonary system [2, 11–13]. The relationship of the microfluctuations to the mean response of the accommodative system is of primary interest, because the physical nature of the process changes depending on the level of accommodative effort. Several authors have reported that the amplitude of the high frequency component increases with the target vergence [6,14–16].

Variability of the location of peaks in the PSD across subjects had been noted by previous authors [3]. The temporal dynamics of the eye’s focusing power have been characterised as non-stationary [15, 17]. Therefore, periodogram analysis is limited in what it can tell us, because the spectral power of the signal varies over time. Iskander et al. [17], identified low-frequency (<2 Hz) components in Zernike defocus that varied in both magnitude and frequency throughout measurements. Though they are incapable of resolving the time-varying frequency characteristics of signals, periodograms can still be of practical use in interpreting the major characteristics of the accommodative system, by giving an indication of the average spectral power.

Autocorrelation and autocovariance functions have an important role in the study of biomedical signals and systems. These functions can help to detect patterns in biomedical processes, which can give an indication of whether there is any neurological processing involved in the regulation of these processes [18]. Long-term correlation and fractal behaviour has been found in the human heartbeat [19], which is known to be closely related to the microfluctuations of accommodation [4].

In this work, we investigate the effect of accommodative effort on the dynamics of monocular accommodative response. We achieve this by examining their frequency content and autocorrelation structure. We look for consistent trends from subject to subject, independent of their refractive error, in order to characterise some features of the signal.

2. Methods

2.1. Data collection

All data was collected with a custom made aberrometer based on a Shack-Hartmann wavefront sensor. The experimental setup is shown in Figure 1. Data was collected from the subject’s ²), independent of its vergence. It was illuminated by an unfiltered green LED (λ = 0.53μm). The subject was aligned using a bite-bar for stabilisation. The probing beam of the aberrometer was set at a power approximately equal to 15 μW (λ = 0.78 μm), and was scanned over a 1° field angle to remove speckle on the raw frames recorded by the Shack-Hartmann. It was centered slightly off axis, so that it was seen by the subject just outside of the green circular background field of diameter 1°. The light level used for the measurements was fixed before the experiment, using a calibrated power meter (Thorlabs®). At 15 μm, the power of the beam was about 38 times lower than the maximum permitted exposure at the 0.78 μm operating wavelength.

 

Fig. 1. Experimental setup.

Download Full Size | PDF

The wavefront aberration from 9 eyes was measured for 3 different accommodative states: near viewing, far viewing, and an intermediate position. 4 trials of 46.24 s length were performed for each state. A further 4 trials were performed under partial cycloplegia, using ophthalmic drops (1% Tropicamide). The far position was first found by the subject who manually translated the target away from the Badal lens L_B of the fixation arm (as seen in Fig. 1). The instruction was to find a comfortable far viewing position. A cylindrical ophthalmic lens could be introduced in the fixation arm at plane ∏ to correct for the subject’s astigmatism. We then introduced a -4 D ophthalmic lens in the plane ∏. This defined the intermediate point for the trials, which corresponded to the same accommodative demand for all subjects (4 D from the subject’s far point). The near viewing position was then found by directing the subject to translate the target towards L_B. We encouraged the subject to find the limit of his/her accommodative range. The near viewing position was then found by directing the subject to translate the target approximately 1 cm towards L_B, which corresponds to an additional 1.5 D demand. Preliminary measurements were taken to ensure that the accommodative response at the near point of viewing was at least 1 D larger than for the intermediate point. We monitored the pupil constriction during measurements to ensure that this requirement was maintained. If the subject reported a sudden drop in the retinal image quality, the trial was discarded.

This research was approved by The National University of Ireland Research Ethics Committee. Informed consent was obtained from all participants. All subjects were treated according to the Helsinki convention.

2.2. Time series preparation

We used a software pupil to estimate ocular wavefronts over a fixed Φ = 3.9 mm pupil. From the modal reconstruction of the Zernike coefficients of the wavefront (up to the fourth radial order), we computed the accommodative response for the discrete acquisition time n, a(n), in dioptres [20]. This is given in Eq. (1) below:

For each natural viewing condition, we define the mean effort as the mean response minus the mean signal obtained with partial cycloplegia.

The aberrometer was calibrated using a point source located at a distance of 938 ± 1 mm away from the measurement plane. It uses a narrow Gaussian beam, of full width at half maximum 0.5 mm in the measurement plane. The signal to noise of the aberrometer is thus robust to ocular aberrations, and can be used at different accommodative levels without any optical compensation of the defocus of the measured beam. We use a simple algorithm to extend the dynamic range of the sensor [21]. The power of the measured wavefronts was reliably measured after removing the divergence of the reference, to an accuracy of 0.001 D. Because of the chromatic effect of the optics of the eye, the measured accommodative response values were biased by approximately -0.5 D compared to the 0.53 μm central wavelength of our stimulus [22–24]. This bias was not removed in our measurements, so the values given throughout this paper actually correspond to the absolute output of our aberrometer in the near infrared.

To prepare the data for analysis, each trial (consisting of 8,000 data points) was first examined for the presence of spurious values caused by the subject blinking. These unwanted values include instances where the eye is in the process of opening or closing, and also the transient after a blink before the process returns to steady-state [15,17]. Data points corresponding to times when the subject was blinking can often be detected by examining the time series. Methods for the removal of these unwanted values based on the time series values have been proposed in previous work [17], however it was felt that to maximise the efficiency of data removal, the individual Shack-Hartmann frames should be checked manually before removing data points. This procedure was carried out manually for each trial that was used in this study. The blank portions of the time series correspond to instances of blinking. Though the blink interval is typically on the order of 250 μs, the transient components of the measured signal associated with the blink can last for over a second. During this time, the measurement can be biased. Our

2.3. Time series analysis

The classical discrete Fourier transform (DFT) approach to the estimation of PSD admits only uniformly sampled data. Therefore, when one attempts to analyse signals that contain missing values (such as an accommodation signal), the DFT cannot be directly applied. One can attempt to get around this problem by using interpolation to fill in the missing values, however this introduces false values into the signal and gives biased results.

The Lomb-Scargle periodogram [25,26] is a least-squares spectral analysis approach, which employs a modification of the classical definition of the periodogram in order to render it capable of dealing with unevenly sampled and/or incomplete data [27]. The Lomb-Scargle periodogram is particularly suited to the representation of frequency components in time series with missing values, as it avoids the possible bias or erroneous results that may arise from replacement of missing data by interpolation techniques [28]. Also, the frequency values at which the periodogram is to be evaluated may be arbitrarily chosen, and they do not have to be evenly spaced. In this study, blinking by the subjects inevitably led to many missing values in our measured accommodation signals, and so the Lomb-Scargle periodogram was the method we chose in order to estimate the PSD throughout this study. The periodogram P_a(ω) was computed as given in Eqs. 2 and 3 below [26]:

where ω is the angular frequency, x_j is the value of the j-th data point, and τ is defined by

Periodograms were obtained for each trial. The periodograms were estimated from trials with 73 % valid data on average. The spectra were evaluated at 1,000 uniformly spaced frequencies in the range 0.01–10 Hz. In least-squares spectral analysis, the number of invalid points in the time series does not impact on the resolution, because the spectrum is evaluated at arbitrary frequencies [26].

The periodograms were transformed to a logarithmic scale in both axes. Using this representation it was possible to fit a piecewise straight line slope model. A straight line slope in a log-log representation implies a power law relationship [29], in this case between spectral power and frequency. Two separate lines are fitted in each case, comprising a lower frequency region slope m ₁ and a higher frequency region slope m ₂, with the breaking point between them determined empirically by visual inspection of the data. Fitting slopes to the periodograms was achieved using a robust linear regression method.

Many types of commonly observed non-stationary processes are known to have stationary increments. The increments of a time series can provide information about the underlying process, and are useful in the removal of some nonstationary trends [30]. The increments (or first difference) x(n) of the discrete accommodation signal a(n) are given by x(n) = a(n) - a(n - 1).

In order to test time-series for stationarity, we made use of the well-known runs test [31]. This test was chosen because as a non-parametric method it does not require any knowledge about the nature of the system, unlike many model-based methods which require an assumption that the data is normally distributed [32]. The runs test can be based on many different sample statistics [31], e.g. sample mean, sample standard deviation, or mean square error. We chose the sample mean as a suitable statistic for testing the raw accommodation time-series. In each test, the time series under examination was divided into 50 blocks in each case. A 99 % confidence interval was used for the test. Alternatively, one can often receive strong indications on whether or not a time series is non-stationary by simply examining the data, and noting the presence of any trends, discontinuities, transients, or impulsive behaviour. We refer to this process as visual inspection. The visual inspection was performed by manually examining plots of the time series for evidence of trends or transient behaviour [30]. We tested the increments of our measured accommodation time series in the same manner as the raw time series, as seen in Figure 5. Given that the increments are generally zero-mean processes, the sample mean was deemed an unsuitable choice of statistic in this case, and so the sample standard deviation was used instead. Though this test does not conclusively prove non-stationarity, it is a useful indicator, particularly when large amounts of data are available [18]. When dealing with non-stationary processes, traditional methods of analysis such as Fourier analysis are of limited use [17, 18]. Also, the autocorrelation cannot be expressed as the function of a single variable, as in Eq. (4). If one were to drop the assumption of stationarity completely, it would be very difficult to draw any meaningful conclusions from the data. Instead, it is often sufficient to replace the usual assumption of stationarity with a more general notion, such as stationarity of the increments [33].

2.4. Autocorrelation of increments

We also analysed the autocorrelation function (ACF) of the increments of the accommodative response. The ACF r_xx of the discrete random process x(n) as a function of sample lag k can be estimated as given in Eq. (4) below [34]:

where N is the sample size. The ACF of a signal may be used as a quantitative measurement of the memory in the system from which the signal arose [18]. In this study, we computed estimates of the ACF over blocks of N = 500 samples. This enabled us to obtain multiple estimates from a single trial. Any blocks containing data previously marked as invalid were discarded. These estimates were then averaged to give a single ACF estimate for each subject at each viewing condition. The estimates rely on the assumption that the process is wide-sense stationary.

3. Results

3.1. Mean accommodative effort

Figure 2 compares the mean accommodative effort of the 9 subjects at the 3 natural viewing conditions. The effort performed by all subjects at the far point was on average 0.6 D (magenta bars), and corresponds to the accommodation lead usually observed for a relaxed accommodative state. At intermediate viewing, the mean effort was 4.1 D (green). The effort at near viewing corresponds to the maximum accommodative effort that the subject is able to maintain consistently, and was on average 5.9 D (red).

3.2. Spectral analysis

Figure 3 shows periodograms for 3 subjects on a log-log scale, for frequencies up to 10 Hz. We computed slopes m ₁ and m ₂ corresponding to the robust linear fit. Each trace was computed from four trials and averaged. The periodograms tended to converge at frequencies above 10 Hz for all viewing conditions. We attributed this effect to the predominance of the measurement noise, and we therefore ignored this interval when performing the linear fit of the PSD in the high frequency region. Table 1 shows the fitted slope values for all nine subjects at each of the viewing conditions, along with the average value across subjects. We performed t-tests with a 5% significance level to test if the slope values differed significantly over viewing conditions. A Bonferroni correction was included to allow for the data not being independent between accommodation conditions. The null hypothesis was that the slope values were equal across viewing conditions. For m ₁, the test results indicated no significant differences between any of the viewing conditions (0.16 < p < 0.76 in all cases). For m ₂, the values at the intermediate point were found to be different from the values at the near point (p = 0.0003), far point (p = 0.0003), and under partial cycloplegia (p = 0.0003). The tests found no significant differences in the value of m ₂ among the other viewing conditions (0.5 < p < 0.96 in these cases).

 

Fig. 2. Comparison of the mean accommodative effort of the 9 subjects at the 3 natural viewing conditions.

Download Full Size | PDF

Table 1. Fitted slopes for the 9 subjects.

View Table

Figure 4 shows the periodogram averaged across 8 subjects, with the 4 viewing conditions overlaid on a single plot. Subject DDB was excluded from Fig. 4, as the PSD estimates for this subject were significantly different in shape from the other 8 subjects. The reasons for this are unknown, but the difference may be related to the fact that this subject was the only one in the sample who had undergone refractive surgery. The subject was formerly myopic, but had a small degree of hypermetropia post-surgery. Note that in this case, the periodograms were computed at 3,000 uniformly spaced frequencies in the range 0.01–86.5 Hz.

 

Fig. 3. Periodograms of the accommodative response for 3 subjects at each of the viewing conditions with fitted slopes. The values m ₁ and m ₂ denote the fitted slopes for the lower and higher frequency regions respectively.

Download Full Size | PDF

We tested all accommodation measurements for the presence of non-stationarity. Both the runs test and visual inspection approaches were used, and the outcomes are shown in Figure 5. The runs test detected non-stationarity in 75 % of the measured time series, and 6 % of the increments.

3.3. Autocorrelation of increments

Figure 6 shows the autocorrelation function of the increments of the accommodation signal, the increments of Zernike defocus, and the increments of Zernike spherical aberration. Figure 7 shows the normalised ACF of increments of Zernike defocus Z ⁰ ₂ for 3 subjects, at each of the 4 viewing conditions. Each plot is averaged over 4 separate trials.

 

Fig. 4. Averaged periodograms of the accommodative response for different accommodative conditions. Each trace represents the average estimated spectral power across 8 subjects.

Download Full Size | PDF

 

Fig. 5. Assessing the temporal stationarity of the accommodation measurements. There was an overall discrepancy between the two testing methods of 17 % for the raw time series, and 6 % for the increments.

Download Full Size | PDF

 

Fig. 6. Illustration of the effects of noise on the autocorrelation of the increments. The spherical aberration term Z ⁰ ₄ is more greatly affected by noise than the Zernike defocus Z ⁰ ₂ term.

Download Full Size | PDF

 

Fig. 7. Normalised ACF of the increments of the Zernike defocus Z ⁰ ₂ for 3 subjects at each of the 4 viewing conditions.

Download Full Size | PDF

4. Discussion

We have presented measurements of ocular accommodative effort consisting of a large amount of data measured for 9 subjects, with a high-speed, low-noise aberrometer. We found good repeatability over separate trials, and good consistency in our results from subject to subject.

When a subject is fixating at the near limit of their accommodative range, the fluctuations of accommodation are greatly reduced. A small relaxation of the accommodative effort leads to larger microfluctuations. This effect is illustrated by Figure 8. For t > 19 s, the near point signal resembles what we typically observe at intermediate viewing. For t < 19 s, the amplitude of the microfluctuations is reduced in a manner similar to what we typically observe at the far point. Changes of behaviour like this within a single trial at the near point were quite commonly observed, and we consider them as artifacts in our measurements. We suggest that measuring accommodation in subjects at the near limit of their accommodative range might be better served by open view conditions.

 

Fig. 8. Time series measurements of accommodation for subject ED at the 4 viewing conditions. During the near point measurement, the subject was unable to fully hold fixation from t = 19 s onwards.

Download Full Size | PDF

The difference in the shape of the periodogram (Fig. 3) for the intermediate viewing is still clear after averaging (Fig. 4). It is interesting to note that the estimated PSD at the near point is to some degree a mixture of the estimated PSD obtained for the far point and the intermediate point, with an average slope of m ₂ = -2.1 in the 2–10 Hz range. This was expected from observations of the mixed behaviour of the accommodation signal at the near point (Figure 8).

Providing a concise statistical analysis of these time series is a difficult task, because of the non-stationarity of the time series. From plots of the raw time series, we consistently observe differences between the intermediate point measurements and measurements at the other viewing conditions. We attempted to quantify this using standard tools.

A simple approach to quantifying the microfluctuations of accommodation for each measurement trial is to compute the root mean square (RMS) fluctuations of the accommodation signal [5–7,14]. This approach does not show much difference between the intermediate point and the near point, as we found mean RMS values of 0.18 D (near), 0.22 D (intermediate), 0.16 D (far), and 0.11 D (partial cycloplegia). The RMS is very sensitive to drifts of the signal, and is consequently not adequate to quantify non-stationary signals, such as those observed at conditions other than the intermediate point. We therefore pursued other methods to quantify our observations, through spectral analysis and assessing the autocorrelation function of the time series increments.

4.1. Accommodative signal at intermediate viewing

At the intermediate point, we computed a systematically more negative slope of the periodogram in the higher frequency slope range, with an average value of m ₂ = -3.2. At the far point (and with partial cycloplegia), the periodogram could be better fitted by a single straight line over the whole 0.2–10 Hz range. We found average slope values of m ₁ = -1.2 and m ₂ = -2.1 for the far point, and m ₁ = -1.1 and m ₂ = -2.0 for the partial cycloplegia condition. These slope values are comparable to previous studies [15, 35,36].

Looking at the periodogram on a log-log plot emphasises the global shape of the PSD, rather than finer details like peaks associated with the heartbeat and breathing. Such peaks are known to vary in frequency over short periods of time and can be better observed using time-frequency analysis [17]. The correlation between peaks in the PSD of ocular wavefront signals and the cardiopulmonary system have been examined in detail by previous authors [12, 13].

A piecewise linear relationship between spectral power and frequency on a logarithmic scale suggests that the signals could have a degree of self-affinity [18, 37]. Self-affinity and self-similarity have been reported in the temporal fluctuations of many biological processes, including heartbeat [38], blood pressure [39], and balancing [40]. Many of these diverse biological processes exhibit power-law behavior in the PSD. This implies that the current value of the biological signal depends not only on its most recent value but also with its long-term history in a scale-invariant, fractal manner [41]. We observe different slopes in the periodogram depending on the mean accommodative response, in a manner that is consistent from subject to subject. Therefore the scaling of the microfluctuations is related to the accommodative state of the eye, and is significantly altered when the accommodative system is in its active range i.e. between the near point and far point.

The ACF plots of Fig. 7 show unambiguously the difference in the structure of the increments in the accommodation signal at intermediate viewing compared to the other viewing conditions. The slower decay in the ACF suggests longer lasting correlation of the increments of the process for intermediate viewing. This suggests there is some memory in the process in the case of intermediate viewing [18], whereas for the other conditions the ACF more closely resembles physiological noise. We note that the most significant ACF values are at lags below 0.1 s, while Fig. 3 indicated that the estimated PSD of the accommodative response under all the viewing conditions tends to converge at values of approximately 10 Hz and above. It is possible that the memory in the system could play some role in stabilising the accommodative response when the subject is viewing targets at their intermediate point. From the stationarity tests performed on the increments of accommodation, we conclude that it is safe to assume that our accommodation measurements are wide-sense stationary in the increments. This was the motivation for the analysis of the autocorrelation of the increments in accommodation (Figs. 6 and 7).

It is apparent from Fig. 6 that the noise on the measured accommodative signal has an important impact on the ACF. The ACF of the increments of the Zernike defocus Z ⁰ ₂ shows a smoother and longer lasting profile than the full accommodative signal, because the noise on the measured spherical aberration (Zernike Z ⁰ ₄) is increased by a factor of √15 in Eq. (1).

4.2. Stabilisation of the accommodative response

Temporal non-stationarity has been identified as a significant concern when examining aberrometer data [15, 17]. Drifts of the accommodative signal contribute to the low-frequency region of the estimated PSD. The relative absence of these components for intermediate viewing suggests that the accommodation signal could be stationary in this case. To quantify the stabilisation of the accommodative signal, the amplitude of the low-frequency components of the estimated PSD have been analysed [5,8,9,42].

For the three most experienced subjects (CML, CEL, and ED), the slope m ₁ at the intermediate point is noticeably lower than for the other subjects with values of -0.4, -0.7, and -0.4, compared to the average of -0.9. Previous studies of biological processes involving voluntary tasks have shown that statistics of the measurements can change as the subjects gain experience and skill at performing the task [43]. This might indicate that experienced subjects accommodate with improved stability. The results of the stationarity tests (Figure 5) would seem to support this observation. For subjects CML, CEL, and ED, only 1 out of 12 trials of the accommodation signals at the intermediate point were found to contain non-stationarity.