Strategies for formula constant optimisation for intraocular lens power calculation

Achim Langenbucher; Nóra Szentmáry; Alan Cayless; Jascha Wendelstein; Peter Hoffmann

doi:10.1371/journal.pone.0267352

Abstract

Background

To investigate modern nonlinear iterative strategies for formula constant optimisation and show the application and results from a large dataset using a set of disclosed theoretical-optical lens power calculation concepts.

Methods

Nonlinear iterative optimisation algorithms were implemented for optimising the root mean squared (SoSPE), the mean absolute (SoAPE), the mean (MPE), the standard deviation (SDPE), the median (MEDPE), as well as the 90% confidence interval (CLPE) of the prediction error (PE), defined as the difference between postoperative achieved and formula predicted spherical equivalent power of refraction. Optimisation was performed using the Levenberg-Marquardt algorithm (SoSPE and SoAPE) or the interior point method (MPE, SDPE, MEDPE, CLPE) for the SRKT, Hoffer Q, Holladay 1, Haigis, and Castrop formulae. The results were based on a dataset of measurements made on 888 eyes after implantation of an aspherical hydrophobic monofocal intraocular lens (Vivinex, Hoya).

Results

For all formulae and all optimisation metrics, the iterative algorithms showed a fast and stable convergence after a couple of iterations. The results prove that with optimisation for SoSPE, SoAPE, MPE, SDPE, MEDPE, and CLPE the root mean squared PE, mean absolute PE, mean PE, standard deviation of PE, median PE, and confidence interval of PE could be minimised in all situations. The results in terms of cumulative distribution function are quite coherent with optimisation for SoSPE, SoAPE, MPE and MEDPE, whereas with optimisation for SDPE and CLPE the standard deviation and confidence interval of the PE distribution could only be minimised at the cost of a systematic offset in mean and median PE.

Conclusion

Nonlinear iterative techniques are capable of minimising any statistical metrics (e.g. root mean squared or mean absolute error) of any target parameter (e.g. PE). These optimisation strategies are an important step towards optimising for the target parameters which are used for evaluating the performance of lens power calculation formulae.

Citation: Langenbucher A, Szentmáry N, Cayless A, Wendelstein J, Hoffmann P (2022) Strategies for formula constant optimisation for intraocular lens power calculation. PLoS ONE 17(5): e0267352. https://doi.org/10.1371/journal.pone.0267352

Editor: Asaf Achiron, Edith Wolfson Medical Center, ISRAEL

Received: January 30, 2022; Accepted: April 6, 2022; Published: May 5, 2022

Copyright: © 2022 Langenbucher et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: The minimal dataset has been uploaded to figshare: https://figshare.com/articles/dataset/Minimal_dataset_Strategies_for_formula_constant_optimization/19656765.

Funding: The author(s) received no specific funding for this work.

Competing interests: The authors have declared that no competing interests exist.

Introduction

The refractive power of intraocular lenses can be calculated using either empirical formulae, or so called theoretical-optical formulae based on deep learning algorithms, or with combinations [1, 2] of these. With the SRK or SRK2 formula [1, 3], the power of the intraocular lens (IOL) for emmetropisation is derived from the corneal power (calculated from the corneal radius using a keratometer index (n_K) of 1.3375), axial length (AL) and a formula constant (A), which adapts the formula to a specific lens design. In addition to this simplistic regression formula, several theoretical-optical formulae have been published based on a pseudophakic eye model, in which the IOL power for emmetropisation is based on: the AL, the corneal power calculated from corneal curvature with a keratometer index, other optional biometric measures such as anterior chamber depth (ACD), central thickness of the cornea (CCT) or the crystalline lens (LT), the horizontal corneal diameter (W2W), the age of the patient, and one or more formula constants which again adapt the generally defined formula to the characteristics of a specific lens design [1, 4]. The classical calculation concepts use formula constants A (SRKT formula [3, 5]), personalised anterior chamber depth pACD (Hoffer Q formula [6–8]), surgeon factor SF (Holladay 1 formula [9]), or constant triplets a0/a1/a2 (Haigis formula [1]). All of these theoretical optical formulae based on the Gernet or Fyodorov (published independently in 1970 [10] and 1975 [11]) are restricted to linear Gaussian optics (paraxial optics). The prediction of the axial position of the thin IOL implant (effective lens position, ELP) is mostly performed empirically using the formula constant [12–14]. In addition to these classical calculation concepts, over the last 2 or decades many IOL power calculation strategies such as the Holladay 2, the Barrett, Kane, T2, or DGS formula have been presented, and their relative merits discussed in many scientific reports [15–17]. Only a few of these modern IOL power calculation formulae (for example, the Castrop formula) have been disclosed by the formula authors [16, 18]. In contrast to all these formulae, IOL power calculations using deep learning (e.g. the Hill RBF calculator) are based on a big data setup instead of an eye model, and they do not require formula constants to consider the specific characteristics of a lens model.

It is well known that formula constants are specific not only to IOL models, but also to the clinical population, the surgical technique, and to the measurement techniques (e.g. refractometry) or the equipment (e.g. the biometer). Therefore formula constants should be customised to the surgeon or the surgical centre [19, 20]. For optimisation of formula constants, in addition to the preoperative biometric measures which are considered in the lens power calculation formula, the power of the IOL inserted in the patient eye and the achieved refraction (spherical equivalent, SEQ) after cataract surgery are required. Special attention should be given to the postoperative refraction, which should be documented in a time interval at least 4 weeks after cataract surgery using manual refractometry using trial glasses in a trial frame. Eyes with a postoperative Snellen decimal visual acuity lower than 0.8 should not be considered as the refraction might be unreliable.

However, there is in general no consensus on the optimisation technique for formula constants [18–20]. If the IOL power calculation formula is fully disclosed and uses only one formula constant, the formula can be reorganised to solve for the formula constant. For each individual eye a formula constant is extracted, combining the preoperative biometric measures, the power of the implanted IOL, and the postoperative refraction. Ultimately, any statistical measure such as the mean or median could be used to identify the proper formula constant for this study population from the formula constant distribution. However, when using this ‘optimised’ formula constant in a back-calculation procedure, due to nonlinearity of the formula, we have to be aware that such a procedure does not necessarily result in a zero mean or median prediction error for the spherical equivalent, the least mean absolute prediction error (MAE), the least root mean squared prediction error (RMS), or the least confidence intervals (CL) in the distribution of the prediction error [19]. With more than 1 formula constant in the calculation concept, a direct back calculation of the formula constants for an individual patient eye is not possible. In the simplest case, if all the formula constants together define the ELP in a multilinear superposition, the formula can be solved for the ELP and the formula constants retrieved with a multivariable linear regression analysis from the ELP (e.g. a0/a1/a2 with the Haigis formula). In IOL calculation concepts with more than one formula constant which do not define the ELP with a multilinear superposition (such as the Castrop formula), optimisation procedures fail, and nonlinear optimisation strategies have to be applied [18].

The purpose of the present study is to present a methodology for nonlinear iterative optimising formula constants:

which could be used for any fully disclosed IOL power calculation concept, and
which is in general capable of minimising any statistical measure (e.g. mean, median, standard deviation, confidence interval, or combinations) of any target parameter (e.g. SEQ prediction error).

This optimisation strategy is applied to a large dataset of measurements made after cataract surgery from one clinical centre using one IOL type, in order to explain the meaning of the formula constant optimisation results more in detail.

Materials and methods

Dataset for formula constant optimisation

In this retrospective study we analysed a dataset containing measurements from 888 eyes from a cataract population from Augen- und Laserklinik Castrop-Rauxel, Castrop-Rauxel, Germany which was transferred to us (490 right eyes and 398 left eyes; 495 female and 392 male). The mean age was 71.2±9.1 years (median: 71 years, range: 47 to 91 years). The study was registered with the local ethics committee (Ärztekammer des Saarlandes, registration number 157/21), and a patient informed consent was not required for this study. The data were transferred to us in an anonymised fashion, which precludes back-tracing of the patient. The anonymised data contained preoperative biometric data derived with the IOLMaster 700 (Carl-Zeiss-Meditec, Jena, Germany) including axial length AL, anterior chamber depth ACD measured from the corneal front apex to the anterior apex of the crystalline lens, lens thickness LT, and the corneal front surface radius measured in the flat (R1) and in the steep meridian (R2). In all cases a Vivinex 1 piece hydrophobic aspherical (aberration correcting) monofocal intraocular lens (Hoya Surgical Optics, Singapore) was inserted. In addition to the refractive power of the inserted lens (PIOL), the postoperative refraction (sphere and cylinder) 5 to 12 weeks after cataract surgery was measured by an experienced optometrist and recorded in the dataset. The dataset included only data with a postoperative Snellen decimal visual acuity of 0.8 (20/25 Snellen lines) or higher to ensure that the postoperative refraction was reliable. The descriptive data on pre-cataract biometry, PIOL and postoperative refraction are summarised in Table 1.

Download:

Table 1. Descriptive statistics of the dataset with mean, standard deviation (SD), median, and the lower (quantile 5%) and upper (quantile 95%) boundary of the 90% confidence interval.

https://doi.org/10.1371/journal.pone.0267352.t001

The anonymised Excel data (.xlsx-format) was imported into MATLAB to Matlab (Matlab 2019b, MathWorks, Natick, USA) for further processing.

Preprocessing of the data

Custom software was written in Matlab. The patient age was derived from the date of cataract surgery and date of birth. The mean corneal radius of curvature Rmean was calculated as Rmean = ½(R1+R2), and the mean corneal power Kmean was derived from R1 and R2 as Kmean = ½((n_K-1)/R1+(n_K-1)/R2) with a keratometer index n_K as indicated in the formula definition. The following lens power calculation formulae were considered in this constant optimisation process:

SRKT formula published by Sanders, Retzlaff and Kraff [3, 5],
Hoffer Q formula published by Hoffer [6–8],
Holladay 1 formula published by Holladay and Prager [9],
Haigis formula [1], as well as the
Castrop formula published by Wendelstein et al. and Langenbucher et al. [16, 18].

The SRKT, Hoffer Q, and Holladay 1 formulae consider the AL and corneal curvature/power data together with one formula constant (A, pACD, and SF, respectively), the Haigis formula considers the AL, ACD, and corneal curvature together with a formula constant triplet a0/a1/a2s, and the Castrop formula considers AL, CCT, ACD, LT and corneal curvature of the front and back surface together with a formula constant triplet C/H/R. For the Haigis formula we used 2 versions: the simplified option (Haigis1) with preset values a1 = 0.4/a2 = 0.1 and customisation of a0, and additionally the option with a customised formula constant triplet a0/a1/a2 (Haigis3) [19]. For simplicity when using the Castrop formula, instead of using the CCT and the curvature data of the corneal back surface, we replaced CCT by 0.55 mm and the corneal back surface curvature was estimated from the corneal front surface curvature with a fixed anterior to posterior curvature ratio of 0.834 as published by Langenbucher et al. [18], and a sum of segments correction according to Cooke [21, 22] was used for the axial length. All formulae included in this analysis were reorganised and solved for the SEQ as a function of preoperative biometrical data and PIOL. The difference between the achieved SEQ (from the postoperative follow-up examination) and the SEQ predicted by the formula was considered as the formula prediction error PE.

Formula constant optimisation

In this context, formula constants were optimised for statistical metrics of the SEQ prediction error PE. For the statistical metrics we used: the sum of squared PE (SoSPE) which minimises the ‘energy’ of the prediction error, the sum of the absolute PE (SoAPE, as typically used in scientific reports on the performance of lens power prediction formulae), zeroing of the mean the PE (MPE), the standard deviation of the PE (SDPE) and a linear combination of both (SDMPE), and zeroing of the median PE (MEDPE), the 90% confidence interval (CLPE) and combination of both (CLMEDPE). For all formulae under test (SRKT, Hoffer Q, Holladay 1, Haigis1, Haigis3, and Castrop) formula constant optimisation programming code was written in Matlab for the above mentioned metrics.

For optimisation of the PE for minimising the SoSPE and SoAPE we used the Levenberg-Marquardt algorithm, which is also quoted in the literature as a damped least-squares method and in general solves non-linear least squares problems. This optimisation technique was first described by Levenberg in 1944 [23] and 2 decades later rediscovered by Marquardt [24]. This algorithm is typically used in least squares curve fitting [25] and in most cases converges faster than simple back-propagation methods. It combines the classical Gauss-Newton algorithm with the gradient descent algorithms and is–in most applications–more robust in terms of finding the global minimum, with the potential drawback that the convergence may be somewhat slower compared to the Gauss-Newton algorithm.

Formula constant optimisation for MPE, SDPE, SDMPE, MEDPE, CLPE, and CLMEDPE was implemented using the interior point methods, which refer to a family of optimisation techniques for solving linear and nonlinear convex optimisation problems. This algorithm is also quoted as the barrier method. First published by Dikin in 1967 [26] and 2 decades later re-discovered by Karmarker [27], this algorithm is very efficient in terms of minimising the number of iterations and the number of function evaluations, but seems to be less robust in searching a global minimum [28]. It can be directly implemented with linear programming techniques, and in many situations shows a better performance compared to the simplex algorithm.

Both optimisation techniques–the Levenberg-Marquardt algorithm [23, 24] and the interior point method [26, 27]–were used with box constraint boundaries. The formula constants used for initialisation were: 118.9 (boundaries 116.0 to 121.0) for the SRKT formula, pACD = 5.4 (boundaries from 4.5 to 6.3) for the Hoffer Q, SF = 1.5 (boundaries from 0.5 to 3.0) for the Holladay, a0 = 1.8 (boundaries from -1.0 to 2.5) for the Haigis1, a0/a1/a2 = 1.8/0.4/0.1 (boundaries from -1.0 to 2.5/0.0 to 0.8/0.0 to 0.3 for the Haigis3, and C/H/R = 0.4/0.0/0.0 (boundaries from 0.25 to 0.45/-0.35 to 0.35/-0.35 to 0.35 for the Castrop formula, respectively. After optimising the formula constants for SoSPE, these constants were used as presets for all other optimisations, and the respective boundaries were set symmetrically to ±0.5, ±0.5, ±0.5, ±0.5, ±0.5/±0.1/±0.1, and ±0.1/±0.35/±0.35, respectively. The box constraints were not reached in any iteration step of any of the optimisations.

The SEQ prediction was back-calculated using the optimised constants for each formula and each optimisation metric, and the prediction error PE was derived and analysed. The performance of the formula constant optimisation process was documented with the same metrics which were used earlier for the constant optimisation (sum of squared PR, sum of absolute PE, mean PE, SD of PE, median PE, and 90% CL of PE).

Results

In general, the convergence of the iterative optimisation process for all formulae was faster when optimising for the sum of squared PE (up to 26 iterations and 152 function evaluations), the mean PE (up to 30 iterations and 148 function evaluations) or the median PE (up to 32 iterations and 177 function evaluations), as compared to the optimisation for the sum of absolute PE (up to 48 iterations and 255 function evaluations) standard deviation (up to 50 iterations and 311 function evaluations), 90% confidence interval (up to 48 iterations and 198 function evaluations) or combinations of standard deviation and mean or 90% confidence interval and median PE (up to 77 iterations and 510 function evaluations). Convergence was achieved in all optimisations for all formulae under test.

The resulting optimised formula constants in terms of minimising the sum of squared prediction error (SoSPE), sum of absolute prediction error (SoAPE, mean prediction error (MPE), standard deviation of prediction error (SDPE), combinations of standard deviation and mean prediction error (SDMPE), median prediction error (MEDPE), 90% confidence interval of prediction error (CLPE), and combinations of 90% confidence interval of prediction error and median prediction error (CLMEDPE) are listed in Table 2.

Download:

Table 2. Optimised formula constants for the SRKT, the Hoffer Q, the Holladay 1, Haigis (with optimised a0 and preset values a1 = 0.4 / a2 = 0.1, Haigis1; and with optimised a0 / a1 / a2 constant triplet, Haigis3), and Castrop formula.

https://doi.org/10.1371/journal.pone.0267352.t002

Except for the Castrop formula where all optimisation strategies yielded consistent results, an optimisation solely for the standard deviation or the 90% confidence interval of the prediction error may in some cases yield formula constants which deviate from those derived from an optimisation for sum of squared or sum of absolute prediction error, mean or median prediction error, or combinations of standard deviation and mean or confidence interval and median prediction error.

In the next step, the optimised formula constants for the 8 different metrics were used to derive the formula predicted refraction and to calculate the prediction error in terms of achieved spherical equivalent minus formula predicted refraction. In Table 3, the mean, standard deviation, median, lower and upper boundary of the 90% confidence interval together with the width of the 90% confidence interval, the mean absolute and the root mean squared prediction error are listed for all formulae under test.

Download:

Table 3. Prediction error (PE) as the difference between achieved and formula predicted spherical equivalent for 8 different statistical metrics of formula constant optimisation and various formulae under test.

https://doi.org/10.1371/journal.pone.0267352.t003

The numbers marked in bold in each row indicate the prediction error corresponding to the optimisation metric for that row, and which is therefore expected to yield the lowest value. For example, in the row corresponding to optimisation for the mean or the median prediction error, the resulting mean or median prediction error would be expected to have the smallest value. In general, except when optimising for the standard deviation or the 90% confidence interval without considering the respective mean or median value, the optimisation metrics provide consistent results. When using the standard deviation or the 90% confidence intervals as the sole metrics for formula constant optimisation, the mean / median prediction error may in some cases have a systematic offset up to 0.4 dpt (underlined numbers) / 0.14 dpt (italic numbers). This systematic offset in mean / median prediction error could be mostly compensated for by using metrics for formula constant optimisation which consider combinations of the standard deviation and mean or 90% confidence interval and median.

Fig 1 shows the cumulative distribution function (CDF) of the prediction error (difference between the achieved SEQ from the formula predicted SEQ) for the 6 lens power calculation formulae (SRKT, Hoffer Q, Holladay 1, Haigis with preset a1/a2 and optimisation of a0 (Haigis1) and with optimisation of the constant triplet a0/a1/a2 (Haigis3), and Castrop formula).

Download:

Fig 1. Cumulative density function (CDF) of prediction error (PE, achieved spherical equivalent–formula predicted spherical equivalent) for the SRKT, Hoffer Q, Holladay1, Haigis (with optimised a0 and preset values a1 = 0.4 / a2 = 0.1 (Haigis1); and with optimised a0 / a1 / a2 constant triplet (Haigis3)), and Castrop formula.

Formula constants were optimised to minimise the sum of squared prediction errors, the sum of absolute prediction errors, the mean prediction error, standard deviation (SD) of prediction error, a combination of mean and standard deviation of prediction error (SD & mean), median prediction error, the 90% confidence interval (CL) of prediction error, and a combination of median and 90% confidence interval of prediction error (90% CL & median). The respective formula constants used for calculating the PE are shown in the graphs. Formula constant optimisation for SD or for 90% CL without considering the mean or median may in some formulae lead to a significant shift in the CDF for PE, as indicated by the dashed cyan and dashed magenta lines, whereas the optimisations for sum of squared PE, sum of absolute PE, mean PE, SD & mean PE, and 90% CL & median PE yield consistent results.

https://doi.org/10.1371/journal.pone.0267352.g001

For each formula under test, formula constant(s) were optimised to minimise the sum of squared prediction errors, the sum of absolute prediction errors, the mean prediction error, standard deviation (SD) of prediction error, a combination of mean and standard deviation of prediction error (SD & mean), median prediction error, the 90% confidence interval (CL) of prediction error, and a combination of median and 90% confidence interval of prediction error (90% CL & median). From the first graph (SRKT) we see that the CDFs for the prediction error are almost consistent except for the constant optimisation for the standard deviation, which ends up with an A constant which shifts the patient refraction moderately in the direction of myopia (dashed cyan line shifted to the left). From the second graph (Hoffer Q) we see that the CDFs for the prediction error are mostly consistent except for the constant optimisations for the standard deviation and for the confidence interval, both of which end up with a pACD constant with a moderate to severe shift of the patient refraction in the direction of hyperopia (dashed cyan / dashed magenta line shifted to the right). It can be seen from the third graph (Holladay 1) that, similar to the Hoffer Q formula, the CDFs for the prediction error are mostly consistent except for the constant optimisation for the standard deviation and the confidence interval, both of which end up with a surgeon factor SF which moderately or severely shifts the patient refraction in the direction of hyperopia (dashed cyan / dashed magenta line shifted to the right). The fourth graph (Haigis1) shows that the CDFs for the prediction error are mostly consistent except for the constant optimisation for the standard deviation, which ends up with an a0 which severely shifts the patient refraction in the direction of hyperopia (dashed cyan line shifted to the right). The fifth graph (Haigis3) shows that–similar to the Haigis1 formula–the CDFs for the prediction error are mostly consistent except for the constant optimisation for the standard deviation, which ends up with a a0/a1/a2 constant triplet which severely shifts the patient refraction in the direction hyperopia (dashed cyan line shifted to the right). And last but not least, from the sixth graph (Castrop) we see that the CDFs for the prediction error all are consistent. This means that all of the optimisation metrics tested in this study result in constant triplets C/H/R that do not induce systematic offsets in the formula predicted refraction.

Discussion

Numerous formulae for calculation of intraocular lens power have been proposed in the last 20 years. In contrast to the basic formulae of Fyodorov [11] or Gernet [12] or the classical formulae of Sanders, Retzlaff and Kraff (SRKT), Hoffer (Hoffer-Q), Holladay (Holladay1) or Haigis (simplified Haigis with 1 optimised constant and Haigis formula with 3 optimised constants) [1, 4–9, 14], most of the formula authors nowadays do not disclose or publish the calculation strategy. At best they offer WEB based applications or software solutions for calculating the lenses. Such software tools do not allow batch calculations on a large set of patient data. Today, classical formulae are increasingly being replaced by ‘modern’ calculation strategies such as the Barrett Universal II, Kane, Pearl, EVO, VRF/VRF-G, Hill RBF or T2 formula in many countries of the world [2, 16]. To compare the prediction performance with other formulae it is necessary to enter the data from preoperative data (biometry), intraoperative data (lens power) and postoperative data (manual refraction after 4 weeks to 6 months) manually, introducing a large risk of transcription errors. Additionally, a systematic optimisation of constants is not possible for undisclosed formulae [15, 19, 20]. Obtaining the appropriate formula constant is mostly achieved by trial and error e.g. by varying the constants in the calculation scheme to eliminate the mean or median prediction error.

Currently there are no generally accepted guidelines or recommendations for formula constant optimisation [15, 20, 29]. Where a formula is fully disclosed, and in the simple case of a formula with a single constant, it would be possible to reorganise the formula to solve for the formula constant (for each clinical case), and from the distribution of the formula constants the clinician could select statistical metrics such as the arithmetic or geometric mean or the median. However, back-calculating the prediction error with such constants does not necessarily yield the best formula performance in terms of least root mean squared or least mean absolute error or zeroing of the mean or median of the prediction error [19]. In the case of undisclosed formulae or formulae with more than 1 constant there is no straightforward technique for back-tracing the appropriate formula constant(s) for an individual case. In the best case, when all formula constants together describe a linear superposition of a parameter for which the formula could be solved, (e.g. the effective lens position ‘d’ in the Haigis formula) a linear regression technique could be used to derive the regression coefficients (for the Haigis formula a0 as intercept, and a1 / a2 as coefficients for the phakic ACD / AL) for the dataset [19]. In some studies, when investigating the performance of a formula and comparing to other formulae, the formula constants are not optimised at all. Instead, the prediction error calculated with a given preset formula constant is zeroed for its mean or median [15, 20, 30] to analyse performance metrics such as mean absolute or root mean squared prediction error, the standard deviation, or confidence intervals of the distribution. However, such techniques do not yield proper results for the formula performance, as statistical metrics such as standard deviation or confidence intervals can change depending on the formula constant used for predicting the refractive outcome. When zeroing the mean or median prediction error, the entire distribution of the prediction error is simply shifted to centre it on a zero prediction error, which is rather different to the distribution of the prediction error when the formula constant is optimised.

A systematic optimisation of formula constants for statistical metrics of the refractive outcome requires more advanced techniques than simple reorganisation of the formula to solve for the formula constant [18]. In engineering, nonlinear iterative techniques are typically used to derive the best formula constant(s) in terms of minimising any classical statistical metrics of the prediction error. The requirements for such nonlinear techniques are: reliable and fast convergence with all datasets after a couple of iterations, robustness in a way that they result in global, rather than local, minima, and a straightforward implementation with a programming code.

In the present paper, we have used classical statistical metrics to show the capability of nonlinear iterative optimisation strategies: the sum of squared prediction error, the sum of absolute prediction error, the mean, the standard deviation of the prediction error distribution, a combination of the standard deviation of the prediction error distribution and the mean prediction error, the median prediction error, the 90% confidence interval of the prediction error, and a combination of the 90% confidence interval of the prediction error and the median prediction error. In general, such nonlinear iterative optimisation strategies are not restricted to classical metrics [28], and any continuous metrics could be used for the merit function. For statistical metrics, sum of squares and sum of absolute prediction error, special techniques which have a fast and reliable convergence to a global minimum of the merit function have been developed. Especially for minimising the sum of squared prediction error, which is equivalent to the minimisation of the root mean squared prediction error, the Levenberg-Marquardt algorithm [23, 24] is well-known in engineering disciplines to show excellent performance and robustness. For minimisation of other statistical metrics, other optimisation techniques such as the interior point methods [26, 27] or the simplex algorithm [25] are used. However, in general, if the merit function to be minimised is a ‘good-natured problem’ with a global convexity [25, 28] and without local minima, any iterative optimisation algorithm will result in a good solution for the formula constant (but with different numbers of iteration cycles and function evaluations). For this study we used the Levenberg-Marquardt algorithm [23, 24] for optimising the constant(s) for the root mean squared and the mean absolute error, and for all other metrics we used the interior point method [26, 27]. The iterative algorithms under test showed stable convergence for all formulae and for all optimisation metrics.

We see from our results that formula constant optimisation for the root mean squared and mean absolute prediction error, and also for the mean and median prediction error yields mostly consistent results for all formulae. Fig 1 shows us that the cumulative distribution functions (CDFs) as a measure of the formula performance with this optimisation are quite coherent, and there are no systematic differences in the shape or the horizontal position of the CDFs. However, looking more in detail at the respective formula constants for the Haigis3 and the Castrop formula, both of which use constant triplets, it is obvious that the individual constants in the triplet can vary significantly (e.g. for the Castrop formula: C/H/R with optimising for SoSPE 0.2814/0.3500/0.0848 vs. 0.2517/0.5014/0.0554 with optimising for SoAPE) without significantly affecting the distributions of the prediction error [19]. All of the constant triplets work quite well when calculating for eyes with ‘normal’ biometrical parameters, but if applying these constant triplets to extreme values of biometric parameters (highly myopic/hyperopic eyes or with extremely flat or steep corneal curvatures), the resulting formula predicted refraction can vary significantly between the different constant triplets.

What can also be seen from Fig 1 is that the standard deviation of the prediction error distribution and the 90% confidence interval of the prediction error as statistical metrics for formula constant optimisation do not yield the best results in the overall performance. In particular, optimising for the standard deviation of the prediction error without considering the mean prediction error shows a systematic offset in the CDF graphs, indicating that with the formula optimised constants having the smallest mean prediction errors, the standard deviation of the prediction error is not smallest. Only with the Castrop formula are all of the optimisation strategies more or less consistent. This is mostly due to the internal structure of the formula, where one of the formula constants (R) is used as a systematic offset value for the formula predicted refraction, mainly to account for different lane distances in the refraction measurements. With this constant the CDF may be centred using the constant optimisation strategy without changing the ‘shape’ of the CDF. However, optimising formula constants for the 90% confidence interval does not yield acceptable results for all of the formulae under test. In particular, the CDF for the Hoffer Q and the Holladay 1 formula is not centred and indicates some systematic offset. This means that optimising for the confidence interval without considering the median prediction error may provide constant(s) with a narrow distribution, but probably with a systematic offset in the prediction error. To overcome this problem we completed our optimisation strategies using a compromise of using combinations of the standard deviation of the prediction error distribution and the mean prediction error to ensure that the CDF is centred and the standard deviation is considered in the optimisation (optimising for SDMPE). Accordingly, we completed our optimisation strategies using a compromise of using combinations of the 90% confidence interval of the prediction error and the median prediction error to ensure that the CDF is centred and the confidence interval is considered in the optimisation (optimising for CLMEDPE). Both optimisations show a good performance in terms of low systematic offset of the prediction error and a narrow CDF.

The results of the classical performance indicators for lens power calculation concepts are summarised in Table 3. The mean error mostly accounts for a systematic shift towards either myopia or hyperopia in the formula predicted refraction, and the standard deviation accounts for the width of the distribution if a (parametric) normal distribution is considered. The median and the lower and upper boundaries of the 90% confidence interval (together with the width of the 90% confidence interval) refer to the (nonparametric) metrics of an arbitrary distribution for the prediction error. Finally, the mean absolute and the root mean squared prediction error are classical measures in the literature to evaluate the performance of formula predicted refractions, but both are significantly affected by systematic offsets in the prediction error. This means that they are not proper measures of the consistency of the results. The listing implies that the optimisation strategies worked quite well: where we optimised for the SoSPE the resulting root mean squared error (RMS) is the lowest in the row for all formulae. When optimising for the SoAPE it is clear that on our dataset the mean absolute error (ABS) is the lowest in the row. Where we have optimised for the mean prediction error (MPE) or the median prediction error (MEDPE), it is can be seen that the resulting mean prediction error (MEAN) or the median prediction error (MEDIAN) is lowest in the row. When optimising for the standard deviation of the prediction error distribution (SDPE) we see that the standard deviation (SD) is lowest in the row, but at the cost of a systematic shift in mean or median formula predicted refraction (MEAN or MEDIAN). Accordingly, when optimising for the 90% confidence interval, the width of the confidence interval (90% CL) is lowest in the row at the cost of the mean or median value (MEAN and MEDIAN) which are systematically shifted, at least in the Hoffer Q and the Holladay formula.

In conclusion, this study describes techniques for optimising formula constants on a clinical dataset using modern nonlinear iterative optimisation strategies. These optimisation strategies are generally capable of dealing with all disclosed lens power calculation formulae having one or more formula constants, and with any statistical metric of any outcome measure. The applicability of this constant optimisation technique is shown for the SRKT, Hoffer Q, Holladay 1, Haigis, and Castrop formula on a clinical dataset with 888 data from 1 clinical centre treated with the Hoya Vivinex intraocular lens. The prediction error, defined as the difference between the postoperatively achieved and the formula predicted spherical equivalent power of refraction, was used as a merit function. Optimisations for the root mean squared and mean absolute prediction error and for the mean and median prediction error yielded mostly consistent performance results, but optimisations for the standard deviation of the prediction error distribution or the confidence interval may result in a narrow distribution but with some systematic offset in the formula predicted refraction which is not appropriate in a clinical setting.

References

1. Shammas HJ. Intraocular lens power calculations. Slack Inc 2004. ISBN-13: 978–1556426520.
2. Melles RB, Kane JX, Olsen T, Chang WJ. Update on intraocular lens calculation formulae. Ophthalmology 2019; 126(9):1334–1335. pmid:30980854
- View Article
- PubMed/NCBI
- Google Scholar
3. Sanders DR, Retzlaff JA, Kraff MC, Gimbel HV, Raanan MG. Comparison of the SRK/T formula and other theoretical and regression formulas. J Cataract Refract Surg 1990; 16(3):341–346. pmid:2355322
- View Article
- PubMed/NCBI
- Google Scholar
4. Savini G, Taroni L, Hoffer KJ. Recent developments in intraocular lens power calculation methods-update 2020. Ann Transl Med. 2020;8(22):1553. pmid:33313298
- View Article
- PubMed/NCBI
- Google Scholar
5. Retzlaff JA, Sanders DR, Kraff MC. Development of the SRK/T intraocular lens implant power calculation formula. J Cataract Refract Surg 1990; 16(3):333–340. pmid:2355321
- View Article
- PubMed/NCBI
- Google Scholar
6. Hoffer KJ. Steps for IOL power calculation. Am Intraocul Implant Soc 1980; 6(4):370. pmid:7440385.
- View Article
- PubMed/NCBI
- Google Scholar
7. Hoffer KJ. Intraocular lens calculation: the problem of the short eye. Ophthalmic Surg 1981; 12(4):269–272. pmid:7254770.
- View Article
- PubMed/NCBI
- Google Scholar
8. Hoffer KJ. The Hoffer Q formula: a comparison of theoretic and regression formulas. J Cataract Refract Surg 1993; 19(6):700–712. pmid:8271165
- View Article
- PubMed/NCBI
- Google Scholar
9. Holladay JT, Prager TC, Chandler TY, Musgrove KH, Lewis JW, Ruiz RS. A three-part system for refining intraocular lens power calculations. J Cataract Refract Surg 1988; 14(1):17–24. pmid:3339543
- View Article
- PubMed/NCBI
- Google Scholar
10. Gernet H, Ostholt H, Werner H. Die präoperative Berechnung intraocularer Binkhorst-Linsen. 122. Vers. d. Ver. Rhein.-Westfäl. Augenärzte. Balve, Verlag Zimmermann 1970: pp. 54–55.
- View Article
- Google Scholar
11. Fyodorov SN, Galin MA, Linksz A. Calculation of the optical power of intraocular lenses. Invest Ophthalmol 1975; 14: 625–628. pmid:1150402
- View Article
- PubMed/NCBI
- Google Scholar
12. Norrby NE, Koranyi G. Prediction of intraocular lens power using the lens haptic plane concept. J Cataract Refract Surg 1997; 23(2):254–259. pmid:9113578
- View Article
- PubMed/NCBI
- Google Scholar
13. Olsen T. J Prediction of the effective postoperative (intraocular lens) anterior chamber depth. J Cataract Refract Surg 2006; 32(3):419–424. pmid:16631049
- View Article
- PubMed/NCBI
- Google Scholar
14. Olsen T, Hoffmann P. C constant: new concept for ray tracing-assisted intraocular lens power calculation. J Cataract Refract Surg 2014; 40(5):764–773. pmid:24767910
- View Article
- PubMed/NCBI
- Google Scholar
15. Aristodemou P, Knox Cartwright NE, Sparrow JM, Johnston RL. Intraocular lens formula constant optimization and partial coherence interferometry biometry: Refractive outcomes in 8108 eyes after cataract surgery. J Cataract Refract Surg 2011; 37(1):50–62. pmid:21183099
- View Article
- PubMed/NCBI
- Google Scholar
16. Wendelstein J, Hoffmann P, Hirnschall N, et al. Project hyperopic power prediction: accuracy of 13 different concepts for intraocular lens calculation in short eyes. Br J Ophthalmol 2021 27: bjophthalmol-2020-318272. pmid:33504489
- View Article
- PubMed/NCBI
- Google Scholar
17. Vega Y, Gershoni A, Achiron A, Tuuminen R, Weinberger Y, Livny E, et al. High agreement between Barrett Universal II calculations with and without Utilization of optional biometry parameters. J Clin Med. 2021;10(3):542. pmid:33540639
- View Article
- PubMed/NCBI
- Google Scholar
18. Langenbucher A, Szentmáry N, Cayless A, Weisensee J, Fabian E, Wendelstein J, et al. Considerations on the Castrop formula for calculation of intraocular lens power. PLoS One. 2021 Jun 2;16(6):e0252102. pmid:34077432
- View Article
- PubMed/NCBI
- Google Scholar
19. Schröder S, Leydolt C, Menapace R, Eppig T, Langenbucher A. Determination of Personalized IOL-Constants for the Haigis Formula under consideration of measurement precision. PLoS One 2016; 11(7):e0158988. pmid:27391100
- View Article
- PubMed/NCBI
- Google Scholar
20. Zhang JQ, Zou XY, Zheng DY, Chen WR, Sun A, Luo LX. Effect of lens constants optimization on the accuracy of intraocular lens power calculation formulas for highly myopic eyes. Int J Ophthalmol 2019; 12(6):943–948. pmid:31236350
- View Article
- PubMed/NCBI
- Google Scholar
21. Cooke DL, Cooke TL. Approximating sum-of-segments axial length from a traditional optical low-coherence reflectometry measurement. J Cataract Refract Surg 2019; 45(3):351–354. pmid:30851808
- View Article
- PubMed/NCBI
- Google Scholar
22. Cooke DL, Cooke TL. A comparison of two methods to calculate axial length. J Cataract Refract Surg 2019 45(3):284–292. pmid:30851805
- View Article
- PubMed/NCBI
- Google Scholar
23. Levenberg K. A method for the solution of certain problems in least squares. Quart Appl Math 1944; 2:164–168.
- View Article
- Google Scholar
24. Marquardt D. An algorithm for least-squares estimation of nonlinear parameters. SIAM J Appl Math 1963; 11:431–441.
- View Article
- Google Scholar
25. Kanzow C, Yamashita N, Fukushima M. Levenberg–Marquardt methods with strong local convergence properties for solving nonlinear equations with convex constraints. J Comp Applied Mathematics 2004. 172 (2): 375–397.
- View Article
- Google Scholar
26. Dikin II. Iterative solution of problems of linear and quadratic programming. Dokl. Akad 1967. Nauk SSSR. 174 (1): 747–748.
- View Article
- Google Scholar
27. Karmarkar N. A new polynomial-time algorithm for linear programming. Proceedings of the sixteenth annual ACM symposium on Theory of computing–STOC 1984. p. 302. https://doi.org/10.1145/800057.808695
28. Boyd S, Vandenberghe L. Convex Optimization. Cambridge: Cambridge University Press 2004. p. 143. ISBN 978-0-521-83378-3. MR 2061575.
29. Langenbucher A, Szentmáry N, Cayless A, Müller M, Eppig T, Schröder S, et al. IOL formula constants: strategies for optimization and defining standards for presenting data. Ophthalmic Res. 2021;64(6):1055–1067. pmid:33530082
- View Article
- PubMed/NCBI
- Google Scholar
30. Hoffer KJ, Savini G. Update on Intraocular Lens Power Calculation Study Protocols: The Better Way to Design and Report Clinical Trials. Ophthalmology. 2020 9:S0161-6420(20)30638-2. pmid:32653457
- View Article
- PubMed/NCBI
- Google Scholar

[ref1] 1. Shammas HJ. Intraocular lens power calculations. Slack Inc 2004. ISBN-13: 978–1556426520.

[ref2] 2. Melles RB, Kane JX, Olsen T, Chang WJ. Update on intraocular lens calculation formulae. Ophthalmology 2019; 126(9):1334–1335. pmid:30980854
View Article
PubMed/NCBI
Google Scholar

[3] View Article

[4] PubMed/NCBI

[5] Google Scholar

[ref3] 3. Sanders DR, Retzlaff JA, Kraff MC, Gimbel HV, Raanan MG. Comparison of the SRK/T formula and other theoretical and regression formulas. J Cataract Refract Surg 1990; 16(3):341–346. pmid:2355322
View Article
PubMed/NCBI
Google Scholar

[7] View Article

[8] PubMed/NCBI

[9] Google Scholar

[ref4] 4. Savini G, Taroni L, Hoffer KJ. Recent developments in intraocular lens power calculation methods-update 2020. Ann Transl Med. 2020;8(22):1553. pmid:33313298
View Article
PubMed/NCBI
Google Scholar

[11] View Article

[12] PubMed/NCBI

[13] Google Scholar

[ref5] 5. Retzlaff JA, Sanders DR, Kraff MC. Development of the SRK/T intraocular lens implant power calculation formula. J Cataract Refract Surg 1990; 16(3):333–340. pmid:2355321
View Article
PubMed/NCBI
Google Scholar

[15] View Article

[16] PubMed/NCBI

[17] Google Scholar

[ref6] 6. Hoffer KJ. Steps for IOL power calculation. Am Intraocul Implant Soc 1980; 6(4):370. pmid:7440385.
View Article
PubMed/NCBI
Google Scholar

[19] View Article

[20] PubMed/NCBI

[21] Google Scholar

[ref7] 7. Hoffer KJ. Intraocular lens calculation: the problem of the short eye. Ophthalmic Surg 1981; 12(4):269–272. pmid:7254770.
View Article
PubMed/NCBI
Google Scholar

[23] View Article

[24] PubMed/NCBI

[25] Google Scholar

[ref8] 8. Hoffer KJ. The Hoffer Q formula: a comparison of theoretic and regression formulas. J Cataract Refract Surg 1993; 19(6):700–712. pmid:8271165
View Article
PubMed/NCBI
Google Scholar

[27] View Article

[28] PubMed/NCBI

[29] Google Scholar

[ref9] 9. Holladay JT, Prager TC, Chandler TY, Musgrove KH, Lewis JW, Ruiz RS. A three-part system for refining intraocular lens power calculations. J Cataract Refract Surg 1988; 14(1):17–24. pmid:3339543
View Article
PubMed/NCBI
Google Scholar

[31] View Article

[32] PubMed/NCBI

[33] Google Scholar

[ref10] 10. Gernet H, Ostholt H, Werner H. Die präoperative Berechnung intraocularer Binkhorst-Linsen. 122. Vers. d. Ver. Rhein.-Westfäl. Augenärzte. Balve, Verlag Zimmermann 1970: pp. 54–55.
View Article
Google Scholar

[35] View Article

[36] Google Scholar

[ref11] 11. Fyodorov SN, Galin MA, Linksz A. Calculation of the optical power of intraocular lenses. Invest Ophthalmol 1975; 14: 625–628. pmid:1150402
View Article
PubMed/NCBI
Google Scholar

[38] View Article

[39] PubMed/NCBI

[40] Google Scholar

[ref12] 12. Norrby NE, Koranyi G. Prediction of intraocular lens power using the lens haptic plane concept. J Cataract Refract Surg 1997; 23(2):254–259. pmid:9113578
View Article
PubMed/NCBI
Google Scholar

[42] View Article

[43] PubMed/NCBI

[44] Google Scholar

[ref13] 13. Olsen T. J Prediction of the effective postoperative (intraocular lens) anterior chamber depth. J Cataract Refract Surg 2006; 32(3):419–424. pmid:16631049
View Article
PubMed/NCBI
Google Scholar

[46] View Article

[47] PubMed/NCBI

[48] Google Scholar

[ref14] 14. Olsen T, Hoffmann P. C constant: new concept for ray tracing-assisted intraocular lens power calculation. J Cataract Refract Surg 2014; 40(5):764–773. pmid:24767910
View Article
PubMed/NCBI
Google Scholar

[50] View Article

[51] PubMed/NCBI

[52] Google Scholar

[ref15] 15. Aristodemou P, Knox Cartwright NE, Sparrow JM, Johnston RL. Intraocular lens formula constant optimization and partial coherence interferometry biometry: Refractive outcomes in 8108 eyes after cataract surgery. J Cataract Refract Surg 2011; 37(1):50–62. pmid:21183099
View Article
PubMed/NCBI
Google Scholar

[54] View Article

[55] PubMed/NCBI

[56] Google Scholar

[ref16] 16. Wendelstein J, Hoffmann P, Hirnschall N, et al. Project hyperopic power prediction: accuracy of 13 different concepts for intraocular lens calculation in short eyes. Br J Ophthalmol 2021 27: bjophthalmol-2020-318272. pmid:33504489
View Article
PubMed/NCBI
Google Scholar

[58] View Article

[59] PubMed/NCBI

[60] Google Scholar

[ref17] 17. Vega Y, Gershoni A, Achiron A, Tuuminen R, Weinberger Y, Livny E, et al. High agreement between Barrett Universal II calculations with and without Utilization of optional biometry parameters. J Clin Med. 2021;10(3):542. pmid:33540639
View Article
PubMed/NCBI
Google Scholar

[62] View Article

[63] PubMed/NCBI

[64] Google Scholar

[ref18] 18. Langenbucher A, Szentmáry N, Cayless A, Weisensee J, Fabian E, Wendelstein J, et al. Considerations on the Castrop formula for calculation of intraocular lens power. PLoS One. 2021 Jun 2;16(6):e0252102. pmid:34077432
View Article
PubMed/NCBI
Google Scholar

[66] View Article

[67] PubMed/NCBI

[68] Google Scholar

[ref19] 19. Schröder S, Leydolt C, Menapace R, Eppig T, Langenbucher A. Determination of Personalized IOL-Constants for the Haigis Formula under consideration of measurement precision. PLoS One 2016; 11(7):e0158988. pmid:27391100
View Article
PubMed/NCBI
Google Scholar

[70] View Article

[71] PubMed/NCBI

[72] Google Scholar

[ref20] 20. Zhang JQ, Zou XY, Zheng DY, Chen WR, Sun A, Luo LX. Effect of lens constants optimization on the accuracy of intraocular lens power calculation formulas for highly myopic eyes. Int J Ophthalmol 2019; 12(6):943–948. pmid:31236350
View Article
PubMed/NCBI
Google Scholar

[74] View Article

[75] PubMed/NCBI

[76] Google Scholar

[ref21] 21. Cooke DL, Cooke TL. Approximating sum-of-segments axial length from a traditional optical low-coherence reflectometry measurement. J Cataract Refract Surg 2019; 45(3):351–354. pmid:30851808
View Article
PubMed/NCBI
Google Scholar

[78] View Article

[79] PubMed/NCBI

[80] Google Scholar

[ref22] 22. Cooke DL, Cooke TL. A comparison of two methods to calculate axial length. J Cataract Refract Surg 2019 45(3):284–292. pmid:30851805
View Article
PubMed/NCBI
Google Scholar

[82] View Article

[83] PubMed/NCBI

[84] Google Scholar

[ref23] 23. Levenberg K. A method for the solution of certain problems in least squares. Quart Appl Math 1944; 2:164–168.
View Article
Google Scholar

[86] View Article

[87] Google Scholar

[ref24] 24. Marquardt D. An algorithm for least-squares estimation of nonlinear parameters. SIAM J Appl Math 1963; 11:431–441.
View Article
Google Scholar

[89] View Article

[90] Google Scholar

[ref25] 25. Kanzow C, Yamashita N, Fukushima M. Levenberg–Marquardt methods with strong local convergence properties for solving nonlinear equations with convex constraints. J Comp Applied Mathematics 2004. 172 (2): 375–397.
View Article
Google Scholar

[92] View Article

[93] Google Scholar

[ref26] 26. Dikin II. Iterative solution of problems of linear and quadratic programming. Dokl. Akad 1967. Nauk SSSR. 174 (1): 747–748.
View Article
Google Scholar

[95] View Article

[96] Google Scholar

[ref27] 27. Karmarkar N. A new polynomial-time algorithm for linear programming. Proceedings of the sixteenth annual ACM symposium on Theory of computing–STOC 1984. p. 302. https://doi.org/10.1145/800057.808695

[ref28] 28. Boyd S, Vandenberghe L. Convex Optimization. Cambridge: Cambridge University Press 2004. p. 143. ISBN 978-0-521-83378-3. MR 2061575.

[ref29] 29. Langenbucher A, Szentmáry N, Cayless A, Müller M, Eppig T, Schröder S, et al. IOL formula constants: strategies for optimization and defining standards for presenting data. Ophthalmic Res. 2021;64(6):1055–1067. pmid:33530082
View Article
PubMed/NCBI
Google Scholar

[100] View Article

[101] PubMed/NCBI

[102] Google Scholar

[ref30] 30. Hoffer KJ, Savini G. Update on Intraocular Lens Power Calculation Study Protocols: The Better Way to Design and Report Clinical Trials. Ophthalmology. 2020 9:S0161-6420(20)30638-2. pmid:32653457
View Article
PubMed/NCBI
Google Scholar

[104] View Article

[105] PubMed/NCBI

[106] Google Scholar

Figures

Abstract

Background

Methods

Results

Conclusion

Introduction

Materials and methods

Dataset for formula constant optimisation

Preprocessing of the data

Formula constant optimisation

Results

Discussion

References