## Abstract

The focusing distance of the eye fluctuates during accommodation. However, the visual role of these accommodation fluctuations is not yet fully understood. The fluctuation complexity is one of the obstacles to this long standing challenge in visual science. In this work we seek to develop a statistical approach that i) accurately describes experimental measurements and ii) directly generates randomized and realistic simulations of accommodation fluctuations for use in future experiments. To do so we use the random walk approach, which is usually appropriate to describe the dynamics of systems that combine both randomness and memory.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

The focusing distance fluctuates when the eye fixates a steady target, mostly because of the changing optical power of the crystalline lens. Accommodation fluctuations are often considered a non-stationary random process, with possibly drifting mean [1] and/or time-varying spectral density [2]. Power-laws approximate the global structure of the accommodation spectral density [3–5], although a distribution of power-law exponents better describe typical spectra [6,7]. Despite this mathematical complexity, factors affecting accommodation fluctuations have been well characterized in literature using spectral analysis of time series. Non-stationarity is probably attributable to the low frequency component of the spectrum ($<0.5$ Hz), which increases in amplitude when the viewing conditions make the eye relatively insensitive to defocus error: with low luminance [8], small pupil size [9,10], low spatial frequency of the sinusoidal target [11] and (for some subjects) when chromatic cues are suppressed [12,13]. Besides their characteristics, challenges remain in understanding the role of fluctuations in the accommodation system [14]. One approach to address these questions is to develop optical eye models that simulate accommodation fluctuations during visual experiments [15–17].

Accommodation fluctuations appear to belong to the class of random processes that are non-stationary (with mean, variance, or power spectrum that vary with time), but that have stationary increments [18]. Leahy et al. tested accommodation fluctuations for stationarity using the runs test [1]. They detected non-stationarity in $75\%$ of time series of the accommodation signal, and in only $6\%$ of the corresponding time series of increments. Therefore, accommodation fluctuations are conveniently specified by the statistics of their stationary increments. As in Leahy et al. [1], we define the signal increment as the difference between signals at two consecutive measurement times. With this definition, signal increment is proportional to instantaneous signal velocity. The stationarity of increments simplifies the statistical description of the temporal dynamics of accommodation. Owing to the stationarity of increments, the autocorrelation of increments (ACI) is a function of one single variable (the temporal lag $\tau$):

In Eq. (1), $\delta s(t)$ is the increment of the accommodation signal $s(t)$, in diopters (D). The signal increment is the difference between successive measurements separated by the temporal resolution of measurement ($\delta t$) : $\delta s(t)=s(t+\delta t)-s(t)$. The symbol $\langle \rangle$ stands for ensemble-averaging over a large number of independent observations. It is a probabilistic definition of averaging.

An analogy between accommodation fluctuations and particle dynamics offers an interesting framework to describe their complex temporal dynamics. The temporal dynamics of particle motion are complex because they can have different statistical characteristics depending on the time scale of observation: random versus directional, stable versus drifting. When studying particle dynamics, we may define the particle position recursively and in that case the statistics of position increments defines the particle trajectory as a random walk. The ACI function of particle position quantifies its tendency to continue in the same direction [19]. Positive $ACI(\tau )$ value indicates that there is on average no change of direction in the particle motion during a $\tau$ duration. When the ACI function is positive over the $[\tau _1;\tau _2]$ interval of temporal lags, the particle motion is said to be persistent during that timescale. Persistence occurs for example in trajectories described by the fractional Langevin equation, as a result of particle mass and persistence of the random force exerted by the thermal bath ("fractional noise") [19,20]. Anti-persistence can occur as a result of elasticity of the surrounding medium or potential energy [21]. The reference for particle dynamics is the Brownian motion, for which position increments at two distinct times are not correlated: the ACI function of particle position is null, except for $\tau =0$.

In general, particle trajectories are described by the mean square displacement, which quantifies the fluctuation range of position with time. For a Brownian motion the mean square displacement is proportional to time, and the coefficient of proportionality is twice the diffusivity (for a random walk along one direction). The mean square displacement is a very illustrative statistical quantity. Consider for example a large number of molecules contained in a tiny drop of dye that we disperse in a glass of water at $t=0$. The law of large number states that the surface area of the diffusing drop, at time $t$, is proportional to the mean square displacement. In this work we apply the particle dynamics approach to study accommodation fluctuations. By analogy with the mean square displacement of particle position, we define the mean square defocus (MSD) of the accommodation signal $s(t)$:

The symbol $\langle \rangle$ stands for ensemble-averaging, over a set of independent observations, at fixed measurement time $t$.

In this work we analyze the ACI and MSD functions of the accommodation signal, for both near and far vision. We describe how the two functions mathematically relate. We also describe how the ACI function can be used as model input to perform simulations of accommodation fluctuations as random walks that accurately reproduce the average spectra of accommodation measurements.

## 2. Methods

#### 2.1 Data collection

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. The dominant eye of eight young and healthy subjects was analyzed using a custom-built aberrometer. The aberrometer used near-infrared light (780 nm) and had a high temporal sampling (173 Hz). It measured Zernike coefficients of ocular wavefronts over a fixed pupil of 3.9 mm diameter. Subjects looked at a 6/12 Snellen O through a Badal lens, at photopic light level (80 cd/m$^2$ luminance). The black letter was printed on a white paper and illuminated with green light (530 nm). We measured time series of the accommodative signal for two viewing distances : subjects first fixated the Snellen O at their far point, then a 4 diopters (D) accommodative demand was added to their far point. Four 46.24 s long measurement trials were performed, for each viewing distance. The raw data of accommodation signals have been analyzed in a previous publication [1]. Subjects reported either emmetropia or slight myopia for their dominant eye. With the calibrated aberrometer, we measured their Zernike spherical refraction (from the measured Zernike defocus coefficient $z_2^0$ [22]) in the -0.70 D to +0.20 D range. Mean Zernike refraction was -0.31 D.

#### 2.2 Data analysis

Each trial led to a time series of 8000 Zernike defocus coefficients, $z_2^0(t_k)$, which we defined using the standard convention for vision science [23]. For $1\leq k \leq 8000$, the accommodation signal was computed in diopters (D) as:

where the pupil diameter for Zernike analysis is $\Phi =3.9$ mm. In Eq. (3), we neglected higher-order aberrations in the definition of the accommodation signal. Incorporating spherical aberration in the calculation increased noise significantly and was detrimental to study the microfluctuations in accommodation. In Ref. [1] (Fig. 6), we showed on one example that incorporating spherical aberration in the definition of accommodation signal decreases the amplitude of its ACI function at non-zero time lags. Such decorrelation is characteristic of increased noise. For $1\leq k \leq 7999$, accommodation increment was computed as the difference between successive measurements:We have estimated the ACI (in D$^2$) at time lag $\tau _n=t_{i+n}-t_i$ ($0\leq \tau _n\leq 5.8$ s, for $0\leq n\leq N-1$, with $N=1000$):

The ACI is also defined as an even function for negative time lags, $ACI(-\tau _n)=ACI(\tau _n)$. The ACI is estimated as a temporal average (summation) instead of the ensemble average of Eq. (1). The equivalence between the two methods of averaging is, by definition, valid when the $\delta s(t)$ random process is ergodic and $N$ large [18]. At zero lag, the $ACI(0)$ value is an estimate of the variance of increments, assuming that they have zero mean. We note that the ACI function scales as the square of the temporal resolution of measurement, $\delta t^2$. We can normalize and analyze the $ACI(\tau _n)/\delta t^2$ function, which represents the autocorrelation function of accommodation velocity.

We have computed the MSD (in D$^2$) at time lag $\tau _n$:

The MSD is estimated as a temporal average (summation) instead of the ensemble average of Eq. (2). The equivalence between Eq. (2) and Eq. (6) depends on the ergodicity of the accommodation signal $s(t)$.

To remove the effect of blinks on our analysis, we have inspected the raw CCD frames acquired by the aberrometer. When the eye was not completely open, we disregarded the corresponding data point. Overall, 13$\%$ of the data were disregarded. In Eq. (5) and Eq. (6), we set to zero each member of the sum for which one term is not valid, and then reduce by 1 unit the normalization factor $N-n$. The duration of measurements without a blink was also analyzed (mean 7.8 s) in order to set the duration over which we analyze the ACI and MSD functions. We choose a 5.8 s interval duration ($N=1000$), which is slightly shorter than the typical duration without a blink.

#### 2.3 Statistical analysis

For each viewing distance, the ACI and MSD were estimated for 256 intervals of duration 5.8 s each (256 observations = 8 subjects $\times$ 4 trials $\times$ 8 intervals per trial). We computed the statistical mean and the standard deviation of ACI and MSD when averaging over the 256 observations, at fixed time lag $\tau _n$.

#### 2.4 Relationship between ACI and MSD

In Appendix A, we show that the increment of MSD at time lag $\tau _n$, $\delta MSD(\tau _n)=MSD(\tau _{n+1})-MSD(\tau _{n})$, is equal to the temporal summation of ACI in the $\pm \tau _{n-1}$ interval of time lags:

For mathematical convenience, the derivation of Eq. (7) is given with the probabilistic definitions of ACI and MSD (Eq. (1) and Eq. (2), respectively) and the hypothesis of stationary increments.

Equation (7) agrees with the well-known statistical properties of Brownian dynamics, for which increments at two different times are uncorrelated: $\sum _{i=-(n-1)}^{i=n-1}ACI (\tau _i)=ACI(0)=\sigma ^2$ ($\sigma ^2$ being the variance of increments). For Brownian dynamics, Eq. (7) gives: $\delta MSD(\tau _n)=\sigma ^2$. The MSD is proportional to time lag $\tau _n$: $MSD(\tau _n)=\sigma ^2 \delta t \tau _n$, with $\delta t=t_{k+1}-t_k$ the temporal resolution of measurements. Using the well-known MSD law of Brownian motion, $MSD(\tau _n)=2K\tau _n$, we identify the diffusion coefficient : $K=\sigma ^2\delta t/2$.

The accuracy of Eq. (7), for our experimental estimates of ACI and MSD, is investigated in the Results section.

#### 2.5 Simulation of accommodation fluctuations as a random walk

Owing to their stationarity, accommodation increments conveniently describe accommodation fluctuations. We modeled accommodation fluctuations as a random walk, with the mean ACI function as model input, for a given viewing distance. Each time series of increments ($\delta s(t_k)$, with $k\leq$1000 and $t_k\leq 5.8$ s) is generated as a multivariate normal random vector (size $1000\times 1$), with zero mean and correlation matrix $\textbf {K}$ (size $1000\times 1000$ ) of components $k_{ij}=ACI(t_i-t_j)$. Time series of simulated accommodation signal are then obtained by temporal summation of increments. Simulations were performed with a MATLAB code that is available as supplementary material, Dataset 1, Ref. [24]. We simulated 10000 time series, and we have compared simulations with measurements for each viewing distance in order to check that the random walk accurately simulates the temporal statistics of measurements. First we have compared the MSD of accommodation simulations and measurements. Then we compared the periodograms of accommodation signal simulations and measurements. As in Leahy et al. [1], we have used the Lomb-Scargle periodogram as an estimate of the power spectral density of the accommodative signal. This approach is robust to missing data caused by the subject blinking.

## 3. Results

Figure 1(a) shows the ACI as a function of time lag, for each viewing distance (solid black line: near vision, solid green line : far vision). The solid line is the mean, and the shaded area shows the mean $\pm$ standard deviation of measurements. We estimated the ACI for intervals of duration 5.8 s, and we averaged over 256 observations for each viewing distance (8 subjects $\times$ 4 trials $\times$ 8 intervals per trial). Figure 1(b) zooms in the ACI function in the $0<\tau <0.5$ s interval of time lags. For clarity, the $ACI(0)$ peak was removed in Fig. 1(b). Figure 1(c) shows the MSD as a function of time lag, for each viewing distance (solid black line: near vision, solid green line : far vision). The solid line is the mean, and the shaded area shows the mean $\pm$ standard deviation of measurements. We measured MSD for intervals of duration 5.8 s, and we averaged over 256 observations for each viewing distance (8 intervals $\times$ 4 trials $\times$ 8 subjects).

For far vision, the ACI peaks at zero lag and quickly drops to zero (Fig. 1(a)) like for a Brownian motion. The MSD is approximately proportional to time lag (Fig. 1(c)), which also confirms the approximation as Brownian dynamics. The dashed green line shows the least square linear fit with zero intercept. The slope is 0.0036 D$^2$/s, and corresponds to $2\times$ the diffusion coefficient for a model of Brownian dynamics.

For near vision, the ACI is positive for short time lags ($\tau < 0.1$ s, Fig. 1(b)): accommodation fluctuations are persistent at that time scale. The ACI is negative for intermediate time lags ($0.1<\tau < 0.5$ s, Fig. 1(b)): accommodation fluctuations are anti-persistent at that time scale. For long time lags ($\tau >1$ s, Fig. 1(a)), the ACI is essentially null and the MSD is well approximated by a linear function of time lag (Fig. 1(c)). In Fig. 1(c), the slope of the linear approximation for near vision (dashed black line) is equal to 0.0036 D$^2$/s, as for far vision. The non-zero intercept, for near vision, distinguishes accommodation fluctuations from Brownian dynamics.

Figure 2(a-b) illustrates Eq. (7) with our estimates of ACI and MSD, for near (a) and far (b) vision. The increase rate of MSD at time lag $\tau$ (solid green line) is well approximated by the temporal sum of the ACI function in the $\pm \tau$ interval (dashed black line). Figure 2(a-b) also emphasizes, for near (a) and far (b) vision, the departure from Brownian dynamics, i.e. constant increase rate of MSD. The dashed red line, in Fig. 2(a-b), corresponds to the 0.0036 D$^2$/s value that approximates the increase rate of MSD at the long time scale (slope of dashed lines in Fig. 1(c)). For far vision, a constant increase rate of MSD is a good approximation for $\tau <0.25$ s approximately. For near vision, the increase rate of MSD oscillates near the 0.0036 D$^2/s$ but does not appear to converge. Moreover, in Fig. 2, the departure from Brownian dynamics appears more pronounced for near than for far vision.

Figure 3 shows simulations of accommodation fluctuations, using the ACI function as model input for near and far vision. In Fig. 3(a-b) we compare one simulation (solid green line) with one measurement (solid black line) of accommodative signal over a time interval of 5.6 s, for near (a) and far (b) vision. Simulation fluctuations and drifts are similar to measurements. Only the mean level of accommodation is not reproduced by the simulation. Figure 3(c-d) compares the MSD of measurements (solid black line, mean of 256 measurements) and simulations (solid green line, mean of 10000 simulations), for near (c) and far (d) vision. Shaded areas delimit the mean $\pm$ standard deviation. Simulations reproduce very well the MSD of measurements. The standard deviation of simulations is lower than the standard deviation of measurements, especially at short time lags. Figure 3(e-f) compares the power spectral density of measurements (solid black line: mean of 256 measurements) and simulations (solid green line: mean of 10000 simulations), for near (e) and far (f) vision. Shaded areas delimit the mean $\pm$ standard deviation.

Figure 4 shows, for near vision, the differences in the accommodation fluctuations of two subjects. The condition of near vision in our study was stimulated with a +4 D lens. For subject ED, near vision corresponded to the middle of his/her accommodative range. Mean accommodative effort was 3.67 D for near vision, but subject ED could maintain a stable accommodative effort up to 6.47 D during another experiment (see Fig. 2 in Ref. [1]). Mean accommodative effort for subject EL at near vision (3.31 D) was much closer to his/her maximal effort (4.65 D). For the three descriptions of accommodation fluctuations in Fig. 4 (PSD in Fig. 4(a), ACI in Fig. 4(b), and MSD in Fig. 4(c)), the data of subject ED (green lines) and EL (red lines) differ significantly and lay on both sides of the mean data (black lines). As discussed in [1], the PSD at low spatial frequencies ($\leq 1$ Hz) is flat for subject ED, who is focusing at the middle of his/her accommodative range. The PSD of subject EL resembles more a straight line, as it was observed for most subjects at both ends of their accommodative range [1]. Remarkably, the ACI function (Fig. 4(b)) of subject ED shows an extended interval with negative values ($ACI(\tau )<0$ for $0.11<\tau <0.51$ s), and higher $|ACI(\tau )|$ values in the $\tau <0.5$ s interval. As a result, and according to the model of Eq. (7), the MSD function (Fig. 4(c)) increases faster at short time lags for subject ED and then reaches a plateau for $\tau >2$ s. This plateau is characteristic of stabilized focus.

## 4. Discussion

For near and far vision, we observed persistence of accommodation fluctuations at short time lags ($ACI(\tau )>0$ for $\tau <0.1$ s, Fig. 1(b)) and anti-persistence at intermediate time lags ($ACI(\tau )<0$ for $0.1<\tau <0.5$ s, Fig. 1(b)). Similar scales of persistence and anti-persistence have been reported in a study of eye motion when fixating a steady target [25]. The analysis of eye motions as a random walk brought new perspectives in modeling the process with physiologically-inspired features such as self-avoidance and potential minimization [26]. A self-avoiding random walk accounts for the persistence of eye motion at short time scale and is a strategy to reduce perceptual fading by avoiding looking at recently visited features of the target. Potential energy minimization accounts for the anti-persistence of eye motion at longer time scale and models the ability to fixate a steady target. Here in our work the simulation of accommodation fluctuations as a random walk uses measurement statistics (mean ACI) as model input, and is therefore more of a descriptive model; nevertheless, a physiologically-inspired random walk may also be an interesting approach, for future research, to model accommodation fluctuations.

We found that the mean square defocus (MSD, Fig. 3(c-d)) and the spectral density (PSD, Fig. 3(e-f)) of simulations show good agreement with measurements. We interpret the good agreement in terms of MSD as a consequence of Eq. (7). The good agreement in terms of spectral density was less predictable, and shows that the ACI function captures sufficient information to model accurately accommodation fluctuations as a random walk.

The ACI function is almost null for long time lags ($\tau >1s$, Fig. 1(a)) and, in agreement with Eq. (7), the MSD increase rate is equal ($\simeq 0.0036$ D$^2$/s) for both near and far vision (see linear approximations, as parallel dashed lines in Fig. 1(c)). In that sense, accommodation fluctuations at the long time scales are similar for near and far vision. At the long time scale, the influence of viewing distance on accommodation fluctuations does not reach consensus in the literature [14]. Some studies report an increase in the amplitude of the low frequency component of accommodation with decreasing viewing distance [4,27–30]. Gambra et al., however, do not report such a trend [31] and are more in line with our result. At the intermediate time scale (0.1-0.5 s), the ACI function of the accommodation signal is negative (Fig. 1(a)). Accommodation fluctuations are anti-persistent at that time scale, and we expect that this feature may well describe the ability of subjects to maintain a stable focus. More experiments with variable depth of focus, as described in the literature [8–14], will allow us to study this hypothesis. At the short time scale (high frequency), there is a general agreement in the literature that accommodation fluctuations increase in amplitude for near vision [10,27–30,32,33]. The comparison of MSD functions (near and far, in Fig. 1(c)) clearly confirms this result. Different results have however been obtained concerning the decrease in amplitude for fluctuations at the subject’s near point (see Ref. [27] for example). With only 4 D of accommodation stimulus and young subjects, we likely did not reach the near point for most subjects. Some of the observed inter-subject differences in the accommodation fluctuations may be related to the proximity between the near point and the accommodative state induced by the 4D stimulus. Figure 4 shows that two subjects with distinct accommodative range had different accommodation fluctuations for near vision, at any time scale. Moreover, we note that significant inter-subject differences in accommodation fluctuations have been attributed to individual pulse rates and breathing rates [34,35]. On the one hand, we could choose not to pool the data of different subjects because of those well documented inter-subject differences. On the other hand, the random walk approach introduces probabilistic quantities (Eqs. (1),2) that require averaging over a large number of independent experiments. Because a subject cannot perform a large number of measurement trials we averaged across subjects, and time (Eqs. (5),6), to estimate those probabilistic quantities.

Time scales of (anti-)persistence are similar for near and far vision. However, for far vision, the effect of persistence is barely visible on the MSD plot (Fig. 1(c)) because the $ACI(0)$ peak dominates the ACI curve (Fig. 1(a)) and imposes near-Brownian dynamics at any time scale: the MSD is proportional to temporal lag. Brownian dynamics have long been associated to biological noise because they describe the motion of particles with temperature. Recent advances in particle tracking have shown that there are many examples of Brownian-like particle transport associated to vital functions at cell and organ scales [19]. Accommodation fluctuations may be another example. Often depicted as plant noise [36], in the sense that they are neither controlled nor optimized, accommodation fluctuations can provide temporal cues to guide the accommodation system [37,38]. Previous authors have investigated the detectability [39] and role as a directional cue [40] of accommodation fluctuations, which were artificially introduced in the visual scene as simple harmonics. Similarly, we plan to investigate the role of accommodation fluctuations by introducing either optical defocus or numerical blur [41] in the visual scene as a computer-generated random walk.

## 5. Conclusion

The description of accommodation fluctuations as a random walk allows us to simulate time series with realistic temporal characteristics, which should be useful for designing visual experiments that investigate their role in vision. Persistence of accommodation increments (or velocity) reveals inspiring similarity with eye motion during fixation of a stationary target.

## Appendix A: derivation of Eq. (7)

We start from the probabilistic definition of the mean square defocus (MSD), which we write for the discrete time $t_n$ :

We express the accommodation signal $s(t_n)$ as a function of its increment $\delta s(t_k)$:

Hence we re-write :

We change the index in the second sum ($i=h-k$) and obtain :

We now use the hypothesis of stationary increments and introduce the autocorrelation function of increments (ACI) as a function of time lag $t_i=t_{k+i}-t_k$:

and we obtain :Extracting the $(k=n-1)$ term in the sum gives:

Extracting the $(i=n-1-k)$ term in the sum gives:

We change the index in sum ($j=n-1-k$) and obtain :

We combine the two single-indexed sums into one sum:

According to Eq. (A1) we have:

We identify $MSD(t_{n-1})$ in Eq. (A2) and obtain:

We introduce the increment: $\delta MSD(t_{n-1})=MSD(t_n)-MSD(t_{n-1})$ and we obtain:

We change the index and obtain Eq. (7):

## Acknowledgement

The authors are grateful to Prof. Chris Dainty and Dr. Luis Diaz-Santana for initiating and inspiring this work.

## Disclosures

The authors declare no conflicts of interest.

## Data availability

All data and MATLAB code underlying the results presented in this paper are available in Dataset 1, Ref. [24]. The code also includes the simulation of accommodation signal as a random walk, using the ACI function as model input.

## References

**1. **C. Leahy, C. Leroux, C. Dainty, and L. Diaz-Santana, “Temporal dynamics and statistical characteristics of the microfluctuations of accommodation: Dependence on the mean accommodative effort,” Opt. Express **18**(3), 2668–2681 (2010). [CrossRef]

**2. **D. Iskander, M. Collins, M. Morelande, and M. Zhu, “Analyzing the dynamic wavefront aberrations in the human eye,” IEEE Trans. Biomed. Eng. **51**(11), 1969–1980 (2004). [CrossRef]

**3. **H. Hofer, P. Artal, B. Singer, J. L. Aragón, and D. R. Williams, “Dynamics of the eye’s wave aberration,” J. Opt. Soc. Am. A **18**(3), 497–506 (2001). [CrossRef]

**4. **A. Mira-Agudelo, L. Lundström, and P. Artal, “Temporal dynamics of ocular aberrations: monocular vs binocular vision,” Ophthalmic Physiol. Opt. **29**(3), 256–263 (2009). [CrossRef]

**5. **M. Zhu, M. J. Collins, and D. R. Iskander, “Microfluctuations of wavefront aberrations of the eye,” Ophthalmic Physiol. Opt. **24**(6), 562–571 (2004). [CrossRef]

**6. **K. M. Hampson and E. A. H. Mallen, “Multifractal nature of ocular aberration dynamics of the human eye,” Biomed. Opt. Express **2**(3), 464–470 (2011). [CrossRef]

**7. **K. M. Hampson and E. A. H. Mallen, “Chaos in ocular aberration dynamics of the human eye,” Biomed. Opt. Express **3**(5), 863–877 (2012). [CrossRef]

**8. **L. S. Gray, B. Winn, and B. Gilmartin, “Effect of target luminance on microfluctuations of accommodation,” Ophthalmic Physiol. Opt. **13**(3), 258–265 (1993). [CrossRef]

**9. **L. S. Gray, B. Winn, and B. Gilmartin, “Accommodative microfluctuations and pupil diameter,” Vision Res. **33**(15), 2083–2090 (1993). [CrossRef]

**10. **L. R. Stark and D. A. Atchison, “Pupil size, mean accommodation response and the fluctuations of accommodation,” Ophthalmic Physiol Opt. **17**(4), 316–323 (1997). [CrossRef]

**11. **M. Day, L. S. Gray, D. Seidel, and N. C. Strang, “The relationship between object spatial profile and accommodation microfluctuations in emmetropes and myopes,” J. Vis. **9**(10), 5 (2009). [CrossRef]

**12. **P. B. Kruger, K. R. Aggarwala, S. Bean, and S. Mathews, “Accommodation to stationary and moving targets,” Optom. Vision Sci. **74**(7), 505–510 (1997). [CrossRef]

**13. **D. A. Atchison, N. C. Strang, and L. R. Stark, “Dynamic accommodation responses to stationary colored targets,” Optometry Vision Sci. **81**(9), 699–711 (2004). [CrossRef]

**14. **W. N. Charman and G. Heron, “Microfluctuations in accommodation: an update on their characteristics and possible role,” Ophthalmic Physiol. Opt. **35**(5), 476–499 (2015). [CrossRef]

**15. **J. J. Esteve-Taboada, A. J. Del Águila-Carrasco, I. Marín-Franch, P. Bernal-Molina, R. Montés-Micó, and N. López-Gil, “Opto-mechanical artificial eye with accommodative ability,” Opt. Express **23**(15), 19396–19404 (2015). [CrossRef]

**16. **E. J. Fernández and P. Artal, “Dynamic eye model for adaptive optics testing,” Appl. Opt. **46**(28), 6971–6977 (2007). [CrossRef]

**17. **C. Leahy and C. Dainty, “A non-stationary model for simulating the dynamics of ocular aberrations,” Opt. Express **18**(20), 21386–21396 (2010). [CrossRef]

**18. **J. W. Goodman, * Statistical optics* (John Wiley & Sons, 2015), chap. 3.

**19. **R. Metzler, J.-H. Jeon, A. G. Cherstvy, and E. Barkai, “Anomalous diffusion models and their properties: non-stationarity, non-ergodicity, and ageing at the centenary of single particle tracking,” Phys. Chem. Chem. Phys. **16**(44), 24128–24164 (2014). [CrossRef]

**20. **B. B. Mandelbrot and J. W. Van Ness, “Fractional Brownian motions, fractional noises and applications,” SIAM Rev. **10**(4), 422–437 (1968). [CrossRef]

**21. **I. Goychuk, “Viscoelastic subdiffusion: from anomalous to normal,” Phys. Rev. E **80**(4), 046125 (2009). [CrossRef]

**22. **J. Martin, B. Vasudevan, N. Himebaugh, A. Bradley, and L. Thibos, “Unbiased estimation of refractive state of aberrated eyes,” Vision Res. **51**(17), 1932–1940 (2011). [CrossRef]

**23. **L. N. Thibos, R. A. Applegate, J. T. Schwiegerling, and R. Webb, “Standards for reporting the optical aberrations of eyes,” J. Refract. Surg. **18**(5), S652–S660 (2002). [CrossRef]

**24. **C. Leroux and C. Leahy, “Matlab tools for analyzing accommodation fluctuations as random walk,” figshare2021, https://doi.org/10.6084/m9.figshare.14761779.

**25. **R. Engbert and R. Kliegl, “Microsaccades keep the eyes’ balance during fixation,” Psychologic Sci. **15**(6), 431 (2004). [CrossRef]

**26. **R. Engbert, K. Mergenthaler, P. Sinn, and A. Pikovsky, “An integrated model of fixational eye movements and microsaccades,” Proc. Natl. Acad. Sci. **108**(39), E765–E770 (2011). [CrossRef]

**27. **K. Toshida, F. Okuyama, and T. Tokoro, “Influences of the accommodative stimulus and aging on the accommodative microfluctuations,” Optom. Vision Sci. **75**(3), 221–226 (1998). [CrossRef]

**28. **M. Day, N. C. Strang, D. Seidel, L. S. Gray, and E. A. H. Mallen, “Refractive group differences in accommodation microfluctuations with changing accommodation stimulus,” Ophthalmic Physiol. Opt. **26**(1), 88–96 (2006). [CrossRef]

**29. **E. Gambra, L. Sawides, C. Dorronsoro, and S. Marcos, “Accommodative lag and fluctuations when optical aberrations are manipulated,” J. Vis. **9**(6), 4 (2009). [CrossRef]

**30. **G. L. v. d. Heijde, A. P. A. Beers, and M. Dubbelman, “Microfluctuations of steady-state accommodation measured with ultrasonography,” Ophthalmic Physiol. Opt. **16**(3), 216–221 (1996). [CrossRef]

**31. **E. Gambra, S. Ortiz, P. Perez-Merino, M. Gora, M. Wojtkowski, and S. Marcos, “Static and dynamic crystalline lens accommodation evaluated using quantitative 3-D OCT,” Biomed. Opt. Express **4**(9), 1595–1609 (2013). [CrossRef]

**32. **P. Yao, H. Lin, J. Huang, R. Chu, and B.-C. Jiang, “Objective depth-of-focus is different from subjective depth-of-focus and correlated with accommodative microfluctuations,” Vision Res. **50**(13), 1266–1273 (2010). [CrossRef]

**33. **S. Plainis, H. S. Ginis, and A. Pallikaris, “The effect of ocular aberrations on steady-state errors of accommodative response,” J. Vis. **5**(5), 7 (2005). [CrossRef]

**34. **M. Muma, D. R. Iskander, and M. J. Collins, “The role of cardiopulmonary signals in the dynamics of the eye’s wavefront aberrations,” IEEE Trans. Biomed. Eng. **57**(2), 373–383 (2010). [CrossRef]

**35. **K. Hampson, I. Munro, C. Paterson, and C. Dainty, “Weak correlation between the aberration dynamics of the human eye and the cardiopulmonary system,” J. Opt. Soc. Am. A **22**(7), 1241–1250 (2005). [CrossRef]

**36. **J. C. Kotulak and C. M. Schor, “Temporal variations in accommodation during steady-state conditions,” J. Opt. Soc. Am. A **3**(2), 223–227 (1986). [CrossRef]

**37. **J. C. Kotulak and C. M. Schor, “A computational model of the error detector of human visual accommodation,” Biological Cybernetics **54**(3), 189–194 (1986). [CrossRef]

**38. **B. Winn, “Accommodative microfluctuations: A mechanism for steady-state control of accommodation,” in Accommodation and Vergence Mechanisms in the Visual System (Springer, 2000), pp. 129–140.

**39. **G. Walsh and W. N. Charman, “Visual sensitivity to temporal change in focus and its relevance to the accommodation response,” Vision Res. **28**(11), 1207–1221 (1988). [CrossRef]

**40. **S. Metlapally, J. L. Tong, H. J. Tahir, and C. M. Schor, “Potential role for microfluctuations as a temporal directional cue to accommodation,” J. Vis. **16**(6), 19 (2016). [CrossRef]

**41. **A. J. Del Águila-Carrasco, I. Marín-Franch, P. Bernal-Molina, J. J. Esteve-Taboada, P. B. Kruger, R. Montés-Micó, and N. López-Gil, “Accommodation responds to optical vergence and not defocus blur alone,” Invest. Ophthalmol. Visual Sci. **58**(3), 1758–1763 (2017). [CrossRef]