**Purpose**:
A new functional regression model is presented to explain the intersubject variability of the circumpapillary retinal nerve fiber layer (RNFL) thickness in healthy subjects.

**Methods**:
To evaluate the functional regression approach we used data from 202 healthy volunteers, divided equally into training samples (TS) and validation samples (VS). Covariates included RNFL, fovea distance, fovea angle, optic disk ratio, orientation and area provided by Fourier-domain–optical coherence tomography, age, and refractive error. Root mean square errors (RMSE) were calculated for each of the 256 sectors and for the 12 clock-hour sectors in the TS and VS and were compared to the RMSE of the previous model and the standard deviation of the raw data.

**Results**:
With the functional regression approach, we were able to explain on average 27.4% of the variation in the TS and 25.1% of the variation in the VS. The new model performed better compared to a multivariate linear regression model. It performed best in the superior-temporal and inferior-temporal clock-hour sectors where the percentage of RMSE reduction ranged between 26.3% and 44.1% for the TS and between 20.6% and 35.4% for the VS.

**Conclusions**:
The new functional regression approach improves on the multivariate linear regression model and allows an even larger reduction of the amount of intersubject variability, while at the same time using a substantially smaller number of parameters to be estimated.

**Translational Relevance**:
The demonstrated reduction of interindividual variation is expected to translate into an improved diagnostic separation between healthy and glaucomatous subjects, but this remains to be demonstrated in further studies.

^{1}However, these measurements currently present a high intersubject variance,

^{2}mainly due to individual anatomical factors. Because it is the most frequent cause of irreversible blindness in industrialized nations, with an estimate of about 112 million people affected by the disease,

^{3}it is crucial to find more accurate methods that can improve diagnosis.

^{4–6}Specifically, a multivariate linear regression approach including retinal vessel density (RVD)—a parameter reflecting the locations and thicknesses of all measurable circumpapillary retinal vessels—as well as optic disc shape descriptors and fovea parameters allowed a significant decrease in intersubject variability of RNFL.

^{6}We reported a reduction of 18% on average, and up to 29% measured in 12 clock-hour sectors, of the coefficient of variation of RNFL thickness. This previous approach considered 256 independent models around the optic disc (OD) with values of RNFL thickness at a given sector as the dependent variable and all putative parameters, including the values of RVD at the respective sector, as independent factors. Although achieving an actual reduction of intersubject variability, this approach had, from a statistical perspective, some serious shortcomings. Model selection was performed for each of the 256 sectors individually, resulting in potential overfitting due to the large number of parameters to be estimated. Furthermore, the apparent correlation of RNFL between neighboring sectors was not taken into account at all. We propose here a functional regression approach that is meant to overcome these weaknesses.

^{7}whereas more recent developments in terms of application have been described by Morris.

^{8}Here we modeled for each subject the entire RNFL curve, composed of 256 thickness values, using only 16 basis functions. The coefficients corresponding to these basis functions were then used as dependent variables in regression models with eight regressors: vessel thickness plus seven other explanatory variables that were used in the previous multivariate approach. Model selection was performed for each of the 16 linear models. The resulting number of coefficients to be estimated was much smaller than in the multivariate approach, while achieving a better fit of the data, as will be demonstrated by analyzing both the training data set as well as an independent validation sample.

^{6}Shortly, OCT examinations of the OD (cube, 200 × 200) and the macula (cube, 512 × 128) were acquired with Fourier-domain OCT (FD-OCT) (Cirrus; Carl Zeiss Meditec, Inc., Dublin, CA, USA). Scans with a quality index lower than 6 (ranging from 0 to 10) or movement artifacts within the measurement circle were excluded. The image data were exported and analyzed with the Cirrus Research Browser (software version 6.0.2.81), which included segmentation of the RNFL of the 3.4-mm peripapillary circle.

^{6}Using the vessel trees generated in the OD-centered images, we considered all vessels within a band of diameter around the OD center extending from 3.28 to 3.64 mm to integrate a 256-sector vessel profile. This discrete profile was based on vessel thickness and position relative to the OD (angle). The thickness of each vessel was attributed to a specific sector (out of 256 sectors, each 1.4° wide) according to its angle. As opposed to the linear regression model previously presented, we considered only a discrete profile representative of retinal vessel location in a circumpapillary profile, that is, without considering a Gaussian convolution.

^{6}

^{9}will be only briefly outlined here. A more detailed description is given in the Supplementary Material. The first step consisted of landmark registration with respect to three characteristic points of the RNFL curves: the minimum near the fovea and two flanking maximum peaks. After aligning individual curves according to these landmark points, 16 functional principal components were computed and used as base functions to approximately model the aligned RNFL curves. The corresponding principal component analysis (PCA) loadings were then considered as dependent variables in a regression model including vessel thickness and seven additional subject-specific characteristics as explicatory variables. To be able to compare our results with those from Pereira et al.,

^{6}the subject-specific variables included age, spherical refractive error (RE), OD area (ODA), OD orientation (ODO), OD ratio (ODR), as well as fovea angle (FA), which is the angle between the fovea and the OD centers and a horizontal line, and fovea distance (FD) from the OD center. The final step consisted of model selection based on the Akaike Information Criterion (AIC) resulting in a model with 49 parameters.

^{6}were applied to the 101 subjects of TS. RNFL thickness profile, estimated with the functional regression approach for the aligned curves, was transformed back to the original coordinate system by linear interpolation, and residuals were calculated at the 256 original sectors. The different models were assessed according to the root mean square error (RMSE) at each sector. Moreover, for each individual the overall RMSE and the RMSE over 12 clock-hour sectors were calculated for both approaches, and differences were tested using pairwise

*t*-tests. The validation data set underwent the same kind of analysis to compare multivariate linear regression and functional data analysis. The validation step, performed in a sample different from the one used to train the model, ensures that differences in performance were not due to overfitting the training data set.

^{6}(red curve) with the RMSE of the simplest possible approach of considering the overall mean curve without any subject-specific information (null model, black curve). Note that in this last case the RMSE is identical to the standard deviation of the sample.

*P*-value

*t*-test < 0.0001), which reduced the RMSE from 18.97 to 15.42 μm, which translates to a reduction of 18.7%. The new model performed best in the superior-temporal and inferior-temporal areas (clock-hour sectors 2–4 and 10–11), where the percentage of reduction ranged between 26.3% and 44.1%, respectively. In these sectors, the reduction was also significantly higher than the reduction obtained the multivariate model (15.4% to 22.0%). Only in sector 8 did the new functional approach performed significantly worse than the multivariate linear model, with a reduction of 11.7% compared to 21.2%.

*P*-value

*t*-test < 0.0001)

^{6}This new model is based on functional regression, a relatively new development in statistical methodology, which allows for functional modeling of the circumpapillary RNFL thickness profile. In this way, the RNFL thickness is treated as a continuous function and correlation between different sectors is properly accounted for, in contrast to the previous approach in which independent linear models for each sector were considered. A further advantage of the new method is that the number of parameters to be estimated is much smaller than in the previous method, which reduces the danger of overfitting the training data set. Each RNFL curve is no longer represented by 256 points but rather by 19 parameters corresponding to three landmark parameters and 16 basis functions. In terms of estimation, we are dealing with a total of 128 parameters, out of which 49 were selected with AIC to obtain one model representing the RNFL of the complete circumpapillary measurement circle, as compared to 2048 parameters in 256 separate models in the multivariate linear regression approach. The new method also overcomes the issue of having a different best-fitting model in each sector, and as seen in Figure 1, the influence of subject-specific covariates is modeled as a smooth function over all sectors in stark contrast to the previous model.

^{10–12}The difference in our present approach is that both RNFL peaks and the temporal RNFL minimum are used to align the RNFL measurement circle in a way that results in different amounts and directions of rotation around the circle for a given subject. This landmark registration procedure reduces, on average, the intersubject variability by about the same amount as the joint influence of vessels and other subject-specific covariates. There are some clock-hour sectors where our landmark registration slightly increases the RMSE, but only in those sectors where even the null model gives already relatively small errors, which was nasally between 1 and 5 o'clock. In those sectors where the null model is doing particularly badly (close to the peaks P1 and P2, see Supplementary Figure S1), landmark registration gives the largest reduction in RMSE.

^{6}gave results from the multivariate linear regression model in the original coordinate system. In general, all factors had fairly similar behavior along the TSNIT profile for both approaches. In accordance with our previously reported results, age did not show a strong effect throughout the complete RNFL profile, with a loss of 0.036 μm per year on average. In the literature, this loss has mostly been found to be between 0.2 and 0.5 μm per year.

^{13–21}The reason for the weaker correlation in our study might be that our samples are composed by mainly young subjects with a mean age around 31 years. Some authors have published evidence that age-dependent decline in RNFL may be steeper after an age of 40 to 50 years,

^{19,20}although this has not been confirmed in a recent longitudinal study.

^{18}

^{6}The difference between the models might be explained by the fact that the peak location is variable, which is now reflected in the model by the alignment but was not so in the multivariate linear model. Our results confirm previous findings of our own group and others that bigger ODA is associated with thicker RNFL over all sectors but is especially strong in the areas around the two peaks.

^{6,22}

^{6}This might reflect the temporal shift in the RNFL profile in subjects in which the distance between fovea and OD center was increased, as previously reported.

^{23}FA presents a highly positive association with RNFL in superior through superior-nasal areas and a negative influence in nasal and inferior sectors, which may represent the rotation of the superior RNFL peak toward temporal sectors in cases in which the fovea is located more inferiorly to the OD. These findings essentially confirm those obtained with our old approach.

^{6}

^{6}However, with the functional regression model the maximum of this negative association is around the peak areas of RNFL, while in the multivariate linear model there was a positive association at the superior peak. The difference between the two models likely is related to the landmark registration, which now corrects for different peak locations. In the old model, these differences had to be accounted for by other factors such as the vessels and also by the refractive error. It has been previously reported

^{6,24}that in myopic eyes both the retinal vessel arcades and the arcuate bundles of the RNFL tend to be displaced to the temporal side.

^{25}Furthermore, the majority of glaucoma patients show RNFL defects at these locations.

^{26}Reducing the interindividual variation may thus be especially valuable in those sectors.

^{27}This is to be expected since usually a factor with known association with RNFL thickness (for example, age) can be used to either correct for its impact on RNFL thickness or to construct an age-corrected normative database—two sides of the same coin. The effect of age on RNFL thickness as measured by OCT has been first described by Bowd et al. in 2002,

^{13}where it explained less than 17% of the total interindividual variance. Since then, the effect of age on RNFL has been confirmed by numerous authors,

^{14–17,20–21,28}but to the best of our knowledge nobody as yet has evaluated whether compensation for age reduces interindividual variation of RNFL thickness or improves diagnostic separation. Nevertheless, the normative databases of OCT devices take into account the effect of age. In our opinion, this discrepancy between clinical use and lack of scientific proof reflects the natural acceptance by human reasoning that correcting for age will improve the diagnostic performance of OCTs. It is obvious for us that an aged person with thus reduced but healthy RNFL must not be compared with normal limits of young people or even a mixed population without taking into account the effect of age on RNFL. Consequently, as yet nobody has tested this hypothesis. If we take into account several physiologic factors, as we do here, that are all associated with RNFL thickness, the same argument should hold: If these factors explain a considerable part of the interindividual variation (and they do), then this is of clinical relevance and should result in a better diagnostic separation.

^{29}and, in the case of our multivariate model, may also influence some of the described correlations.

^{24,30}In terms of age distribution, both data sets (training and validation) are composed of relatively young subjects (under 40 years old), which may mask or minimize some meaningful impacts that are currently not explicit in the analysis. This fact, however, does not invalidate any of the results presented, and it does not contradict any of the previous reports correlating age and RNFL thickness.

^{20,21,28,31}Moreover, it would be desirable to redefine OD descriptors. The fact that three descriptors (orientation, ratio, and area) are considered may bring some competitiveness between parameters, specifically between ratio and orientation, both dependent on the major axis defined by the OD contour.

**I. Pereira**, None;

**E. Pablik**, None;

**F. Schwarzhans**, None;

**H. Resch**, None;

**G. Fischer**, None;

**C. Vass**, None;

**F. Frommlet**, None

*Br J Ophthal*. 2009; 93: 139–143.

*J Glaucoma*. 2008; 17: 333–340.

*Ophthalmology*. 2014; 121: 2081–2090.

*Br J Ophthalmol*. 2014; 98: 538–543.

*PloS One*. 2015; 10: e0120378.

*Invest Ophthalmol Vis Sci*. 2015; 56: 5290–5298.

*Functional Data Analysis*. 2 ed. New York, NY: Springer-Verlag; 2005.

*Annu Rev Stat Appl*. 2015; 2: 321–359.

*J Stat Softw*. 2012; 51: 1–28.

*Invest Ophthalmol Vis Sci*. 2014; 55: 1419–1426.

*Invest Ophthalmol Vis Sci*. 2014; 55: 7332–7342.

*Br J Ophthalmol*. 2016; 100: 531–536.

*J Opt Soc Am A Opt Image Sci Vis*. 2002; 19: 197–207.

*Br J Ophthalmol*. 2009; 93: 1057–1063.

*Arch Ophthalmol*. 2012; 130: 312–318.

*Korean J Ophthalmol*. 2012; 26: 163–168.

*Invest Ophthalmol Vis Sci*. 2014; 55: 5134–5143.

*Transl Vis Sci Technol*. 2016; 5: 1.

*PLoS One*. 2017; 12: e0179320.

*Ophthalmology*. 2007; 114: 921–926.

*J Glaucoma*. 2013; 22: 532–541.

*Brit J Ophthalmol*. 2005; 89: 489–492.

*Invest Ophthalmol Vis Sci*. 2010; 51: 3515–3523.

*Invest Ophthalmol Vis Sci*. 2013; 54: 5481–5488.

*Ophthalmology*. 2012; 119: 2261–2269.

*Invest Ophthalmol Vis Sci*. 2013; 54: 7338–7343.

*J Glaucoma*. 2010; 19: 158–166.

*Invest Ophthalmol Vis Sci*. 2013; 54: 8095–8103.

*Invest Ophthalmol Vis Sci*. 2012; 53: 4990–4997.

*Invest Ophthalmol Vis Sci*. 2010; 51: 4075–4083

*Invest Ophthalmol Vis Sci*. 2000; 41: 741–748.