#### Abstract

To demonstrate the importance of material properties of the cornea in intraocular pressure (IOP) readings via standard Goldmann applanation tonometry.

A realistic finite element model of the cornea was developed for the simulation of Goldmann applanation tonometry. A virtual cornea population was generated by randomly sampling material properties, central corneal thickness (CCT), and IOP for comparison with 181 clinical cases. The effect of material properties and CCT on IOP prediction in the virtual population was determined via computational simulation.

The results show that corneal biomechanical properties (as characterized in this study by the stiffness parameter E_{init}) are as important as the CCT in influencing measured (Goldmann) IOP.

This study supports the contention that the observed large scatter in standard correlations of clinical measurements of IOP versus CCT can be largely accounted for by plausible individual variations in corneal biomechanical stiffness properties.

From ESI-USA R&D, Champaign (Kwon); and the Department of Civil and Environmental Engineering, University of Illinois at Urbana-Champaign (Ghaboussi, Pecknold, Hashash), Illinois.

This research was partially funded with a grant from the Grainger Foundation, Urbana, Illinois (Ghaboussi).

The authors have no proprietary or financial interest in the materials presented herein.

**AUTHOR CONTRIBUTIONS**

*Study concept and design (T.K., J.G., D.A.P., Y.H.); data collection (T.K., J.G.); analysis and interpretation of data (T.K., J.G., D.A.P.); drafting of the manuscript (T.K., J.G., D.A.P.); critical revision of the manuscript (T.K., J.G., Y.H.); statistical expertise (J.G.); obtained funding (J.G.); administrative, technical, or material support (T.K., J.G., D.A.P., Y.H.)*

Correspondence: Tae-Hyun Kwon, PhD, ESI-USA R&D, 1834 E Amber Ln, #204, Urbana, IL 61802. Tel: 217.419.0182; Fax: 217.344.2619; E-mail: taekwonkr@gmail.com

Measurement of intraocular pressure (IOP) in the human eye is the primary screening method for early detection of those at risk for glaucoma. The most common method for IOP measurement is applanation (Goldmann) tonometry, in which an instrument is used to flatten the central part of the cornea. A number of studies^{1–4} have been carried out on human eyes in vivo. Ehlers et al^{2} established linear correlations between corneal thickness and the error of applanation tonometry. Johnson et al^{3} reported a special case of one patient with 0.9-mm corneal thickness. The patient’s actual IOP, later determined by cannulation, was 11 mmHg, whereas her applanation tonometer reading was 35 mmHg. There have been numerous studies to determine the relationship between central corneal thickness (CCT) and IOP measurement in normal, glaucomatous, and ocular hypertensive eyes. Whitacre et al^{4} also pointed out that measuring the corneal thickness is necessary to properly interpret the results of Goldmann applanation tonometry. Doughty and Zaman,^{1} in a thorough analysis of the literature published from 1968 to 1999, found a statistically significant correlation between apparent IOP and CCT.

All of these studies show systematic errors in IOP tonometer readings that correlate with increased CCT. Such errors have significant implications for glaucoma screening. However, all of the results show a high degree of scatter. The purpose of this study is to demonstrate, via a series of finite element computational simulations, that the variations in other biomechanical characteristics and properties of the cornea, that are not controlled in these various studies, are likely responsible for the observed scatter.

This study uses an improved geometrical and biomechanical model of the cornea. The corneal biomechanical model is first calibrated using experimental data.^{5} Corneal biomechanical and geometrical parameters and true IOPs are then randomly sampled within ranges established from published data to produce a virtual population of corneas. Numerical finite element simulation of Goldmann applanation tonometry is then carried out on each member of the virtual population to determine (simulated) Goldmann IOP. Comparisons are made with the true IOPs to determine the effect of the corneal biomechanical and geometric properties on errors in IOP.

### Materials and Methods

A finite element model of the human cornea, including the effect of the aqueous humor in the anterior chamber, is used in this study. The model is briefly described here, and differences with previous simplified analytical models are noted. A more detailed description is available elsewhere.^{6,7}

#### Corneal Geometry

The corneal surfaces are assumed to be axisymmetric (with axis of symmetry along the optical axis), but not necessarily spherical. The geometry of the anterior and posterior surfaces is expressed by the following equation^{8}:

^{8}

Once the anterior and posterior surfaces are defined, the aspheric finite element model of the cornea is developed using 4-noded axisymmetric solid elements. Depending on the parameters describing the anterior and posterior surfaces, the thickness of the cornea may vary from the optical axis to the limbus. Previous studies, particularly Liu and Roberts^{9} and Orssengo and Pye,^{10} approximated the cornea as a thin spherical shell of constant thickness.

#### Corneal Biomechanical Properties

The nonlinear constitutive behavior of the cornea has been documented in the literature.^{11–14} Anderson et al^{15} conducted inflation tests on porcine corneas and showed significant stiffening in a collagen-regulated phase. Bryant and McDonnell^{5} also noted that, in each of five representative corneas, the knee of inflation curves occurred around the normal IOP (12 to 20 mmHg).

In this study, the cornea is modeled as a nonlinear transversely isotropic material, in which the elastic properties depend on the stress level in the cornea, stiffening as stresses increase. The constitutive model is derived from Fung’s exponential strain energy function W as^{16}

_{0}is the initial (linear) elastic tangent stiffness matrix. The constitutive relation (secant stiffness D

_{s}, which relates 2nd Piola-Kirchhoff stress σ to Green strain ɛ) is obtained by differentiating the strain energy function with respect to strain.

A second differentiation yields the Jacobian matrix D_{t}.

The initial tangent stiffness matrix D_{0} contains five independent elastic constants for transverse isotropy; these have been reduced to a single independent elastic constant E_{init} (the initial Young’s modulus in the plane of isotropy) by using correlations proposed for the human cornea by Bryant and McDonnell.^{5} Then, the initial in-plane Young’s modulus E_{init} and the nonlinearity parameter α are determined by calibration with in vitro inflation tests.^{5} In Bryant and McDonnell’s in vitro inflation tests, each cornea in turn was mounted on an artificial anterior chamber and its apical displacement was measured during incremental pressure increases.^{5} Separate inflation tests are simulated with the finite element model of the cornea. Calibrated finite element results are shown in Figure 1 and compared with the experimental results for five of these tests.

Figure 1. Calibration of Nonlinear Material Model in Finite Element Simulations (lines) with the Experimental Results (circles) by Bryant and McDonnell. |

Note that there are substantial variations between the inflation responses of the different individual corneas. The calibrated values of E_{init} range from 0.12 to 0.35 MPa for the five tests; calibrated values of α ranged from 417 to 1000 (1/MPa). The results also appear to indicate that the greater the initial slope of the pressure-displacement curve, the more rapidly the slope continues to increase with increasing pressure. Therefore, it is found that the material constant α, which controls the extent of the nonlinearity, can be related to E_{init} by a linear interpolation function with sufficient accuracy.^{6} The resulting model, which is a function only of the single material constant E_{init}, is referred to herein as the “Extended Fung Model.”

Previous cornea models have often assumed that the cornea material is linear elastic (no stress dependence) and isotropic. In this regard, it should be noted that the constant Young’s modulus used in previous models of this type represents an average value over the applanation process, starting from the (initial) physiological stress level in the cornea corresponding to the IOP, and is not directly comparable to the parameter E_{init} used in the current model, which represents the stiffness of the cornea at zero physiological stress. Thus, the constant Young’s modulus values used previously in linear elastic cornea models are much higher than the parameter E_{init} used in the model described herein. At higher IOP values, and thus higher physiological stress, corneal stiffness is greater in response to applanating force.

#### Initial Physiological Stress Due to IOP

Equation () describes the geometry of the cornea when it is already in a state of initial stress produced by the IOP. To take account of this initial stress state in the simulation of applanation tonometry, a two-step numerical procedure is used herein.

In the first step of the computational procedure, the IOP is applied using the initial linear properties (D_{0}), which are artificially increased by a suitable factor (a factor of 100 is used). The purpose of this first computational step is to determine a set of in situ stresses σ_{1} in the cornea, which equilibrate the IOP, without producing appreciable deformations. The deformed shape of the cornea at the end of the first computational step is therefore close to that described by Eq.(2) as desired.

Fictitious strains ɛ_{1} in this first step are computed using the artificial initial stiffness. The strains ɛ_{1}^{iter.} that would be consistent with the actual nonlinear constitutive model and the computed initial stresses σ_{1} are determined at each material reference point via an iterative (Newton-Raphson) numerical process.

The applanation of the cornea is then modeled in the second computational step via a series of displacement controlled load increments, now using the actual corneal constitutive model, with total strains corrected by the strains ɛ_{1} and ɛ_{1}^{iter.} determined from the first step as summarized in the following pseudo code:

If step equals 1, then

- σ
_{1}*= kD*ɛ_{0}, where k = 100_{1} - store (ɛ
_{1}, σ_{1}) at the end of the first step - σ
_{1}→ ɛ_{1}^{iter.}by Newton-Raphson iteration.

- ɛ
_{n}_{+}_{1}*=*ɛ+ɛ_{n}−ɛ_{1}^{iter.}_{1} - σ
_{n}_{+}_{1}*= D*ɛ_{s}(_{n}_{+}ɛ_{1})_{n}_{+}_{1} *d*σ_{n}_{+}_{1}*= D*ɛ_{t}(_{n}_{+}ɛ_{1})d_{n}_{+}._{1}

#### Finite Element Simulation of Goldmann Applanation Tonometry

In the full three-dimensional finite element model of the cornea^{6} that included the anterior and posterior chambers separated by the lens, the interaction between the aqueous humor and cornea was simulated with hydrostatic fluid elements. The aqueous humor was assumed to be incompressible with fluid mass density of 1000 kg/m^{3}. The fluid elements provided the coupling between the deformation of the cornea and the pressure exerted by the aqueous humor on the anterior surface of the cornea. Aqueous humor secretion was simulated with the fluid link element by means of a functional relationship between the aqueous humor outflow rate (q) and the pressure difference (Δp) across the fluid link element, defined by *q* = (*0.28μl*/min/*mmHg*) ×Δ*p*. Applanation simulations showed that the IOP was only slightly elevated above its initial value due to displacement of aqueous humor during the applanation. Because the IOP elevation was negligible in practical terms as Goldmann had postulated, the full three-dimensional model is simplified by using the current axisymmetric cornea where the IOP is applied through fluid elements without modeling of aqueous humor secretion.

Schwartz et al^{17} estimated the pressure drop across the tear meniscus between the cornea and the edge of the applanator. The following calculation is based on their equations. The pressure drop (Δd) across the tear surface is illustrated schematically in Figure 2.

*R*is the radius of the fluid surface, 1.52 mm (applanation radius); and

_{t}*r*is the tear radius at the edge of the wetted area.

_{t}*R*is the radius of curvature of the anterior surface of the cornea, taken as 7.77 mm from the corneal geometry.

_{c}Then,

The calculated pressure drop (=4.7 mmHg) is assumed to be the same for all corneas modeled in this study because the corneal radii of curvature are not considered as variables. The applanation procedure is simulated by contact analysis with a rigid surface representing the applanator. The applied force (expressed as a pressure) is calculated when the cornea is applanated over the contact radius (1.52 mm) and the pressure drop 4.7 mmHg is subtracted to obtain the (simulated) measured IOPs. Unless otherwise stated, the values of simulated IOP in this study are corrected by the pressure drop Δd=4.7 mmHg to take into account the effect of the tear meniscus. Additional details of the modeling procedure can be found elsewhere.^{6,7}

This contact analysis results in the central portion of the cornea being flattened as required. Previous models^{9,10} that have used a closed form (thin shell) analytical expression for the corneal deflection under applanation implicitly assume that the contact pressure between applanator and cornea is uniform. Velten et al^{18} showed that the assumption of uniform contact stress produced significantly different results from a displacement controlled contact analysis. The final deformed shape of an average cornea with CCT of 550 μm is shown in Figure 3.

Figure 3. Deformed Shape of Applanated Cornea with Undeformed Shape Superimposed (dotted Lines). Central Corneal Thickness = 550 |

#### Virtual Population of Cornea Finite Element Models

A virtual population of cornea finite element models is generated for simulation. It is intended that this virtual population be roughly comparable to the 181 cases reported in a clinical study by Bhan et al,^{19} with a reasonably similar scatter in the correlation of the apparent IOP with CCT. The virtual population of 181 different cornea models is generated by randomly sampling three variables: CCT, (true) IOP, and the material properties (E_{init}). It is recognized that there may be variables other than those examined herein that affect IOP measurement. For example, cornea radius of curvature is believed to have a minor effect on applanation reading as pointed out by Ehlers et al.^{2} It is assumed that the three variables chosen in this study are the ones with the most significant effect on tonometric IOP readings.

The statistical distribution of CCT is taken as 551.53±49 μm, corresponding to Bhan et al’s clinical measurements.^{19} Well-documented data are also available on the statistical variation of (measured) IOP. Pooled data from large epidemiological studies indicate that the mean population IOP is 16 mmHg, with a standard deviation (SD) of 3.0 mmHg.^{20} The data on IOP typically show a non-Gaussian distribution. However, a normal distribution is used here based on the assumption that it will satisfactorily represent the variation of IOP for the purposes of this study.

As described previously, E_{init} is the only independent material parameter in the Extended Fung Model formulated herein. Also, it is observed that E_{init} ranged from 0.12 to 0.35 MPa when calibrated using the experimental data of Bryant and McDonnell.^{5} This substantial variation is used to establish a mean and SD of E_{init} for the virtual population. The middle value (0.23 MPa) is selected for the mean of E_{init} and 0.04 MPa for the SD in a way that the range of calibrated E_{init} (0.12∼0.35 MPa) is covered by mean±3 SD.

The three variables are assumed to have normal distributions and to be uncorrelated; the mean CCT is 551.52±49 μm, IOP is 16±3.0 mmHg, and E_{init} is 0.23±0.04 MPa.

Bhan et al^{19} report the IOP measurements by Goldmann applanation tonometry from 181 normal corneas of 94 patients attending a general eye clinic over a 6-month period. The same size sample of 181 cornea finite element models is generated with the statistical distribution of the three variables (CCT, IOP, and E_{init}). The comparison of the results of finite element simulations with the clinical measurements, and the effect of each variable on the IOP measurements with Goldmann applanation tonometry, are reported. The finite element simulations of Goldmann applanation tonometry are performed using Abaqus software package 6.3.1 (Hibiit; Karlsson & Sorensen Inc, Pawtucket, Rhode Island).

### Results

#### IOP Statistics in the Virtual Population and Clinical Study

Figure 4 shows apparent IOP as a function of CCT for the virtual population of 181 cornea models, determined from the finite element simulation of Goldmann applanation tonometry. Figure 5 shows the digitized version of the IOP measurements with Goldmann applanation tonometry by Bhan et al^{19} in a clinical study of 181 normal corneas of 94 patients. The trend lines determined from regression analysis of the data are also shown in Figures 4 and 5, respectively.

Figures 4 and 5 show a reasonably similar distribution of IOP versus CCT for the virtual population and the clinical study. Considering the idealizations that have been made in the modeling of the corneal geometry and biomechanical properties and that only three variables (E_{init}, CCT, and true IOP) have been sampled to generate the virtual population, the similarity between Figures 4 and 5 is considered to be satisfactory. Numerical values of the statistics of the two samples are given in Table 1. The range of CCT values and the trend lines in IOP versus CCT are similar, although the clinical measurements show a somewhat higher degree of scatter. In the virtual cornea population, the scatter is due solely to unaccounted variation in true IOP and corneal material properties as represented by E_{init}, whereas in the clinical population, there may be additional variables^{21} involved.

Table 1: Statistics of Intraocular Pressure Versus Central Corneal Thickness Plots |

#### Effect of E_{init} and True IOP on Measured IOP

Variation of calculated IOP as a function of the other two variables, E_{init} and true IOP, in the virtual population is shown in Figures 6 and 7, respectively. Both figures show a considerable scatter in the data, similar to Figures 4 and 5.

Figure 6. Apparent Intraocular Pressure (IOP) in the Virtual Population of 181 Finite Element Cornea Models as a Function of Cornea In-Plane Young’s Modulus E |

Figure 7. Apparent Intraocular Pressure (IOP) in the Virtual Population of 181 Finite Element Cornea Models as a Function of True IOP (IOPT). |

Figure 6 demonstrates that the actual range of randomly sampled E_{init} values in the virtual population of 181 corneas is approximately 0.11 to 0.31 MPa. This range is roughly the mean ±2 SD, and essentially accounts for the scatter in the IOP versus CCT correlation for the virtual population (see Fig 4) because the correlation would be a single curve without any scatter when E_{init} is fixed as a constant. It also may be inferred from a comparison of Figure 4 with Figure 6 that the variation in corneal biomechanical properties (as represented in the virtual population by the parameter E_{init}) is as important as CCT in determining IOP by applanation tonometry.

To emphasize further the dependence of the scatter in the plot of IOP versus CCT on the variability in E_{init} and true IOP, two subsets of the virtual population are selected that include only the finite element cornea models whose E_{init} and true IOP fall within the *mean*±SD(σ) and *mean*±½SD(σ), respectively. The data for these two subsets are shown in Figures 8 and 9.

Figure 8. Data for the Subset of the Virtual Population with E |

Figure 9. Data for the Subset of the Virtual Population with E |

These figures show that the scatter in measured IOP diminishes as the range of E_{init} and true IOP narrows. It is also interesting to note that the slope of the trend lines is reduced as the range of the variation of the other two variables narrows. Numerical values of the slope of the trend line and the SD from the trend line are given in Table 2.

Table 2: Numerical Values of the Statistics for Different Subsets of the Virtual Population |

### Discussion

It has been demonstrated by means of realistic finite element simulations that the observed large scatter in standard correlations of clinical measurements of IOP versus CCT can be largely accounted for by reasonable individual variations in corneal biomechanical stiffness properties. Although previous studies by others have explored this claim via simple analytical models, the present study has used what is believed to be a more realistic geometrical and biomechanical model of the cornea as well as a more accurate simulation of Goldmann applanation tonometry. An important feature of the corneal biomechanical model used herein, which has been calibrated using experimental results, is that the cornea stiffens under increasing physiological stress (due to IOP). Thus, the true IOP becomes a significant variable influencing the error in measured IOP.

From the results of the finite element simulations on a virtual population of 181 corneas, it is found that the corneal biomechanical properties (as characterized in this study by the stiffness parameter E_{init}) are as important as the CCT in influencing measured (Goldmann) IOP. Plausible individual variations in corneal biomechanical properties may account for the observed large scatter in standard correlations of clinical measurement of IOP versus CCT. In addition, it is important to note that these findings are supported by the finite element simulations rather than experimental measurements in this study.

Unfortunately, no methods are currently available for assessing corneal biomechanical properties in vivo, and thus individual variations in these properties cannot be accounted for at present. It is evident that the determination of corneal biomechanical properties has important implications for improved design of individualized laser surgery procedures as well.

The authors have developed a new method for using the results of applanation tonometry for determining corneal biomechanical properties in vivo, as well as correcting measured IOPs. This method will be presented in a forthcoming publication.

### References

- Doughty MJ, Zaman ML. Human corneal thickness and its impact on intraocular pressure measures: a review and meta-analysis approach.
*Surv Ophthalmol*. 2000;44:367–408. doi:10.1016/S0039-6257(00)00110-7 [CrossRef] - Ehlers N, Bramsen T, Sperling S. Applanation tonometry and central corneal thickness.
*Acta Ophthalmol (Copenh)*. 1975;53:34–43. doi:10.1111/j.1755-3768.1975.tb01135.x [CrossRef] - Johnson M, Kass MA, Moses RA, Grodzki WJ. Increased corneal thickness simulating elevated intraocular pressure.
*Arch Ophthalmol*. 1978;96:664–665. - Whitacre MM, Stein RA, Hassanein K. The effect of corneal thickness on applanation tonometry.
*Am J Ophthalmol*. 1993;115:592–596. - Bryant MR, McDonnell PJ. Constitutive laws for biomechanical modeling of refractive surgery.
*J Biomech Eng*. 1996;118:473–481. doi:10.1115/1.2796033 [CrossRef] - Kwon TH.
*Minimally Invasive Characterization and Intraocular Pressure Measurement via Numerical Simulation of Human Cornea*[doctoral thesis]. Urbana, IL: University of Illinois at Urbana-Champaign; 2006. - Kwon TH, Ghaboussi J, Pecknold DA, Hashash YM. Effect of cornea material stiffness on measured intraocular pressure.
*J Biomech*. 2008;41:1707–1713. doi:10.1016/j.jbiomech.2008.03.004 [CrossRef] - Liou HL, Brennan NA. Anatomically accurate, finite model eye for optical modeling.
*J Opt Soc Am A Opt Image Sci Vis*. 1997;14:1684–1695. doi:10.1364/JOSAA.14.001684 [CrossRef] - Liu J, Roberts CJ. Influence of corneal biomechanical properties on intraocular pressure measurement: quantitative analysis.
*J Cataract Refract Surg*. 2005;31:146–155. doi:10.1016/j.jcrs.2004.09.031 [CrossRef] - Orssengo GJ, Pye DC. Determination of the true intraocular pressure and modulus of elasticity of the human cornea in vivo.
*Bull Math Biol*. 1999;61:551–572. doi:10.1006/bulm.1999.0102 [CrossRef] - Hjortdal JO. Regional elastic performance of the human cornea.
*J Biomech*. 1996;29:931–942. doi:10.1016/0021-9290(95)00152-2 [CrossRef] - Jue B, Maurice DM. The mechanical properties of the rabbit and human cornea.
*J Biomech*. 1986;19:847–853. doi:10.1016/0021-9290(86)90135-1 [CrossRef] - Woo SL, Kobayashi AS, Schlegel WA, Lawrence C. Nonlinear material properties of intact cornea and sclera.
*Exp Eye Res*. 1972;14:29–39. doi:10.1016/0014-4835(72)90139-X [CrossRef] - Nyquist GW. Rheology of the cornea: experimental techniques and results.
*Exp Eye Res*. 1968;7:183–188. doi:10.1016/S0014-4835(68)80064-8 [CrossRef] - Anderson K, El-Sheikh A, Newson T. Application of structural analysis to the mechanical behaviour of the cornea.
*J R Soc Interface*. 2004;1:3–15. doi:10.1098/rsif.2004.0002 [CrossRef] - Fung YC, Fronek K, Patitucci P. Pseudoelasticity of arteries and the choice of its mathematical expression.
*Am J Physiol*. 1979;237:H620–H631. - Schwartz NJ, Mackay RS, Sackman JL. A theoretical and experimental study of the mechanical behavior of the cornea with application to the measurement of intraocular pressure.
*Bull Math Biophys*. 1966;28:585–643. doi:10.1007/BF02476865 [CrossRef] - Velten K, Gunther M, Oberacher-Velten I, Lorenz B. Finite-element simulation of corneal applanation.
*J Cataract Refract Surg*. 2006;32:1073–1074. doi:10.1016/j.jcrs.2006.02.052 [CrossRef] - Bhan A, Browning AC, Shah S, Hamilton R, Dave D, Dua HS. Effect of corneal thickness on intraocular pressure measurements with the pneumotonometer, Goldmann applanation tonometer, and Tono-Pen.
*Invest Ophthalmol Vis Sci*. 2002;43:1389–1392. - Colton T, Ederer F. The distribution of intraocular pressures in the general population.
*Surv Ophthalmol*. 1980;25:123–129. doi:10.1016/0039-6257(80)90086-7 [CrossRef] - Whitacre MM, Stein R. Sources of error with use of Goldmann-type tonometers.
*Surv Ophthalmol*. 1993;38:1–30. doi:10.1016/0039-6257(93)90053-A [CrossRef]

Statistics of Intraocular Pressure Versus Central Corneal Thickness Plots

IOP | ||

Mean | 16.18 | 15.60 |

Standard deviation | 3.53 | 4.83 |

IOP increase per 10 | 0.25 | 0.28 |

Correlation coefficient (R) | 0.366 | 0.221 |

Numerical Values of the Statistics for Different Subsets of the Virtual Population

SD from the trend line | 3.27 | 1.75 | 0.92 |

IOP increase per 10 | 0.25 | 0.19 | 0.18 |