Optimization of Linear Filtering Model to Predict Post-LASIK Corneal Smoothing Based on Training Data Sets

Laser vision correction is a rapidly growing field for correcting nearsightedness, farsightedness as well as astigmatism with dominating laser-assisted in situ keratomileusis (LASIK) procedures. While the technique works well for correcting spherocylindrical aberrations, it does not fully correct high order aberrations (HOAs), in particular spherical aberration (SA), due to unexpected induction of HOAs post-surgery. Corneal epithelial remodeling was proposed as one source to account for such HOA induction process. This work proposes a dual-scale linear filtering kernel to model such a process. Several retrospective clinical data sets were used as training data sets to construct the model, with a downhill simplex algorithm to optimize the two free parameters of the kernel. The performance of the optimized kernel was testedon new clinical data sets that were not previously used for the optimization.


Introduction
Historically, eyeglasses and contact lenses have been used to alleviate refractive problems such as nearsightedness, farsightedness, and astigmatism.With the advent of excimer lasers [1] specially designed for laser-assisted in situ keratomeliusis (LASIK) and photorefractive keratectomy (PRK), patients started to enjoy a new type of vision correction that is free of eyeglasses.With wavefront-guided LASIK [2], the correction of ocular aberrations is no longer limited to the so-called low-order aberrations, i.e., the spherocylindrical error that can be corrected with traditional eyeglasses.This new technology enables the correction of higher-order aberrations (HOAs) that are beyond the spherocylindrical error, most notably spherical aberration and coma.Thus, super sharp vision is attainable in theory with the wavefront-guided LASIK.
Unfortunately, the human cornea is not a piece of plastic [3].With LASIK, it involves first cutting a flap on the corneal stroma, lifting it to the side, then delivering the UV laser pulses to remove tissue, and finally putting back the flap, which heals shortly after surgery.The precise design of an ablation target may cut the corneal stroma as needed to achieve a desired shape immediately after surgery.However, the biomechanical process and the corneal epithelial remodeling after surgery change the surface of the cornea, resulting in deviations from the original optical design of the ablation shape.Therefore, the post-operative induction of HOAs, especially spherecal aberration (SA), is currently among the most serious challenges for laser vision correction technology.Among several possible root causes of SA induction, the postoperative cornea remodeling was found the most important [4].The main effect of the cornea remodeling is the smoothing of epithelial anterior surface, when the epithelium tends to grow thicker at the center and fill in the dips of the cornea surface, created by refractive surgery [5].The epithelial smoothing causes some spherocylindrical regression after refractive surgery, which can be corrected by a linear adjustment of the intended refractive correction.It also leads to the induction of high-order aberrations, which are increasingly strong for high myopia and hyperopia cases [4,6].
Among the HOAs induced, spherical aberration is the most significant.In general, the amount of the SA induction tends to increase with post-surgery time.Several months after surgery when the cornea stabilizes, the in-duced SA shows a statistically significant trend versus the magnitude of the treated refraction.Figure 1 shows the post-operative SA over a 6 mm diameter as a function of the pre-operative manifest refraction in spherical equivalent (MRSE).The regression slope of the induction is remarkably consistent between different data sets.
The purpose of this study is to find a corneal smoothing model to represent the corneal change post-surgery using an optimization algorithm, based on retrospectively available clinical data.The kernel is then tested with other clinical data sets that were not previously used for the optimization.This well tested kernel can then be used to "reverse" the biological corneal smoothing effect by a mathematical deconvolution process.An improved treatment algorithm can then be designed, which hopefully will remove the induced spherical aberration.

Modeling of Post-Operative Corneal Smoothing
Various models can capture geometric changes to the surface of the human cornea occurring after the surgery.
We considered an optimized linear filter (OLF) model, which describes post-operative smoothing of the corneal ablation.This model is characterized by a small set of parameters determined by a model optimization based on retrospective clinical data.The post-operative epithelial smoothing process can be simulated by means of a simple mathematical model.This model defines the shape of the post-operative cornea surface as a convolution of the ablation target profile with a linear smoothing filter as where h stands for the elevation maps of the corneal surface for pre-operative and post-operative situations, respectively,  denotes a convolution operation, T(x,y) is the ablation target profile and K(x,y) is the linear smoothing filter kernel.A simple squared Butterworth low-pass filter [7] has been proposed.A squared Butterworth filter of the first order takes a form with the square term of the spatial frequency as where , and s is a parameter representing the scale of the kernel.With limited success using the squared Butterworth filter as defined by Equation ( 2), we began consideration of dual-scale and triple-scale OLFs that has a somewhat similar shape as the squared Butterworth filter.Our tests show that a dual-scale OLF model has the advantage of faster convergence and proper account of biological change of the epithelial cells than the triple-scale model.Therefore, we have used a dual-scale OLF kernel that is defined as is the radial distance from the coordinate origin, s 2 and s 4 are two unknown free parameters to be determined.For our application, r, s 2 and s 4 all have dimensions in mm. Figure 2 shows the cross-section of the kernel and its power spectrum.

Optimization of the Kernel Parameters
"obs" stands for observation, i.e., clinical outcome.
Several different optimization algorithms have been tested.Our choice of the downhill simplex method [8,9] works well with our model and the data sets.With four clinical data sets (two parts of Data Set 1, low myopia and high myopia, and Data Sets 3 and 5) used for optimization, we found that s 2 = 0.0334 mm and s 4 = 0.464 mm give the minimum σ as defined in Equation (4).By using Equation (3) in Equation ( 1) using the pre-operative and post-operative wavefront data as well as the treatment targets for various previously treated eyes, the two unknown parameters s 2 and s 4 can be obtained by minimizing the difference between the simulated postoperative wavefront error and the observed post-operative wavefront error.This minimization is a least-squares type which minimizes the regression slopes of the postoperative spherical equivalent (SE) and post-operative SA as a function of the pre-operative SE for all eyes as     With these two kernel parameters, application of the treatment parameters for those eyes using Equation (1) enables us to obtain the simulated clinical outcome, which includes spherocylindrical error (wavefront spherical equivalent, or WSE) and SA.The WSE is measured in diopters (D) and the SA is measured in microns (µm) over a 6 mm diameter.Figure 3 shows the comparison between the observed and simulated post-operative outcome for two data sets that were used for the optimization.Both the post-operative WSE and SA as a function of the pre-operative WRSE are plotted.With no surprise, the regression slopes of the simulated eyes agree well with those of the observed eyes.where slopeSE and slopeSA are the regression slopes of the post-operative SE versus pre-operative SE and postoperative SA versus pre-operative SE, respectively.δ stands for 95% confidence interval of the observed slope.Subscript "simu" stands for simulation and subscript Once the parameters s 2 and s 4 are determined, the OLF kernel can be determined based on Equation (3).To obtain a new target shape that is capable of removing the post-operative induction of spherical aberration, a deconvolution process of Equation ( 1) can be employed as where F(·) stands for a Fourier transform, * denotes a complex conjugate, T current is the current treatment target with induction of post-operative SA, T new is the new target that is expected to remove the post-operative SA, and K INV is the inverse kernel of K(x,y).This is the typical Wiener filtering technique [9].The SNR is used to prevent noise amplification and oscillation at the edge.A constant value of 0.1 was used for practical purpose.

Verification of the Model with New Test Data
The effect of post-LASIK central corneal thickening caused by epithelial smoothing has been observed previouslyin the literature [4,5,[10][11][12] and provided at least partial explanation for regression after refractive surgery for myopia.The OLF obtained in this study confirms the central corneal thickening phenomena, which gives biological support of the kernel.
For the first verification, we used the kernel optimized with Data Sets 1, 3 and 5 and applied it to Data Sets 2 and 4, as depicted in Figure 4.It is interesting to see that the regression slopes of the simulated eyes agree well with those of the observed eyes, even though these data sets were not used for the optimization.This result is expected because the model is supposed to simulate the post-operative corneal smoothing process, which should not be different for different data sets.
For the second verification, recall that we used SE and SA as two aberration parameters for optimization in Equation ( 4).More convincingly, we used the same kernel to obtain similar regression slopes for the secondary spherical aberration, which is not a parameter used in the optimization, as shown in Figure 5. Again, we have good matches between the observed and the simulated slopes.This result comes as expected, as the induction of HOAs from the corneal smoothing is primarily rotationally symmetric.Secondary spherical aberration is the most important rotationally symmetric after sphere and primary spherical aberrations.

Some regression plots show a constant offset between
Observed: y = 0.1146x -0.0436 R² = 0.1323 Simulated: y = 0.0894x -0.2302 R² = 0.7088 the simulated and the observed trend lines.These offsets for post-operative SE or SA trend are about the same for all pre-operative MRSE values, indicating that they do not depend on ablation depth.They may be caused by the creation of the LASIK flap [10][11][12][13].Depending on the choice of microkeratome and individual surgeon technique, the flap-induced aberrations may differ from site to site or surgeon to surgeon [14].

Discussion
Search of the smoothing kernel to model the postoperative induction of spherical aberration is a relatively new field of study.Based on the observation of post-operative refractive regression, a simple low-pass squared Butterworth filter was proposed [7].This kernel is defined by a single free parameter, which characterizes the scale of smoothing.Unfortunately, this model does not provide a satisfactory fitting for the regression slopes for both post-operative low-order refraction and high order aberrations simultaneously.Optimization for both refraction and spherical aberration leads to diverged outcomes.
The proposed OLF with two free parameters models a dual-scale smoothing.Looking at the cross-section of the kernel in logarithmic scale in Figure 2, the sharp core corresponds to the short scale diffusion process and the wide wings correspond to the long scale smoothing process.These two separated processes can be linked to the post-operative corneal change in low-order and highorder aberrations, respectively.Consequently, the OLF yields a good match for both low-order aberrations (WSE) and high order aberrations (SA) observed clinically.Furthermore, this match can be extended to different data sets and different aberration types (secondary spherical aberration).
Looking at the inverse kernel K INV , as depicted in Figure 6, we found that the peak of the power spectrum of the inverse kernel corresponds to the size of the superficial cells of the epithelium.This may not be just a coincidence, as the movement of the epithelial cells, especially the superficial cells, it's attributed to the mechanism of the post-operative corneal smoothing.As the smoothing kernels generally smooth high curvature areas, the effect of the inverse kernel works exactly the opposite, sharpening areas that have high curvature changes.The link of the peak of the power spectrum of the inverse kernel to the size of the superficial cells of the epithelium rovides another layer of support of our model.p Similar to the neural network technique where a data set is used to train the network and a new data set is used to test the network, we used combination of data sets to optimize the free parameters of the proposed kernel and used new data sets to test the effectiveness of the model.The successful test of the model strengthens is usefulness, which can be important when costly clinical trials are decided.

Figure 1 .
Figure 1.Post-LASIK spherical aberration (SA) as a function of the pre-operative manifest refraction in spherical equivalent (MRSE) for two data sets.

Figure 2 .
Figure 2. Cross-section of the optimized linear filter.Left panel, linear scale of the center of the kernel; middle panel, logaithmic scale of the entire kernel; right panel, power spectrum of the kernel.r

Figure 6 .
Figure 6.Cross-section of the inverse kernel of the OLF.Left panel, the inverse kernel; right panel, the power spectrum.