A Transient Diffusion Model of the Cornea for the Assessment of Oxygen Diffusivity and Consumption

Abstract

purpose. Oxygen diffusivity and consumption in the human cornea have not been directly measured yet; current models rely on properties measured in vitro in rabbit corneas. The aim of this study was to present a mathematical model of time-dependent oxygen diffusion that permits the estimation of corneal consumption and diffusivity.

methods. The current oxygen diffusion model was extended to include the temporal domain and was used to simulate in vivo noninvasive measurements of tear oxygen tension in human corneas.

results. The new model reproduced experimental data successfully, provided values for corneal diffusivity and consumption, and described the relationship between oxygen consumption and oxygen tension in the cornea. Estimated values were three times higher than those reported previously in in vitro rabbit experiments.

conclusions. This model allowed for the further investigation of oxygen transport in the cornea, including a better mathematical description and a determination of the transport properties of the cornea and the specific oxygen uptake rate of the tissue. The combination of this model and tear oxygen tension measurements can be useful in determining the individual oxygen uptake rate and exploring the relationship between oxygen transport and corneal abnormalities.

Extended duration of contact lens wear has been associated with a number of corneal disorders. One of the possible causes of corneal disorder from extended contact lens wear is that, under the resultant hypoxia, lactic acid production rates increase and the accumulation of lactate leads to stromal edema. Subsequent increases in osmotic pressure and tissue acidosis have also been related to altered endothelial morphology and function. ^{1 2} The cornea is an avascular tissue whose health depends on oxygen transport through its atmospheric boundary, but contact lenses provide additional resistance to this process. ³ Further, the oxygen transmissibility (Dk/t) of the lens is insufficient to evaluate the changes in oxygen supply to the cornea under contact lens wear because the relationship between these two variables is not linear. ⁴ Other terms derived from mathematical models, such as the equivalent oxygen pressure (EOP) and biological apparent oxygen transmissibility (BOAT), have been introduced to better assess the performance of contact lenses with regard to the physiological changes in oxygen supply they induce. ^{5 6} However, these measures rely heavily on the accuracy of the mathematical models and tissue properties used. 

Transport of oxygen in the cornea has been intensively studied since the beginning of mainstream use of contact lenses in the early 1970s. The current mathematical model of oxygen transport in the cornea was established by Fatt et al. at approximately the same time. ^{5 7 8 9 10} Since its introduction, this model has seen some improvements by different researchers ^{11 12} and has been extended to a two-dimensional model to determine oxygen transport through and around intracorneal lenses. ¹³ Its main limitation, though, is its use of tissue properties derived from in vitro experiments on rabbit corneas. This limitation arises because of the experimental difficulty in obtaining analogous values for in vivo human corneas. 

With the use of a phosphorescence-based technique and an optical sensor, Bonanno et al. ³ published some novel measurements of tear film oxygen tension during contact lens wear on human eyes in vivo. Their measurements at the cornea-contact lens interface provide abundant and useful data to characterize oxygen diffusion in the cornea. The concentration underneath the contact lens depends not only on lens transmissibility but on corneal tissue properties (oxygen diffusivity and metabolic consumption of oxygen). Their experiments give open-eye (OE) and closed-eye (CE) steady state oxygen tension at the lens-cornea interface and the transient response of oxygen tension at this location. To make best use of this experimental data, the current diffusion models must be extended, making the diffusion equation time dependent and not steady state. 

The aim of this article is to present a mathematical model of time-dependent oxygen diffusion in the cornea. The experimental data of Bonanno et al. ³ were used to validate the model and to determine the oxygen consumption and diffusivity of the cornea by fitting the model to three different cases of contact lens transmissibility. This model, in conjunction with tear film oxygen tension measurement, could then be used for multiple purposes, such as patient-specific estimation of corneal oxygen uptake, comparison of contact lens performance, or study of oxygen transport in the cornea and its relation to different corneal abnormalities. 

Methods

Solute transport by diffusion is governed by equation 1 , which is obtained combining Fick’s law and conservation of mass.  

\[D\ \frac{{\partial}^{2}c_{\mathrm{O}_{2}}}{{\partial}x^{2}}{-}Q{=}\ \frac{{\partial}c_{\mathrm{O}_{2}}}{{\partial}t}\]

It is a common practice to relate the concentration of oxygen to its partial pressure according to Henry’s law, \(c_{\mathrm{O}_{2}}\) = Po ₂× k, which leads to the established equation for oxygen transport in the cornea (equation 2) , whereDis the diffusivity (mm²/s),xis the position (mm), Po ₂is the partial pressure of oxygen (mm Hg),Qis the consumption of oxygen (nmol/mm³· s),kis the Henry’s solubility constant (nmol/mm³· mm Hg), andtis time (s).  

\[D\ \frac{{\partial}^{2}\mathrm{P}\mbox{\textsc{\mathrm{o}}}_{2}}{{\partial}x^{2}}\ {-}\ \frac{Q}{k}{=}\ \frac{{\partial}\mathrm{P}\mbox{\textsc{\mathrm{o}}}_{2}}{{\partial}t}\]

The consumption of oxygen is a result of corneal cell metabolism, which depends on a number of factors. In this study, the assumption was made that oxygen consumption depends only on the partial pressure of oxygen. An appropriate relationship should be continuous, yield zero consumption when Po ₂ is zero, and yield increasing consumption with increasing Po ₂until the saturation level is reached. A Michaelis-Menten–type relationship for oxygen uptake kinetics was used for this purpose (equation 3 ; see 1 2 3 4 Fig. 5 ), where Q* is the oxygen consumption rate at the saturated oxygen tension p* (in this case, 155 mm Hg) and a is a constant that determines the shape of the Q vs. Po ₂curve.  

\[Q{=}{\beta}{\times}Q{\ast};\ {\beta}{=}\ \frac{\mathrm{P}\mbox{\textsc{\mathrm{o}}}_{2}(a{+}p{\ast})}{p{\ast}(a{+}\mathrm{P}\mbox{\textsc{\mathrm{o}}}_{2})}\]

The cornea was divided into three layers, epithelium (0.05 mm), stroma (0.45 mm) and endothelium (0.005 mm). Each layer’s diffusivity and consumption interrelationship was constrained by the ratios measured by Fatt et al. ^{5 7 8} :  

\[Q{\ast}{=}Q_{\mathrm{stroma}}{=}0.11{\times}Q_{\mathrm{epithelium}}{=}0.02{\times}Q_{\mathrm{endothelium}}\]

 

\[D{=}D_{\mathrm{stroma}}{=}1.59{\times}D_{\mathrm{epithelium}}{=}5.66{\times}D_{\mathrm{endothelium}}\]

Three different contact lenses were used in the simulations: a high transmissibility hydrogel-silicone lens (Dk = 99 × 10⁻¹¹ cm² · mlO₂/s · mL · mm Hg; t = 0.1 mm), a thin hydrogel lens (Dk = 8.4 × 10⁻¹¹ cm² · mL O₂/s · mL · mm Hg; t = 0.06 mm), and a thick hydrogel lens (Dk = 8.4 × 10⁻¹¹ cm² · mL O₂/s · mL · mm Hg; t = 0.2 mm). 

The partial differential equation solver pdepe in MATLAB (The MathWorks Inc., Natick, MA) was used to solve the diffusion equation. Bonnano’s experiments were simulated with the mathematical model by changing the boundary³ condition at the anterior surface from CE (45 mm Hg) to OE (155 mm Hg) after 5 minutes of lens wear with the eye closed; during this time the CE steady state condition was reached. An example of the transient response of the whole simulation is shown in Figure 1 . 

An iterative least squares optimization process was developed to find the values of Q*/k, D, and a that fit best the data of Bonanno et al. This process minimized the difference between the tear film Po ₂of the transient model and the exponential fit proposed by Bonanno.³  

Results

The mathematical model reproduced the experiments successfully, and the optimization process provided realistic values for the optimized variables for the corneal stroma: Q*/k = 1.87 mm Hg/s, D = 2.81 × 10⁻³ mm²/s, and a = 20. The calculated steady state Po ₂at the tear film was within the SD of the experiments by Bonanno et al. ³ (Fig. 2) , and the model followed the experimental data closely for each contact lens. 

Values for Q* and Dk in converted units are reported in Table 1for direct comparison with the experimental results from Fatt et al. ¹⁰ Oxygen solubility in water at 35°C of k = 1.0268 · 10⁻² μmol/mL · kPa was used for this purpose. ¹⁴  

The level of agreement between the model and the data of Bonanno was not the same for every lens. Except for the CE steady state oxygen tension under the thick hydrogel, the predicted values lie within the SD of the experimental data. For the high-transmissibility lens, OE and CE steady state values were close to the average experimental data (Fig. 2) . The model overestimated the CE steady state for the thin hydrogel by 1.8 mm Hg but underestimated the CE steady state for the thick hydrogel by 2.0 mm Hg. In all cases, the transient response was accurately reproduced, as evidenced by the delayed increase of tear Po ₂ (approximately 7 seconds) with the thick hydrogel after eye opening (Fig. 3C)and the time delay needed to reach steady state (Fig. 3) . 

Additional postlens tear oxygen tension measurements have been published recently by Bonanno et al. (Bonanno JA, et al. IOVS 2008;49:ARVO E-Abstract 4846). The transient model presented here agrees with the new OE and CE steady state oxygen tension values over several lens transmissibilities. Unfortunately, the overall transient behavior could not be compared because these data were not explicitly presented 

Discussion

Oxygen diffusivity and consumption have yet to be directly measured in the human cornea. The aim of this study was to present a mathematical model of time-dependent oxygen diffusion in the cornea. Such a model in conjunction with in vivo experimental data permits the estimation of human oxygen diffusivity and consumption. The model successfully reproduces experimental results for transient oxygen tension after CE contact lens wear and steady state oxygen tension, which implies a high fidelity in the model and corneal tissue properties. 

Previous mathematical models used to infer knowledge about oxygen transport in the human cornea have relied on tissue properties measured in vitro with a diffusion chamber and rabbit corneas. ^{10 11 12 13} The noninvasive in vivo experimental measures by Bonanno ³ provide sufficient data to determine oxygen consumption and diffusivity of the human cornea. Although the relative merits of diffusion and tear oxygen experiments are debatable, the preference to acquire human in vivo data rather than rabbit in vitro data would seem under dispute. If the tissue properties from previous models are used, the oxygen tension is overestimated at the tears for CE and OE conditions. The initial time response is also faster, yet more time is needed to reach steady state (Fig. 4A) . Nevertheless, only small differences can be observed in the steady state oxygen distribution through the cornea (Fig. 4B) . The threefold increase in oxygen consumption has, therefore, little influence on the resultant oxygen distribution except at the epithelium and the anterior stroma. The calculated diffusivity is also approximately three times higher; hence, more oxygen is consumed, but it is also transported more quickly, leading to similar oxygen tension distribution. The oxygen supply from the aqueous humor explains why the difference is smaller at the endothelium and the posterior stroma. 

The computational model presented here relies on experimental data and boundary conditions chosen for the anterior and posterior surfaces. Boundary conditions and corneal geometry were set using mean values from the experimental population, but these parameters will vary in individual cases, leading to results with implications that may differ from those of the main population. For this reason, we believe that patient-specific modeling would be the ideal way to predict abnormalities related to altered oxygen transport in the cornea. Although there are some inconsistencies in the available experimental results (Bonanno JA, et al. IOVS 2008;49:ARVO E-Abstract 4846), ³ transient oxygen tension measures at the tears have the potential to provide accurate data for these models. 

The boundary cell layers of the cornea (epithelium and endothelium) have lower diffusivity than the stroma because of their tight structure and protective function, resulting in a decrease in oxygen tension at both boundaries (Fig. 4B) . This decrease is enhanced because oxygen consumption by these layers is increased as a result of higher cell density and cell activity. Therefore, for a realistic idea of the distribution of oxygen tension in the cornea, the different layer properties must be taken into account. Freeman ⁸ and Freeman and Fatt ⁹ previously measured an important difference in oxygen diffusivity and consumption between the layers of the rabbit cornea. In the present study, the assumption was made that the consumption/diffusivity ratio between the layers was the same as what had been measured in rabbits (equations 4and 5 ). This assumption was made for two reasons. First, calculating the “average” diffusivity and consumption for the whole cornea would discard the layer information, and the model would then be inadequate to predict the oxygen tension distribution in the cornea. Second, it is not possible to calculate the diffusivity and consumption of each layer based on tear oxygen tension. Intracorneal measurements would be necessary for this purpose. 

The oxygen consumption of corneal cells depends on the nutrient availability at every location. Different relations have been proposed in the literature to define corneal oxygen consumption as a function of oxygen tension. ^{10 11 12 15} An appropriate relationship should give zero consumption at zero Po ₂and should increase Q as Po ₂increases. In an initial optimization of the model, the consumption law was not taken into account. Only Q/k and D were optimized, and the Michaelis-Menten constant, a, was set constant. This model reproduced individual experiments well (for a given lens), but it could not reach a good solution for the three lenses simultaneously. When the consumption law was included in the optimization loop, the model converged to a successful solution. By this iteration, we found that the oxygen consumption law was of key importance in reaching the solution and was a fundamental part of the oxygen transport model in the cornea. The resultant curve (Fig. 5)was close to the logarithmic relationship between oxygen flux to the cornea and postlens Po ₂proposed by Brennan. ¹⁶ A similar relationship was found by Bibby et al. ¹⁷ when they measured intervertebral disc (IVD) cell metabolism. The IVD is the largest avascular tissue in the human body; its nutrient transport and oxygen demands are similar to those of the cornea. 

As stated earlier, corneal oxygen consumption depends on oxygen availability, but it may also depend on other factors such as the availability of additional nutrients (glucose) or the tissue pH. Harvitt and Bonanno ¹⁸ found that corneal oxygen consumption increases under acute pH reductions in rabbit corneas in vitro. Although other researchers have reported decreased corneal oxygen uptake after prolonged use of contact lenses, ^{19 20} this may correspond to a long-term adaptive change in cellular metabolism induced by an altered oxygen distribution. To develop the existing oxygen transport models, it is important to improve our understanding on how corneal cell metabolism behaves under the extreme environment found in the cornea during contact lens use. 

A standard oxygen diffusion model has been improved with the addition of the temporal dimension; this has allowed us to estimate corneal diffusivity and oxygen consumption directly from noninvasive in vivo measures in humans. We found that values of oxygen consumption and diffusivity are approximately three times those that have been proposed. The oxygen consumption law has also been characterized in detail. We think the combination of this model with tear film oxygen tension measurements can be useful in determining individual oxygen uptake rates and in further investigating the relationship between oxygen transport and corneal abnormalities. 

 Supported by NCCR Co-Me of the Swiss National Science Foundation and the Swiss Innovation Promotion Agency, CTI (6801.3;6 FHS-LS).

 Submitted for publication June 23, 2008; revised August 6, 2008; accepted January 22, 2009.

 Disclosure: X. Larrea, None; P. Büchler, None

 The publication costs of this article were defrayed in part by page charge payment. This article must therefore be marked “advertisement” in accordance with 18 U.S.C. §1734 solely to indicate this fact.

 

Figure 5.

 

Relation between oxygen consumption rate and oxygen tension after optimization. This curve is the plot of equation 3 , with a = 20 and p* = 155.

Figure 1.

 

Simulation of the experiments of Bonanno et al. ³ using the mathematical model. Subjects were asked to close their eyes for 5 minutes. When they opened their eyes, oxygen tension was measured at the tear film (between the lens and the cornea) for at least 40 seconds (○) with the phosphorescence decay of a porphyrin-protein complex (oxygen sensitive dye) bound to the posterior surface of the contact lens.

Figure 2.

 

Steady state tear oxygen tension for open-eye and closed-eye conditions computed by the model and compared with the mean oxygen tension measured by Bonanno et al. ³ for each lens. Error bars, SD.

Table 1.

 

Oxygen Consumption and Permeability of Human Corneal Stroma

	Transient Model (human)	Fatt et al. 10 (rabbit)
Q*stroma	5.75 · 10−5 \batchmode \documentclass[fleqn,10pt,legalpaper]{article} \usepackage{amssymb} \usepackage{amsfonts} \usepackage{amsmath} \pagestyle{empty} \begin{document} \(\frac{ml\ O_{2}}{ml\ {\cdot}\ s}\) \end{document}	2.24 · 10−5 \batchmode \documentclass[fleqn,10pt,legalpaper]{article} \usepackage{amssymb} \usepackage{amsfonts} \usepackage{amsmath} \pagestyle{empty} \begin{document} \(\frac{ml\ O_{2}}{ml\ {\cdot}\ s}\) \end{document}
Dk stroma	8.62 · 10−10 \batchmode \documentclass[fleqn,10pt,legalpaper]{article} \usepackage{amssymb} \usepackage{amsfonts} \usepackage{amsmath} \pagestyle{empty} \begin{document} \(\frac{cm^{2}\ {\cdot}\ ml\ O_{2}}{s\ {\cdot}\ ml\ {\cdot}\ mmHg}\) \end{document}	3.0 · 10−10 \batchmode \documentclass[fleqn,10pt,legalpaper]{article} \usepackage{amssymb} \usepackage{amsfonts} \usepackage{amsmath} \pagestyle{empty} \begin{document} \(\frac{cm^{2}\ {\cdot}\ ml\ O_{2}}{s\ {\cdot}\ ml\ {\cdot}\ mmHg}\) \end{document}

 Values shown here were derived from the mathematical model and compared with experimental values from rabbit corneas reported by Fatt et al. ¹⁰ Values for epithelium and endothelium can be easily calculated from equations 4and 5 .

Figure 3.

 

Representative results for the tear oxygen tension after 5 minutes of CE lens wear (A–C) and the steady state oxygen tension distribution through corneal thickness (D–F), for contact lens transmissibility of 99 U (A, D), 14 U (B, E), and 4.2 U (C, F): ×10⁻¹¹ cm² · mL O₂/s · mL · mm Hg. (gray) Lens regions.

Figure 4.

 

Comparison of oxygen tension at the tears (A) and oxygen tension distribution through the cornea (B) using tissue properties from Table 1(Dk/t = 99 × 10⁻¹¹ cm² · mL O₂/s · mL · mm Hg). Simulations with other lens transmissibilities led to similar results.