# Addressing Parameter Variability in Corneal Biomechanical Models: A Stepwise Approach for Parameters' Optimization

José González-Cabrero, Carmelo Gómez, Manuel Paredes, Francisco Cavas

► To cite this version:

José González-Cabrero, Carmelo Gómez, Manuel Paredes, Francisco Cavas. Addressing Parameter Variability in Corneal Biomechanical Models: A Stepwise Approach for Parameters' Optimization. Biomimetics, 2025, 10 (10), pp.683. 10.3390/biomimetics10100683. hal-05489458

HAL Id: hal-05489458

https://ut3-toulouseinp.hal.science/hal-05489458v1

Submitted on 2 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Article

# Addressing Parameter Variability in Corneal Biomechanical Models: A Stepwise Approach for Parameters' Optimization

José González-Cabrero <fnref n="1" />,<fnref n="2" />,<fnref n="3" /> ID, Carmelo Gómez <fnref n="1" />,<fnref n="2" /> ID, Manuel Paredes <fnref n="3" /> ID and Francisco Cavas <fnref n="1" />,<fnref n="2" />,* ID

## Abstract

Biomechanical modeling of the cornea is crucial for understanding the progression of some ocular diseases and optimizing surgical treatments. However, hyperelastic nonlinear material models, such as those used for corneal tissue, often yield highly variable parameter sets in the scientific literature, influenced by factors like the chosen optimization intervals and differences between tensile and inflation test curve optimization, both of which are addressed in this study. This variability complicates the understanding of corneal mechanical properties. In this research, the aim is to optimize and calibrate the key parameters of the corneal material model, particularly focusing on c₁, c₂, k₁ and k₂, using the Holzapfel-Gasser-Ogden (HGO) hyperelastic model, and a novel methodology is proposed that separately estimates the isotropic and anisotropic components in a stepwise manner, addressing the issue of multiple parameter sets fitting experimental curves similarly. This approach helps to standardize corneal material models and improve the reliability of parameter estimations. Moreover, accurate biomechanical characterization within this framework contributes not only to clinical applications but also to biomimetics, inspiring the design of artificial corneal substitutes and bioengineered materials.

**Keywords:** corneal biomechanics; Holzapfel model; inflation test; hyperelastic material; parameters estimation

## 1. Introduction

The human eye receives over 90% of external sensory information, highlighting the critical importance of ocular health for quality of life [1]. The cornea, contributing approximately two-thirds of the eye's refractive power [2], plays a vital role in this process. This avascular structure, with an average thickness of 550 µm centrally and 750 µm peripherally, is primarily composed of stroma, which accounts for 90% of its thickness and consists of collagen fibers embedded in a proteoglycan matrix and arranged in a region-specific orientation [3].

The cornea's mechanical behavior is commonly modeled using anisotropic hyperelastic models [4,5], essential to optimize surgical outcomes and understanding diseases like ectasia, characterized by structural weakening of collagen fibers [6]. Hyperelastic models, such as the Gasser-Holzapfel-Ogden (GHO) and Holzapfel-Gasser-Ogden (HGO) models, have become standard in corneal studies due to their capacity to describe anisotropic

<figure><img src="image_2.png" /></figure>

check for
updates

Academic Editor: Yu Pan

Received: 3 September 2025

Revised: 1 October 2025

Accepted: 8 October 2025

Published: 10 October 2025

**Citation:** González-Cabrero, J.; Gómez, C.; Paredes, M.; Cavas, F. Addressing Parameter Variability in Corneal Biomechanical Models: A Stepwise Approach for Parameters' Optimization. *Biomimetics* **2025**, *10*, 683. https://doi.org/10.3390/biomimetics10100683

**Copyright:** © 2025 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

<footnote n="1">¹ Department of Structures, Construction and Graphical Expression, Technical University of Cartagena, 30202 Cartagena, Spain; jose.gonzalezc@edu.upct.es (J.G.-C.); carmelo.gomez@upct.es (C.G.)</footnote>

<footnote n="2">² Bioengineering and Applied Computational Simulation Research Group, Technical University of Cartagena, 30202 Cartagena, Spain</footnote>

<footnote n="3">³ ICA, Université de Toulouse, UPS, INSA, ISAE-SUPAERO, MINES-ALBI, CNRS, 3 rue Caroline Aigle, 31400 Toulouse, France; paredes@insa-toulouse.fr</footnote>

<footnote>* Correspondence: francisco.cavas@upct.es; Tel.: +34-968338856</footnote>
properties and fiber dispersion [7–10]. These models assume near-incompressibility and provide parameters representing isotropic stiffness, fiber stiffness, and non-linear deformation. However, parameter variability across studies, influenced by differences in experimental techniques, tissue properties, and mathematical simplifications, complicates direct comparisons and the reliability of conclusions between studies [11–14].

Tensile and inflation tests are the primary experimental methods to estimate corneal parameters. However, the two approaches yield inconsistent results due to differences in testing conditions and the mathematical models applied. For example, studies using tensile tests, such as those by [11–13], have reported wide variability in parameter values, while inflation tests, as demonstrated by [14,15], highlight the sensitivity of material parameters to experimental conditions.

Experimental comparisons have shown that strip extensometry and inflation produce markedly different responses [16] and that multiple loading scenarios are critical for accurate parameter identification [17]. Moreover, recent advances such as Brillouin microscopy have enabled non-invasive biomechanical characterization [18], while age- and IOP-related collagen recruitment studies further emphasize the variability in corneal mechanical response [19]. Recent ex vivo and animal studies also demonstrate discrepancies between uniaxial and biaxial testing outcomes [20], highlighting the difficulty of establishing universally valid parameter sets.

Variations in geometric measurements, non-uniformity of corneal properties, and model-specific assumptions contribute to these discrepancies. Furthermore, high correlations between experimental data and model predictions (high correlation coefficients $R^2$) do not always ensure unique or physically valid parameter sets. This issue arises when parameter compensations occur, leading to multiple sets that fit experimental data similarly. Consequently, it becomes necessary to establish robust methods for parameter estimation that minimize variability and improve reliability.

This study aims to address these challenges by (1) developing a method to characterize corneal material parameters of the HGO model through a novel stepwise optimization approach that separates the estimation of the isotropic ground-matrix components from the anisotropic fiber parameters to reduce variability and ensure physical consistency in the estimated parameters; (2) proposing parameter bounds and optimization algorithms to demonstrate the high variability inherent in the characterization problem, driven by the large number of parameter combinations that yield similarly high correlation with the experimental curve.

Beyond its clinical relevance, this research also fits within a biomimetic perspective. The cornea, combining transparency, anisotropy, and mechanical resilience, serves as an inspiring model for bioengineered materials. Advances in contact lenses, tissue substitutes, and artificial corneal constructs seek to replicate these properties, highlighting how biomechanical modeling can guide the design of biomimetic solutions that emulate their structural and functional efficiency. In this context, an accurate characterization of the cornea's hyperelastic properties is essential to ensure that such designs achieve physiological-like behavior.

## 2. Materials and Methodology

In this research, a detailed methodology is implemented to optimize the material parameters of corneal tissue ($c_1, c_2, k_1$ and $k_2$) using tensile and inflation tests within the Holzapfel–Gasser–Ogden (HGO) model framework [7]. This model is adopted because it combines physiological realism with simplicity: it separates the isotropic ground matrix ($c_1, c_2$) from collagen–fiber reinforcement and reproduces the cornea’s characteristic exponential stiffening ($k_1, k_2$). Other models, such as the five-parameter Mooney–Rivlin
formulation [21], can also capture pronounced nonlinearity; however, they make it more
difficult to distinguish which parameters govern the initial, matrix-dominated segment of
the curve and which govern the later, fiber-driven non-linear response.

Both tensile and inflation tests are mathematically defined to simulate the cornea's
mechanical behavior under physiological conditions. The methodology is summarized in
the scheme shown in Figure 1.

<figure><img src="image_3.png" /><figcaption>**Figure 1.** The four steps that comprise the material parameter estimation analysis and the contributions of the study.</figcaption></figure>
Python version 3.11.9 (Python Software Foundation) is used for data analysis, pa- rameter optimization, and curve fitting. Core libraries were NumPy for array algebra, SciPy (differential_evolution, best/1/bin) for parameter optimization, and Matplotlib for plotting.

In the optimization process, a differential evolution algorithm is employed to iteratively adjust the material parameters ($c_1, c_2, k_1$ and $k_2$) to minimize the difference between the predicted stress values ($\sigma_1, \sigma_2, \sigma_3$) and the experimental data.

This approach, known as least-squares minimization, aims to reduce the total error in the predicted stress values. Penalties are also introduced to avoid physically meaningless solutions, such as negative stress values ($\sigma_1 \le 0$), ensuring that the parameters remain physically valid throughout the optimization process. The coefficient of determination, $R^2$ is used to evaluate the quality of the fit between the model's predictions and the experimental data; it measures how much of the variance in the experimental data is captured by the model, with a value close to 1 indicating a strong fit.

The optimization analysis is conducted by implementing three different parameter bounds to demonstrate that various parameter sets can achieve high experimental fits with $R^2$ values close to 1, while still maintaining physical validity in order to explain the variability of parameter sets commonly reported in the literature. The goal is to analyze which parameter interval is preferable, ranging from highly restrictive, moderately restrictive, to non-restrictive bounds. Only different intervals are set for $c_1$ and $c_2$, while maintaining for $k_1$, and $k_2$. The optimization process is iterated within these intervals until the set that minimizes the error is found.

Experimental data for the tensile and inflation tests' curves are obtained from published investigations. To mitigate inconsistencies and enhance parameter stability, a stepwise optimization approach is proposed in Section 3.3. First, $c_1$ is optimized for small deformations dominated by isotropic behavior. Once $c_1$ is established, $k_1$ and $k_2$ are optimized to capture the anisotropic non-linear behavior at larger deformations, while $c_2$ is neglected. This sequential approach simplifies the optimization process and reduces parameter variability, addressing a major limitation in prior studies. By first fitting $c_1$ on the small-strain, matrix-dominated segment and then fitting $k_1$ and $k_2$ on the non-linear, fiber-dominated segment, the estimation is split into lower-dimensional subproblems, which removes many compensating parameter combinations and thereby reduces parameter variability even when $R^2$ remains high.

## 2.1. Biomechanical Modeling of Corneal Material

In hyperelastic material modeling, the local deformation is described by the deformation gradient $F = \frac{\partial x}{\partial X}$, which maps a material point from its reference position $X$ to its current deformed position $x$. Its determinant $J = \det(F)$ is the local volume ratio; for incompressible response, $J = 1$ (volume preservation). The Cauchy stress tensor is written as

$$ \sigma = -p\mathbf{I} + \frac{2}{J}F\frac{\partial\psi}{\partial C}F^T, \qquad (1) $$

where $\sigma$ is the Cauchy stress tensor, $p$ is a Lagrange multiplier that enforces incompressibility and represents the hydrostatic/volumetric part of the stress, and $I$ is the identity tensor (unit tensor) that produces an equal normal contribution in every direction. The symbol $\psi$ denotes the strain-energy density function. The tensor $C = F^T F$ is the right Cauchy-Green deformation tensor, which measures strain in the reference configuration. The derivative $\frac{\partial \psi}{\partial C}$ is the tensor of partial derivatives of the energy with respect to $C$. In practice, $\psi$ is often split into a volumetric term and a deviatoric (shape-changing) term. In this work, incompressibility is enforced by the Lagrange multiplier $p$, so the volumetric contribution
is represented by the hydrostatic term $-pI$ in Equation (1); therefore, $\psi$ in Equation (1) denotes only the deviatoric strain-energy density. The deviatoric component of $\psi$ in the HGO model is as follows:

$$ \psi_{\text{deviatoric}} = c_1 (\bar{I}_1 - 3) + c_2 (\bar{I}_2 - 3) + \frac{k_1}{2k_2} \sum_{i=4,6} (\exp[k_2(\bar{I}_i - 3)^2] - 1), \quad (2) $$

where $c_1 + c_2$ defines shear modulus, $k_1$ describes fiber stiffness, and $k_2$ accounts for non-linearity in large strains. The invariants are $I_1 = \text{tr}(C)$, $I_2 = \frac{1}{2}((\text{tr}(C))^2 - \text{tr}(C^2))$, $I_4 = \mathbf{m}_0\mathbf{C}\mathbf{m}_0$ and $I_6 = \mathbf{n}_0\mathbf{C}\mathbf{n}_0$, where $\mathbf{m}_0 = [1, 0, 0]$ and $\mathbf{n}_0 = [0, 1, 0]$ are initial fiber directions.

Assuming incompressibility and following the theory of continuum mechanics, the Cauchy stress tensor can be expressed in terms of the invariants and the left Cauchy-Green deformation tensor $\mathbf{b} = \mathbf{F}\mathbf{F}^T$ through the following expression:

$$ \sigma = -pI + \frac{2}{J} \left[ \left( \frac{\partial \psi}{\partial I_1} + I_1 \frac{\partial \psi}{\partial I_2} \right) \mathbf{b} - \frac{\partial \psi}{\partial I_2} \mathbf{b}^2 + \frac{\partial \psi}{\partial I_4} \mathbf{m} \otimes \mathbf{m} + \frac{\partial \psi}{\partial I_6} \mathbf{n} \otimes \mathbf{n} \right], \quad (3) $$

where $\mathbf{m} = \mathbf{F}\mathbf{m}_0$ and $\mathbf{n} = \mathbf{F}\mathbf{n}_0$ are the orientation vectors in the preferred directions in the deformed configuration, and $I$ is the identity matrix.

## Tensile Test

During the tensile test, the tissue is clamped at both ends, and a tensile force is applied along the length in the 1-direction, resulting in elongation characterized by the stretch ratio $\lambda_1 > 1$. To maintain volume due to incompressibility, the tissue contracts in directions 2 and 3 perpendicularly. The deformation gradient tensor is given by

$$ \mathbf{F} = \begin{bmatrix} \lambda_1 & 0 & 0 \\ 0 & \lambda_2 & 0 \\ 0 & 0 & \lambda_3 \end{bmatrix}, \quad (4) $$

where $J = 1 = \lambda_1\lambda_2\lambda_3$ ensures incompressibility. For a tensile load in the 1-direction, with the tissue confined to the 1–2 plane, $\lambda_3 = (\lambda_1\lambda_2)^{-1}$ and $\lambda_2 = \frac{1}{\sqrt{\lambda_1}}$. The Cauchy stress tensor is expressed as follows:

$$ \sigma = \begin{bmatrix} \sigma_1 & 0 & 0 \\ 0 & \sigma_2 & 0 \\ 0 & 0 & \sigma_3 \end{bmatrix} \quad (5) $$

For uniaxial tensile tests, stress is applied only in the 1-direction, so $\sigma_1 = \sigma$; $\sigma_2 = \sigma_3 = 0$. Figure 2 illustrates the virtual tensile setup where the specimen is clamped at both ends and loaded along direction 1.

Tensile strips are typically cut along the inferior-superior or nasal-temporal anatomic axes; these are widely considered the principal corneal fiber directions in the central stroma of the cornea; to generalize across published datasets, this orientation is assumed for all cases. Accordingly, two preferred in-plane fiber families are defined and taken orthogonal in the reference configuration: $\mathbf{m}_0 = [1, 0, 0]$ aligned with the loading direction and $\mathbf{n}_0 = [0, 1, 0]$ aligned transversely. The local axis {1, 2, 3} is also indicated in the figure.

The schematic Figure 3 represents the stress-stretch curves used in this research for three different corneal specimens under tensile testing obtained from [13,14], highlighting varying degrees of non-linearity in their mechanical responses. Additionally, it illustrates the collagen fiber alignment during the stretching process in three stages.
<figure><img src="image_4.png" /><figcaption>Figure 2. Schematic image of a tensile test.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 3. Schematic image of three tensile test curves.</figcaption></figure>

## 2.2. Inflation Test

The inflation test is critical to understand the cornea's biomechanics, as it closely replicates in vivo conditions by applying pressure similar to intraocular pressure (IOP). Unlike uniaxial tensile tests, it provides insights into the realistic mechanical response of corneal tissue under physiological conditions. Note that an axisymmetric thin-membrane assumption is adopted to match the experimental inflation setup and as a test simplification, implying equibiaxial loading ($\sigma_1 = \sigma_2 = \sigma$; $\sigma_3 = 0$). The cornea, however, is not a perfect spherical cap and exhibits regional variation and fiber dispersion; these aspects are acknowledged in Section Limitations of the Study. The deformation gradient tensor for the inflation test is

$$F = \begin{bmatrix} \lambda_1 & 0 & 0 \\ 0 & \lambda_2 & 0 \\ 0 & 0 & \lambda_3 \end{bmatrix}, \qquad (6)$$

where $\lambda_1$, $\lambda_2$, and $\lambda_3$ are the stretches in the meridional $e_\phi$, circumferential $e_\theta$, and radial $e_R$ directions, respectively. Due to incompressibility, $\lambda_1 = \lambda_2$ and $\lambda_3 = 1/\lambda_1^2$, with $J = \lambda_1^2\lambda_3 = 1$.

Assuming a thin-membrane model, the cornea is treated as sufficiently thin to neglect through-thickness stresses and bending effects in the radial direction $\sigma_R$, resulting in a
biaxial stress state in the meridional $\sigma_\phi$ and circumferential $\sigma_\theta$ directions. The stress tensor is as follows:

$$ \sigma = \begin{bmatrix} \sigma & 0 & 0 \\ 0 & \sigma & 0 \\ 0 & 0 & 0 \end{bmatrix}, \qquad (7) $$

where $\sigma$ represents the stress in both directions. This setup allows studying corneal behavior under applied pressure, as shown in Figure 4.

<figure><img src="image_6.png" /><figcaption>Figure 4. Schematic figure of an inflation test. (a) Applied pressure. (b) Meridional $e_\phi$, circumferential $e_\theta$ and radial directions $e_R$. Fiber directions towards meridional $\mathbf{n}_0$ and circumferential $\mathbf{m}_0$.</figcaption></figure>

The equations to convert pressure-apical displacement data into stress-stretch relationships are detailed in Appendix A according to [22]. These include the conversion of pressure units, calculation of the radius of curvature and stretch, and the application of Laplace's law to determine corneal stress. The equations account for geometric changes in the cornea during deformation, assuming it behaves as a spherical shell. Figure 5 shows the transformation of the curves [10] before and after applying the normalization process.

<figure><img src="image_7.png" /><figcaption>Figure 5. Three corneal curves of the inflation test. Conversion from pressure-apical displacement (a) into stress-stretch (b).</figcaption></figure>

## 2.3. Parameter Optimization for Tensile and Inflation Tests

The four parameters ($c_1, c_2, k_1$ and $k_2$) are fitted simultaneously using the differential evolution algorithm (strategy best/1/bin), implemented in Python via SciPy (from scipy.optimization import differential_evolution). This algorithm efficiently handles complex, non-linear problems typical of hyperelastic material modeling. The optimizer minimizes a least-squares error between model predicted and experimental stresses, augmented with large penalties to enforce test physics: penalties avoid negative stress, and non-monotonic axial stress as stretch increases.
The algorithm runs up to 10,000 iterations with a tolerance of 1 × 10⁻⁶ on the objective. Parameter bounds are implemented as box constraints; their restrictiveness is varied to demonstrate that different intervals can yield different parameter sets, with further discussion provided in Section 3. Goodness of fit is reported with the coefficient of determination R² = 1 - SS<sub>Res</sub>/SS<sub>Tot</sub> (where SS<sub>Res</sub> is the residual sum of squares (predicted vs. experimental), and SS<sub>Tot</sub> is the total variance in experimental data).

Section 3.3 introduces a novel stepwise approach seeking to stabilize parameter identification by separating the estimation of the isotropic part ($c_1$) from the anisotropic part ($k_1$ and $k_2$).

# 3. Results

## 3.1. *Parameter Optimization: Tensile Test*

In the optimization process, three bound intervals are tested (Interval 1 (slightly restrictive), Interval 2 (non-restrictive), and Interval 3 (highly restrictive)) to demonstrate how bound choices affect the resulting parameter sets.

The rationale for the bounds selection is explanatory (to show how bounds affect identifiability), keeping the fiber terms effectively unrestricted ($k_1$ in [0, 100] MPa and $k_2$ in [0, 1000]) while varying only the isotropic bounds. Interval 1 ("slightly restrictive") uses positive-only ranges for $c_1$ and $c_2$ (e.g., [0, 50] MPa); whether the upper limit is 50 or 100 MPa is immaterial, provided the range is sufficiently large to calculate without restriction. Interval 2 ("non-restrictive") allows $c_1$ and $c_2$ in [-100, 100] MPa, admitting negative $c_2$ values, sometimes reported to compensate $c_1$; any solutions with $c_1 + c_2 \le 0$ are deemed physically inadmissible even if the statistical fit is high. Interval 3 ("highly restrictive") adopts narrow literature-based ranges around typical $c_1$ and $c_2$ values to illustrate how tight priors can over-constrain the fit and promote parameter compensation in $k_1$ and $k_2$.

Interval 1 (Slightly Restrictive) provides the most reliable results, with parameter bounds of [0, 50] for $c_1$ and $c_2$, and [0, 100] and [0, 1000] for $k_1$ and $k_2$, respectively. This interval ensures physical validity, avoiding negative shear stress values ($c_1 + c_2$). The algorithm consistently assigns $c_2$ near zero, simplifying the isotropic component to a single parameter ($c_1$), which effectively captures the material's initial stiffness. This interval results in high R² values while maintaining interpretability, with parameter values reflecting the stiffness and fiber behavior across different non-linearities of the corneal tissue.

Interval 2 (Non-Restrictive), with bounds of [-100, 100] for $c_1$ and $c_2$, leads to unrealistic results, such as negative shear stress, rendering the parameters physically invalid despite high R² values. This demonstrates that overly broad bounds can cause overfitting, producing results that are statistically accurate but lack meaningful representation of material properties.

Interval 3 (Highly Restrictive) tightly constrains $c_1$ and $c_2$ to [0.005, 0.05] and [-0.05, -0.005], respectively, while maintaining the same bounds for $k_1$ and $k_2$ as in Interval 1. Although this interval produces physically valid results, it forces some parameters to their boundary values, leading to compensatory adjustments in $k_1$ and $k_2$. While useful in specific scenarios where parameter ranges are well-defined, such tight constraints may limit the optimization process and risk over-constraining the system.

Results of the estimated HGO material parameters for the optimization intervals are shown in Table 1.
<table><caption>**Table 1.** Estimated HGO material parameters from the tensile tests for three different optimization intervals of the three curves.</caption><thead><tr><th>Interval</th><th>Curve</th><th>c<sub>1</sub><br>(MPa)</th><th>c<sub>2</sub><br>(MPa)</th><th>k<sub>1</sub><br>(MPa)</th><th>k<sub>2</sub><br>(-)</th><th>R<sup>2</sup></th></tr></thead><tbody><tr><td rowspan="3">Interval 1 (Slightly Restrictive):<br>[0, 50], [0, 50], [0, 100], [0, 1000]</td><td>Curve 1 (Highly non-linear)</td><td>0.39</td><td>0</td><td>0.2721</td><td>24.656</td><td>0.99</td></tr><tr><td>Curve 2 (Moderately non-linear)</td><td>1.57</td><td>0</td><td>0.8743</td><td>8.254</td><td>0.99</td></tr><tr><td>Curve 3 (Slightly non-linear)</td><td>3.74</td><td>0</td><td>0.6872</td><td>9.378</td><td>0.97</td></tr><tr><td rowspan="3">Interval 2 (Non-Restrictive):<br>[-100, 100], [-100, 100], [0, 100],<br>[0, 1000]</td><td>Curve 1 (Highly non-linear)</td><td>12.57</td><td>-12.89</td><td>0.1911</td><td>26.914</td><td>0.99</td></tr><tr><td>Curve 2 (Moderately non-linear)</td><td>39.89</td><td>-40.43</td><td>0.3106</td><td>10.726</td><td>0.99</td></tr><tr><td>Curve 3 (Slightly non-linear)</td><td>65.97</td><td>-66.05</td><td>0.0426</td><td>0.750</td><td>0.99</td></tr><tr><td rowspan="3">Interval 3 (Highly Restrictive):<br>[0.005, 0.05], [-0.05, -0.005], [0,<br>100], [0, 1000]</td><td>Curve 1 (Highly non-linear)</td><td>0.05</td><td>-0.005</td><td>0.3203</td><td>23.504</td><td>0.99</td></tr><tr><td>Curve 2 (Moderately non-linear)</td><td>0.05</td><td>-0.005</td><td>1.393</td><td>6.318</td><td>0.99</td></tr><tr><td>Curve 3 (Slightly non-linear)</td><td>0.05</td><td>-0.005</td><td>2.0899</td><td>4.571</td><td>0.98</td></tr></tbody></table>

Figure 6 presents the stress-stretch curves for Interval 1 (Slightly Restrictive). Plot (a) corresponds to the highly non-linear curve, (b) represents the moderately non-linear curve, and (c) shows the slightly non-linear curve. It illustrates the correlation between the experimental stress and the predicted stress. The high R² values (ranging from 0.97 to 0.99) demonstrate an excellent fit between the experimental and calculated curve, confirming the accuracy of the optimization process.

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Stress-stretch tensile curves (highly non-linear (a), moderately non-linear (b), and slightly non-linear (c)).</figcaption></figure>

The three datasets fitted in this section correspond to the literature tensile curves shown previously in Figure 3 obtained from [13,14]. Stresses in directions 2 ($\sigma_2$) and 3 ($\sigma_3$) are expected to remain zero, as required by the boundary conditions of a uniaxial tensile test, confirming that the parameter set obtained from Interval 1 produces physically meaningful and reliable results across all directions.

## 3.2. Parameter Optimization: Inflation Test

Results from the inflation test are visually represented in Figure 5b, where three corneas (Cornea 1, Cornea 2, and Cornea 3) are displayed. These curves reflect a typical stress-stretch behavior, showing increasing stiffness as stretch progresses. It can be observed that the sudden non-linear behavior, or “knee,” occurs at different points along the curves, with some corneas showing an earlier rapid increase in stiffness, while others present a more gradual transition.

Cornea 1 exhibits a gradual increase in stress at lower strains, followed by a sharper rise later on. This indicates that Cornea 1 shows delayed fiber stiffening, as it is dominated by a lower initial stiffness, influenced mainly by a lower $c_1$. Cornea 2, on the other hand, displays a more intermediate response, with a smoother balance between early-stage stiffness and later non-linear behavior. Cornea 3 presents the sharpest increase in stress almost from the beginning of the test, which suggests that Cornea 3 has a higher initial
stiffness (likely due to a higher $c_1$) and reaches the non-linear fiber-dominated region earlier than the others.

As in Section 3.1, the parameter optimization for the inflation test is carried out for three different intervals, from which the material parameters are optimized using the differential evolution algorithm. A detailed analysis of the set of parameters obtained for each interval and the strong correlation with the experiments are shown in Table 2, with $R^2$ values close to 1.

<table><caption>**Table 2.** Estimated HGO material parameters from the inflation tests for three different optimization intervals of the three curves.</caption><thead><tr><th>Interval</th><th>Curve</th><th>c<sub>1</sub><br>(MPa)</th><th>c<sub>2</sub><br>(MPa)</th><th>k<sub>1</sub><br>(MPa)</th><th>k<sub>2</sub><br>(-)</th><th>R<sup>2</sup></th></tr></thead><tbody><tr><td rowspan="3">Interval 1 (Slightly Restrictive):<br>[0, 50], [0, 50], [0, 100], [0, 1000]</td><td>Curve 1</td><td>0</td><td>0.006</td><td>0.0243</td><td>223.66</td><td>0.99</td></tr><tr><td>Curve 2</td><td>0</td><td>0</td><td>0.0601</td><td>229.26</td><td>0.99</td></tr><tr><td>Curve 3</td><td>0</td><td>0</td><td>0.0717</td><td>434.87</td><td>0.99</td></tr><tr><td rowspan="3">Interval 2 (Non-Restrictive):<br>[-100, 100], [-100, 100], [0, 100],<br>[0, 1000]</td><td>Curve 1</td><td>-0.60</td><td>0.6054</td><td>0.0054</td><td>349.09</td><td>0.99</td></tr><tr><td>Curve 2</td><td>2.905</td><td>-3.792</td><td>1.393</td><td>28.87</td><td>0.99</td></tr><tr><td>Curve 3</td><td>1.494</td><td>-1.632</td><td>0.2917</td><td>216.23</td><td>0.99</td></tr><tr><td rowspan="3">Interval 3 (Highly Restrictive):<br>[0.005, 0.05], [-0.05, -0.005], [0,<br>100], [0, 1000]</td><td>Curve 1</td><td>0.0109</td><td>-0.005</td><td>0.0251</td><td>221.07</td><td>0.99</td></tr><tr><td>Curve 2</td><td>0.005</td><td>-0.05</td><td>0.1092</td><td>166.47</td><td>0.99</td></tr><tr><td>Curve 3</td><td>0.005</td><td>-0.05</td><td>0.1204</td><td>333.22</td><td>0.99</td></tr></tbody></table>

Interval 1 (Slightly Restrictive) emerges as the most suitable for parameter estimation due to its ability to provide physically meaningful values without encountering negative or unrealistic results. However, as shown in Table 2, all curves in this interval exhibit zero values for both $c_1$ and $c_2$, effectively neglecting the isotropic contribution to the cornea's material behavior. This simplification suggests that the model primarily captures anisotropic stiffening due to fiber alignment, which may not fully represent the cornea's initial stiffness, as the isotropic matrix plays a critical role in early deformation responses. This limitation highlights the need for further refinement in optimizing $c_1$, as addressed in Section 3.3.

In contrast, Interval 2 (Non-Restrictive) allows for a broader exploration of parameter values but yields physically invalid results. For example, Cornea 2 shows $c_1 = 2.905$ MPa and $c_2 = -3.792$ MPa, resulting in a negative shear modulus ($c_1 + c_2$), which is inconsistent with the physical interpretation of the corneal material. Despite achieving high $R^2$ values, these results lack physical coherence and cannot reliably describe the cornea's mechanical properties.

Interval 3 (Highly Restrictive), while yielding physically valid results, imposes excessively tight constraints that limit the optimization process. Parameters like $c_1$ are consistently forced to boundary values, such as $c_1 = 0.005$ MPa, restricting the model's ability to explore the cornea's true mechanical behavior. Although this interval may be useful in specific scenarios, its restrictive nature makes it less suitable for general applications.

Figure 7 illustrates the stress-stretch curves for the inflation test using parameters from Interval 1, with results presented for highly non-linear (a), moderately non-linear (b), and slightly non-linear (c) cases. The model accurately captures the anisotropic mechanical behavior of the cornea under inflation conditions. In plots (a, b, c), the predicted stress aligns closely with the experimental data, achieving high $R^2$ values ranging from 0.9967 to 0.9993. The circumferential and meridional stresses are equivalent, consistent with the equibiaxial tension assumption, while the radial stress remains near zero, as expected under thin membrane conditions. These findings confirm that the parameter set derived from Interval 1 provides a reliable and consistent representation of the cornea's mechanical response.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Stress-stretch inflation curves (highly non-linear (a), moderately non-linear (b), and slightly non-linear (c)).</figcaption></figure>

### 3.3. *Novel Methodology for Parameter Optimization in Corneal Inflation Tests*

During the optimization process for the inflation test, it became clear that the values of $c_1$ were unrealistically low, being calculated as zero in some configurations. This is physically incorrect because $c_1$ is responsible for the isotropic stiffness of the corneal matrix, which plays a crucial role in the initial response to deformation. To address this issue, a strategy is implemented to refine the model and ensure more accurate values for $c_1$.

The initial portion of the stress-stretch curves, where the response is approximately linear and dominated by the isotropic matrix, is isolated and optimized separately. This step allowed the model to focus on adjusting $c_1$, ensuring that the early stiffness of the cornea is accurately captured. The exact cutoff for this linear segment can vary between curves, but a practical guideline is that collagen fibers begin to bear a significant share of the load at roughly 2% strain [20]. Therefore, the 0–2% strain range can be treated as approximately linear and dominated by the isotropic matrix.

After isolating the initial phase of the curves, strict limits close to zero are imposed on $k_1$, and $k_2$ to refine the optimization of $c_1$. It is important to remember that $c_2$ can be neglected to simplify the material model, reducing it to a system with only three parameters to optimize. The final values for $c_1$ are in Table 3, demonstrating a clear increasing trend from Cornea 1 to Cornea 3. This validates the method, as the graphs show a corresponding increase in initial stiffness, aligning with the expected mechanical behavior.

<table><caption>**Table 3.** Isotropic parameter $c_1$ estimated values for the three curves.</caption><thead><tr><th>Parameter</th><th>Curve</th><th>c<sub>1</sub><br>(MPa)</th></tr></thead><tbody><tr><td rowspan="3">Interval [0, 50]</td><td>Curve 1</td><td>0.0198</td></tr><tr><td>Curve 2</td><td>0.0361</td></tr><tr><td>Curve 3</td><td>0.0624</td></tr></tbody></table>

Once $c_1$ is fixed, the anisotropic parameters $k_1$ and $k_2$ are optimized to capture the non-linear behavior of the cornea. The resulting values reflecting the influence of fiber alignment and stiffening in the later stages of the stress-strain response are shown in Table 4:

<table><caption>**Table 4.** Anisotropic parameters $k_1$ and $k_2$ estimated values, with fixed $c_1$ value.</caption><thead><tr><th>Interval<br>(MPa)</th><th>Curve</th><th>c<sub>1</sub><br>(MPa)</th><th>k<sub>1</sub><br>(MPa)</th><th>k<sub>2</sub><br>(-)</th><th>R<sup>2</sup></th></tr></thead><tbody><tr><td rowspan="3">Interval 1 (Slightly Restrictive):<br>[0, 50], [0, 50], [0, 100], [0, 1000]</td><td>Curve 1</td><td>0.0198</td><td>0.0143</td><td>273.48</td><td>0.99</td></tr><tr><td>Curve 2</td><td>0.0361</td><td>0.0299</td><td>311.86</td><td>0.99</td></tr><tr><td>Curve 3</td><td>0.0624</td><td>0.0226</td><td>691.26</td><td>0.99</td></tr></tbody></table>
Upon analyzing $k_1$, there appears to be no clear trend across the curves (it increases from Curve 1 to Curve 2 but then decreases in Curve 3). This pattern suggests that $k_1$ does not solely govern the non-linear behavior of the cornea and likely interacts with both $c_1$ and $k_2$, complicating its interpretation. In contrast, $k_2$ shows a clear increasing trend from Curve 1 to Curve 3, playing a more dominant role in the non-linear stiffening behavior observed at higher strains. This makes sense, as $k_2$ is closely linked to fiber stiffening. The overlap of $k_1$'s effects with both $c_1$ and $k_2$ complicates the optimization, reinforcing the importance of first optimizing $c_1$ to separate isotropic behavior from fiber-driven non-linearity. Once $c_1$ is fixed, $k_1$ and $k_2$ can be optimized more effectively.

In summary, $c_1$ increases steadily from Curve 1 to Curve 3, reflecting greater initial stiffness. $k_1$ follows a non-monotonic pattern, with no clear trend, while $k_2$ consistently rises, capturing the sharp increase in stiffness due to fiber alignment at larger strains, aligning with typical non-linear stiffening seen in the cornea during inflation tests.

## 4. Discussion

The significant variability observed in corneal biomechanical studies, particularly in tensile and inflation tests, arises from various causes. In this research, the authors highlight that these discrepancies mainly occur because different parameter sets can still make physical sense and achieve strong correlations with experimental data. This introduces challenges in defining a universally accepted set of material parameters for corneal tissues.

In tensile tests, a wide range of parameter values is reported across different studies (see Table 5). The parameters $c_1$, $c_2$, $k_1$, and $k_2$ show considerable variability, reflecting the anisotropic stiffening behavior of the cornea under uniaxial loading conditions. For example, $c_1$ values range from 0.0025 MPa [13] to 15.11 MPa [10] (Hoeltzel et al., 1992), illustrating a great difference in the isotropic stiffness parameter. Similarly, $k_1$ spans from 0.0338 MPa [15] to 103.5 MPa [12] (M. Nambiar et al., 2022), showing a wide range in fiber alignment's contribution to corneal stiffness. $k_2$ values vary from 5.74966 [13] to 1280.07 [15], highlighting a massive disparity in scale in the non-linear stiffening behavior of collagen fibers across the studies. These variabilities may also contribute to geometrical differences in the corneal specimens, testing different corneal regions (central vs. peripheral), variations in experimental conditions, and the methods used to model the anisotropic behavior of the tissue.

<table><caption>**Table 5.** Corneal parameter data in the scientific literature collected from tensile tests.</caption><thead><tr><th>Reference</th><th>c<sub>1</sub><br>(MPa)</th><th>c<sub>2</sub><br>(MPa)</th><th>k<sub>1</sub><br>(MPa)</th><th>k<sub>2</sub><br>(-)</th></tr></thead><tbody><tr><td>HGO [11]</td><td>0.3312</td><td>(-)</td><td>60.4949</td><td>114.6006</td></tr><tr><td>HGO Pandolfi [11]</td><td>0.3312</td><td>(-)</td><td>2.9211</td><td>20.7395</td></tr><tr><td>HGO [12]<br>[13]</td><td>0.046<br>0.0025</td><td>(-)<br>(-)</td><td>103.5<br>4.115088</td><td>149.5<br>5.74966</td></tr><tr><td>M. Nambiar et al., 2022 [12]</td><td>0.021</td><td>(-)</td><td>87.93</td><td>15.2</td></tr><tr><td>Wollensak et al., 2003 [10]</td><td>0.11</td><td>-0.01</td><td>0.055</td><td>14</td></tr><tr><td>Bryant et al., 1994 [10]</td><td>2.11</td><td>-2.109</td><td>0.2</td><td>50</td></tr><tr><td>Hoeltzel et al., 1992 [10]</td><td>15.11</td><td>-15.1</td><td>0.175</td><td>15</td></tr><tr><td>Zeng et al., 2001 [10]</td><td>0.0101</td><td>-0.01</td><td>0.055</td><td>16</td></tr><tr><td>Wollensak et al., 2003 (MR)<br>[10]</td><td>1.11</td><td>-1</td><td>(-)</td><td>(-)</td></tr><tr><td>Hoeltzel et al., 1992 (MR) [10]</td><td>20.01</td><td>-20</td><td>(-)</td><td>(-)</td></tr></tbody></table>

To aid comparison across sources, Table 5 (tensile) and Table 6 (inflation) use (-) to indicate values not reported. The differences across studies reflect model choices. Some studies use an HGO formulation with a two-parameter Mooney-Rivlin isotropic base ($c_1$
and $c_2$), others report only a single isotropic parameter ($c_1$), and some fit their curves with a pure two-parameter Mooney-Rivlin model without fiber terms $k_1$ and $k_2$. For Table 6, note that [15] present a re-fit in which $k_1$ and $k_2$ are recalculated after fixing $c_1$ to the mean across specimens (0.004 MPa).

<table><caption>**Table 6.** Corneal parameter data in the scientific literature collected from inflation tests.</caption><thead><tr><th>Reference</th><th>Specimen</th><th>c<sub>1</sub> (MPa)</th><th>c<sub>2</sub> (MPa)</th><th>k<sub>1</sub> (MPa)</th><th>k<sub>2</sub>(-)</th><th>c<sub>1</sub> (Mean) (MPa)</th><th>k<sub>1</sub> (New) (MPa)</th><th>k<sub>2</sub> (New)</th></tr></thead><tbody><tr><td rowspan="10">[15]</td><td>1</td><td>0.00291</td><td>(-)</td><td>0.0338</td><td>291.28</td><td>0.004</td><td>0.0339</td><td>294.41</td></tr><tr><td>2</td><td>0.00357</td><td>(-)</td><td>0.0219</td><td>362.33</td><td>0.004</td><td>0.0215</td><td>368.2</td></tr><tr><td>3</td><td>0.04691</td><td>(-)</td><td>0.0324</td><td>1080.51</td><td>0.004</td><td>0.1054</td><td>409.55</td></tr><tr><td>4</td><td>0.03744</td><td>(-)</td><td>0.0217</td><td>1280.07</td><td>0.004</td><td>0.0758</td><td>518.76</td></tr><tr><td>5</td><td>0.00275</td><td>(-)</td><td>0.068</td><td>435.38</td><td>0.004</td><td>0.0684</td><td>442.72</td></tr><tr><td>6</td><td>0.00906</td><td>(-)</td><td>0.1031</td><td>406.51</td><td>0.004</td><td>0.1049</td><td>376.44</td></tr><tr><td>7</td><td>0.01605</td><td>(-)</td><td>0.1602</td><td>774.71</td><td>0.004</td><td>0.1607</td><td>664</td></tr><tr><td>8</td><td>0.00425</td><td>(-)</td><td>0.1614</td><td>667.85</td><td>0.004</td><td>0.1609</td><td>667.02</td></tr><tr><td>9</td><td>0.00328</td><td>(-)</td><td>0.1482</td><td>605.98</td><td>0.004</td><td>0.1496</td><td>607.08</td></tr><tr><td>Mean</td><td>0.01402</td><td>(-)</td><td>0.0834</td><td>656.07</td><td>0.004</td><td>0.0979</td><td>483.13</td></tr><tr><td rowspan="3">[14]</td><td>A</td><td>0.05</td><td>(-)</td><td>25</td><td>2490</td><td>(-)</td><td>(-)</td><td>(-)</td></tr><tr><td>B</td><td>0.05</td><td>(-)</td><td>60</td><td>2490</td><td>(-)</td><td>(-)</td><td>(-)</td></tr><tr><td>C</td><td>0.05</td><td>(-)</td><td>130.9</td><td>2490</td><td>(-)</td><td>(-)</td><td>(-)</td></tr><tr><td rowspan="3">[10]</td><td>Healthy (Inflation)</td><td>0.32</td><td>-0.25</td><td>0.055</td><td>750</td><td>(-)</td><td>(-)</td><td>(-)</td></tr><tr><td>Keratoconus (periphery) (Inflation)</td><td>0.32</td><td>-0.25</td><td>0.005</td><td>400</td><td>(-)</td><td>(-)</td><td>(-)</td></tr><tr><td>Keratoconus (center) (Inflation)</td><td>0.32</td><td>-0.25</td><td>0.005</td><td>200</td><td>(-)</td><td>(-)</td><td>(-)</td></tr></tbody></table>

Inflation tests present distinct challenges. Variability in inflation studies is also very high (see Table 6). Apart from the previously mentioned possible sources of variability, an additional source arises from converting pressure-apical displacement curves into stress-stretch curves, which relies heavily on geometrical measurements such as the corneal base radius, thickness, and apex height. If these are measured inaccurately, they can drastically alter the resulting curves. Consequently, variation in material parameters in inflation tests stems from differences in these geometric inputs, in addition to model assumptions and experimental conditions. For instance, $c_1$ values in inflation tests range from 0.00291 MPa [15] to 0.32 MPa [10]. Similarly, $k_1$ values vary between 0.0217 MPa [15] and 0.055 MPa [10], and $k_2$ ranges from 291.28 [15] to 750 [10].

Tensile and inflation test parameters are not directly comparable. The mechanical behavior of the cornea under uniaxial tensile loading (tensile test) is characterized by higher stresses and strains compared to the inflation test, which simulates physiological conditions with lower stresses but captures the biaxial response of the cornea. Tensile test curves reach stress values up to 7.5 MPa and stretches up to 1.18, indicating a stronger response to uniaxial deformation. In contrast, inflation test curves display much lower stress levels, peaking around 0.05 MPa, with stretches reaching only 1.045. These differences explain why optimized parameters obtained from tensile and inflation tests are not interchangeable. For example, $k_2$ values in inflation tests are typically much larger than those in tensile tests.

### Limitations of the Study

This study is subject to some limitations that must be acknowledged. First, the datasets used were limited in number and derived from previously published experimental studies, which may not fully capture the biological variability of corneal tissue. In addition, the study relied on representative smoothed stress-stretch curves from literature instead of raw
experimental data, which limits the ability to evaluate the robustness of the optimization
process against real experimental noise and variability.

Regarding the modeling assumptions, the inflation test was treated as a spherical
shell with uniform thickness, neglecting regional variations in geometry and thickness
distribution. While this simplification is common in corneal biomechanics, it inevitably
introduces discrepancies with physiological reality.

A particularly relevant limitation concerns the transformation of pressure–apical dis-
placement curves into stress–stretch curves in inflation experiments. This conversion is
highly sensitive to the assumed geometric input parameters, such as corneal thickness,
apex height, and base radius. Small variations or inaccuracies in these values can signif-
icantly alter the computed stress–stretch relationships, thereby affecting the optimized
material parameters.

The present analysis also adopts a simplified, general model with fixed geometric in-
puts rather than subject-specific anatomy; this choice facilitates comparison across datasets
but does not capture inter-individual variability in contour, stiffness, and thickness. Finally,
the stepwise sequential approach was applied and validated only for the HGO model;
although it may extend to models that clearly separate isotropic and anisotropic contribu-
tions with few parameters, its applicability to highly coupled multi-parameter laws (e.g.,
five-parameter Mooney–Rivlin) is limited.

To mitigate these limitations, future work may use raw experimental data rather than
smoothed curves; include subject-specific geometry from imaging (apex height, base radius,
thickness maps); replace the axisymmetric membrane with a biaxial, non-axisymmetric for-
mulation; add fiber-dispersion and allow non-orthogonal preferred directions; and validate
the estimated parameters on independent datasets and in finite-element simulations under
physiological intraocular pressure.

**5. Conclusions**

The objectives of this study are to optimize the mechanical parameters of the cornea using both tensile and inflation tests and to address the variability in parameters reported across studies. It is established that different testing methods (tensile and inflation) yield significantly different results in the HGO model parameters *c*₁, *k*₁, and *k*₂, excluding *c*₂ to simplify the model. To achieve this, the parameters were optimized across different optimization ranges for three tensile curves with varying degrees of nonlinearity and three inflation curves.

A key factor contributing to the variability across studies includes differences in the
specific corneal zones that are tested, and initial geometric measurements like base, height,
and radius in the case of inflation tests. Mathematical simplifications in virtual tests also
play a role in this variability. However, it is identified and demonstrated as the primary
source of variation that multiple sets of parameters can highly correlate with experimental
stress-stretch curves while still making physical sense. This presents a challenge in reducing
parameter variability.

This study proposes a stepwise approach to mitigate this variability. First, *c*₁ should be optimized based on initial deformations from tensile tests, where isotropic stiffness dominates. Once *c*₁ is fixed, the anisotropic parameters *k*₁ and *k*₂, which capture fiber alignment and non-linear stiffening, may be optimized. This method simplifies the optimization process and reduces variability in parameter estimation.

Through a comparative analysis, the study demonstrates that tensile and inflation
test data are not interchangeable. Each test yields different parameter values, as the stress-
stretch scales differ significantly.
**Author Contributions:** Conceptualization J.G.-C., C.G., and F.C.; methodology J.G.-C., C.G., and F.C.; software J.G.-C. and C.G.; validation M.P.; formal analysis J.G.-C. and M.P.; investigation J.G.-C., C.G., and F.C.; resources M.P.; data curation J.G.-C., M.P., and C.G.; writing—original draft preparation J.G.-C., C.G., M.P., and F.C.; writing—review and editing, M.P. and F.C.; visualization M.P. and F.C.; supervision M.P. and F.C.; project administration J.G.-C., C.G., and F.C.; funding acquisition C.G. and F.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This publication is part of the research project 22665/PI/24, funded by grant FS/10.13039/ 100007801, and was also carried out in collaboration with the Association for the Integration of the Disabled in the Comarca del Mar Menor (AIDEMAR, collaboration protocol UPCT-AIDEMAR).

**Institutional Review Board Statement:** This study was conducted according to the guidelines of the Declaration of Helsinki and approved by the Institutional Ethics Institutional Review Board Statement: Committee of the Technical University of Cartagena (CEI21_001). Informed consent was obtained from all the subjects involved in this study.

**Informed Consent Statement:** Not applicable. The study is based exclusively on previously published data, and no new data involving human participants were collected.

**Data Availability Statement:** Data will be made available upon request.

**Acknowledgments:** DeepSeek version 3.1 was used for translating the original manuscript from Spanish to American English.

**Conflicts of Interest:** The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

List of Acronyms

<table><tr><td>HGO</td><td>Holzapfel-Gasser-Ogden</td></tr><tr><td>a</td><td>Base radius of the cornea</td></tr><tr><td>c<sub>i</sub>(i = 1,2)</td><td>Material constant of the Mooney-Rivlin model</td></tr><tr><td>h</td><td>height of the dome</td></tr><tr><td>h<sub>0</sub></td><td>Initial height of the dome</td></tr><tr><td>k<sub>1</sub></td><td>Hyperelastic material parameter to consider the influence of collagen fibrils' stiffness</td></tr><tr><td>k<sub>2</sub></td><td>Hyperelastic material parameter to consider the influence of collagen fibrils on large extensions</td></tr><tr><td>J</td><td>Jacobian of the deformation gradient tensor</td></tr><tr><td>b</td><td>Left Cauchy-Green deformation tensor</td></tr><tr><td>C</td><td>Right Cauchy-Green deformation tensor</td></tr><tr><td>e<sub>R</sub></td><td>Unidirectional vector in radial direction</td></tr><tr><td>e<sub>φ</sub></td><td>Unidirectional vector in meridional direction</td></tr><tr><td>e<sub>θ</sub></td><td>Unidirectional vector in circumferential direction</td></tr><tr><td>F</td><td>Deformation gradient tensor</td></tr><tr><td>m<sub>0</sub></td><td>Initial preferential orientation of collagen fibril family 1</td></tr><tr><td>m</td><td>Deformed preferential orientation of collagen fibril family 1</td></tr><tr><td>n<sub>0</sub></td><td>Initial preferential orientation of collagen fibril family 2</td></tr><tr><td>n</td><td>Deformed preferential orientation of collagen fibril family 2</td></tr><tr><td>P</td><td>Pressure</td></tr><tr><td>R</td><td>Radius of curvature</td></tr><tr><td>R<sup>2</sup></td><td>Coefficient of determination</td></tr><tr><td>SS<sub>Res</sub></td><td>Sum of squared errors between the predicted and experimental values</td></tr><tr><td>SS<sub>Tot</sub></td><td>Sum of square deviations of the experimental values and their average</td></tr></table>
<table><tr><td>t</td><td>Cornea's thickness</td></tr><tr><td>t<sub>0</sub></td><td>Initial Cornea's thickness before pressurization</td></tr><tr><td>x</td><td>Deformed vector</td></tr><tr><td>X</td><td>Reference vector</td></tr><tr><td>λ<sub>i</sub> (i = 1, 3)</td><td>Stretch ratio component</td></tr><tr><td>σ<sub>i</sub> (i = 1, 3)</td><td>Cauchy stress component</td></tr><tr><td>ψ</td><td>Strain energy density</td></tr><tr><td>ψ<sub>deviatoric</sub></td><td>Deviatoric strain energy density</td></tr><tr><td>I&#773;<sub>i</sub> (i = 1, 2, 4, and 6)</td><td>Invariants of the Right Cauchy Stress Tensor</td></tr></table>

## Appendix A

To derive the stress-stretch relationship according to [22], the applied pressure is converted into the appropriate units in Pa.

$$P_{(Pa)} = P_{(mmHg)} \times 133.322. \qquad (A1)$$

The cornea is approximated as a spherical shell, assuming that the applied pressure generates an apical displacement and homogenous surface deformation. The radius of curvature *R* is calculated as

$$R = \frac{(a^2 + h^2)}{2h}, \qquad (A2)$$

$$h = h_0 + w, \qquad (A3)$$

where *h*₀ is the initial height of the dome, *w* is the apical displacement, and *a* is the base radius of the cornea. This formula assumes the cornea deforms as a spherical cap.

The stretch *λ* is calculated from the deformed radius of curvature *R* and the initial radius of curvature before pressurization *R*₀.

$$\lambda = \frac{\left( R \sin^{-1}(a/R) \right)}{\left( R_0 \sin^{-1}(a/R_0) \right)}. \qquad (A4)$$

The stress in the cornea can be approximated using Laplace's law, which relates the applied pressure *P*<sub>(Pa)</sub>, deformed radius *R*, and the deformed corneal thickness *t*:

$$\sigma = \frac{P_{(Pa)}R}{2t}, \qquad (A5)$$

$$t = \frac{t_0}{\lambda^2}, \qquad (A6)$$

where *t*₀ is the initial thickness of the cornea before pressurization.

## References

1. Burton, M.J.; Ramke, J.; Marques, A.P.; Bourne, R.R.A.; Congdon, N.; Jones, I.; Tong, B.A.M.A.; Arunga, S.; Bachani, D.; Bascaran, C.; et al. The Lancet Global Health Commission on Global Eye Health: Vision beyond 2020. *Lancet Glob. Health* **2021**, *9*, e489–e551. [CrossRef] [PubMed]

2. Enaholo, E.S.; Musa, M.J.; Zeppieri, M. The Spherical Equivalent. In *StatPearls*; StatPearls Publishing: Treasure Island, FL, USA, 2024.

3. Alkanaan, A.; Barsotti, R.; Kirat, O.; Khan, A.; Almubrad, T.; Akhtar, S. Collagen fibrils and proteoglycans of peripheral and central stroma of the keratoconus cornea: Ultrastructure and 3D transmission electron tomography. *Sci. Rep.* **2019**, *9*, 19963. [CrossRef] [PubMed]
4. Liu, T.; Shen, M.; Li, H.; Zhang, Y.; Mu, B.; Zhao, X.; Wang, Y. Changes and quantitative characterization of hyper-viscoelastic biomechanical properties for young corneal stroma after standard corneal cross-linking treatment with different ultraviolet-A energies. *Acta Biomater.* **2020**, *113*, 438–451. [CrossRef] [PubMed]

5. Wang, J.; Liu, X.Y.; Bao, F.J.; Lopes, B.T.; Wang, L.Z.; Eliasy, A.; Abass, A.; Elsheikh, A. Review of ex-vivo characterization of corneal biomechanics. *Med. Nov. Technol. Devices* **2021**, *11*, 100074. [CrossRef]

6. Santodomingo-Rubido, J.; Carracedo, G.; Suzaki, A.; Villa-Collar, C.; Vincent, S.J.; Wolffsohn, J.S. Keratoconus: An updated review. *Contact Lens Anterior Eye* **2022**, *45*, 101559. [CrossRef] [PubMed]

7. Holzapfel, G.A.; Gasser, T.C.; Ogden, R.W. A new constitutive framework for arterial wall mechanics and a comparative study of material models. *J. Elast. Phys. Sci. Solids* **2000**, *61*, 1–48. [CrossRef]

8. Gómez, C.; Piñero, D.P.; Peredes, M.; Alió, J.L.; Cavas, F. Iterative methods for the biomechanical evaluation of corneal response: A case study in the measurement phase. *Appl. Sci.* **2021**, *11*, 10819. [CrossRef]

9. Pandolfi, A. Cornea modelling. *Eye Vis.* **2020**, *7*, 2. [CrossRef] [PubMed]

10. Pandolfi, A.; Manganiello, F. A model for the human cornea: Constitutive formulation and numerical analysis. *Biomech. Model. Mechanobiol.* **2006**, *5*, 237–246. [CrossRef] [PubMed]

11. Liu, T.; Shen, M.; Huang, L.; Xiang, Y.; Li, H.; Zhang, Y.; Wang, Y. Characterization of hyperelastic mechanical properties for youth corneal anterior central stroma based on collagen fibril crimping constitutive model. *J. Mech. Behav. Biomed. Mater.* **2020**, *103*, 103575. [CrossRef] [PubMed]

12. Nambiar, M.H.; Liechti, L.; Studer, H.; Roy, A.S.; Seiler, T.G.; Büchler, P. Patient-specific finite element analysis of human corneal lenticules: An experimental and numerical study. *J. Mech. Behav. Biomed. Mater.* **2023**, *147*, 106141. [CrossRef] [PubMed]

13. Lanchares, E. Modelado Biomecánico de los Componentes Refractivos del Ojo Humano y Tratamientos Refractivos Asociados. Ph.D. Thesis, University of Zaragoza, Zaragoza, Spain, 2010.

14. Ariza-Gracia, M.Á. Methods for Characterising Patient-Specific Corneal Biomechanics. Ph.D. Thesis, University of Zaragoza, Zaragoza, Spain, 2018.

15. Kok, S.; Botha, N.; Inglis, H.M. Calibrating corneal material model parameters using only inflation data: An ill-posed problem. *Int. J. Numer. Methods Biomed. Eng.* **2014**, *30*, 1460–1475. [CrossRef] [PubMed]

16. Elsheikh, A.; Anderson, K. Comparative study of corneal strip extensometry and inflation tests. *J. R. Soc. Interface* **2005**, *2*, 177–185. [CrossRef] [PubMed]

17. Studer, H.; Riedwyl, H.; Büchler, P. Importance of multiple loading scenarios for the identification of material coefficients of the human cornea. *Comput. Methods Biomech. Biomed. Eng.* **2011**, *15*, 93–99. [CrossRef] [PubMed]

18. Scarcelli, G.; Besner, S.; Pineda, R.; Yun, S.H. Biomechanical characterization of keratoconus corneas ex vivo with Brillouin microscopy. *Investig. Ophthalmol. Vis. Sci.* **2014**, *55*, 4490–4495. [CrossRef] [PubMed]

19. Foong, T.Y.; Hua, Y.; Amini, R.; Sigal, I.A. Who bears the load? IOP-induced collagen fiber recruitment over the corneoscleral shell. *Exp. Eye Res.* **2023**, *230*, 109446. [CrossRef] [PubMed]

20. Hatami-Marbini, H.; Emu, M.E. Biomechanical properties of porcine cornea; planar biaxial tests versus uniaxial tensile tests. *J. Mech. Behav. Biomed. Mater.* **2025**, *166*, 106955. [CrossRef] [PubMed]

21. Muñoz-Villaescusa, C.; Núñez-Chongo, O.d.l.C.; Cárdenas-Díaz, T.; Batista-Leyva, A.J.; Cavas-Martínez, F. Experimental Determination of Corneal Elastic Constants and Their Use in Biomechanical Modeling. *Appl. Sci.* **2021**, *11*, 11292. [CrossRef]

22. Roberts, C.J.; Dupps, W.J.; Downs, J.C. *Biomechanics of the Eye*; Kugler Publications: Amsterdam, The Netherlands, 2018.

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.