Selected parameters of the corneal deformation in the Corvis tonometer

Introduction Contemporary ophthalmology knows many methods of measuring intraocular pressure, namely the methods of non-contact and impression applanation tonometry. In non-contact applanation tonometers, e.g. the Corvis, the corneal flattening is caused by an air puff. Image registration of the corneal deflection performed by a tonometer enables to determine other interesting biomechanical parameters of the eye, which are not available in the tonometer. The measurement of new selected parameters is presented in this paper. Material and method Images with an M × N × I resolution of 200 × 576 × 140 pixels were acquired from the Corvis device in the source recording format *.cst. A total of 13'400 2D images of patients examined routinely in the Clinical Department of Ophthalmology, in District Railway Hospital in Katowice, Poland, were analysed in accordance with the Declaration of Helsinki. A new method has been proposed for the analysis of corneal deflection images in the Corvis tonometer with the use of the Canny edge detection method, mathematical morphology methods and context-free operations. Results The resulting image analysis tool allows determination of the response of the cornea and the entire eyeball to an air puff. The paper presents the method that enables the measurement of the amplitude of curvature changes in the frequency range from 150 to 500 Hz and automatic designation of the eyeball movement direction. The analysis of these data resulted in 3 new features of dynamics of the eye reaction to an air puff. Classification of these features enabled to propose 4 classes of deformation. The proposed algorithm allows to obtain reproducible results fully automatically at a time of 5 s per patient using the Core i5 CPU M460 @ 2.5GHz 4GB of RAM. Conclusions The paper presents the possibility of using a profiled algorithm of image analysis, proposed by the authors, to measure additional cornea deformation parameters. The new tool enables automatic measurement of the additional new parameters when using the Corvis tonometer. A detailed clinical examination based on this method will be presented in subsequent papers.


Introduction
Contemporary ophthalmology knows many methods of measuring intraocular pressure. These include non-contact and impression applanation tonometry. They are based on the Imbert-Fick law. Knowing the physical strength necessary to flatten a sphere and the surface of the flattening, it is possible to determine the pressure inside the sphere. Non-contact applanation tonometers, e.g. the Corvis, use this rule. In this tonometer, the flattening of the cornea is caused by an air puff - Figure 1. The device measures the cornea thickness, deflection amplitude, applanation length, corneal deflection speed, intraocular pressure. In addition, the information about places of applanation, pachymetry and other data are provided - Figure 2. The locations of applanation (temporary cornea flattening during a convex-concave transition and vice versa) enable to determine two applanation times t 1 and t 2 . In Ocular Response Analyzer (ORA) the air pressure values P 1 and P 2 are read in applanation times t 1 and t 2 , and then the value (P 1 + P 2 )/2 is calculated. This value after linear calibration is the value of the intraocular pressure IOP. The Corvis tonometer also records sequence of images of the corneal deflection (resulting from the air puff). The sequence of images enables to get much more information using the methods of image analysis and processing. Such analysis allows to obtain three-dimensional reconstruction of corneal deflection in one of time ranges. On this basis, it is possible to determine the locations of applanation, the corneal deformation rate or the entire eyeball response to the air puff.
Analysis of images from the Corvis tonometer in the original software is based on designation of the outer and inner edges of the cornea profile. The edge detection methods of Roberts, Sobel and Prewitt can be all applied here. The Canny filter or local thresholding dependent on local brightness changes can be also used with equally good results. Additionally, the visible cross-sections of the iris which are often the brightest visible area are problematic here. All these methods based on a typical standard approach, enable to detect the corneal deformation in only one selected place (usually the main axis of fixation). From a biomechanical point of view, other locations of the Figure 1 Schematic diagram showing methodology of the measurement method and corneal response. The cornea is subjected to an air-puff stimulus in the Corvis tonometer. In consequence, an image of the corneal deflection in the line arranged on the main axis is created every 23 μs. Individual deflection phases of greater importance are shown on the right. corneal deformation or determination of the eyeball reaction to the air puff are also of interest. These parameters are not adequately described in the known publications in terms of both image analysis and processing as well as biomechanics.
The papers published so far are related to the analysis of corneal biomechanical parameters for the data derived mainly from the ORA device [1,2]. In [1], measurements of several populations indicate that corneal hysteresis, a biomechanical measure, varied over a dynamic range of 1.8 to 14.6 mmHg and was only weakly correlated with corneal thickness (r 2 = 0.12). This is related to the observation that some subjects with relatively thick corneas have less-than-average corneal hysteresis. For the cornea itself, the effect of the patient's age [3][4][5][6], glaucoma [7] or wound healing [8] on its biomechanical parameters is also considered. Predictive numerical simulation of corneal biomechanical parameters was considered in an equally interesting way in [9]. In this paper, twentyfour corneal buttons were tested under posterior inflation conditions while monitoring their behaviour using non-contact methods. The measurement methods of the ultrastructure of the corneal stroma were presented in [10]. There are also very interesting studies profiled to the cornea analysis of Brazilian patients [11] or children with the use of the ORA tonometer [12,13]. A group of authors in [14] presented a correlation between hysteresis generated in a tonometer and pachymetry. By contrast, in [15] the same hysteresis as in [16] was presented in correlation with the glaucoma damage. The paper [17] deals with these issues (hysteresis), but in the Reichert ORA device. There exists also a well-known group of papers concerning the analysis of biomechanical parameters of the cornea in keratoconic eyes [18] or LASIK [19]. Moreover, there are descriptions of other tests of the cornea which contribute to the understanding of the biomechanical properties such as the comparison of the corneal strip extensometry presented in [20].
Last year has appeared series of papers dealing with application of Corvis noncontact tonometer to examine both IOP value as well as biomechanical parameters of the cornea, related to the corneal deformation [21][22][23][24]. However, Corvis tonometer gives more possibilities of corneal measurements after application of additional algorithms and procedures of data processing.
Therefore, there is a need to propose a dedicated automatic algorithm enabling measurements performed at out of axis places of the corneal deformation or determination of the eyeball response to an air puff during pressure measurement. The proposed dedicated algorithm and the results obtained from its use are discussed in detail later in the paper.

Pre-processing
Pre-processing of input images L GRAY (m,n,i), where m-row m∈(1,M), n-column n∈(1,N), i another 2D image i∈(1,I), is related to their automatic reading from the source file with the extension *.cst. In the file, *.cst is a sequence of 2D images preceded by a header which contains basic information about the patient. A sequence of images L GRAY (m,n,i) thus obtained is then subjected to median filtering with a mask h 1 sized M h1 × N h1 × I h1 = 3 × 3 × 3 pixels. The mask size was chosen arbitrarily, taking into account the size of distortions and artefacts entering the optical path. In case of larger mask size h 1, median filtering was causing the deletion of the corneal contour, therefore the size M h1 × N h1 × I h1 = 3 × 3 × 3 pixels was optimal. The filtered image L M (m,n,i) undergoes further preliminary transformations. These include designation of the outer edge of the cornea. Two different independent approaches to determining the outer limit of the cornea were implemented.
The first one involves analysis of brightness for each column of the image L M (m,n = const,i = const). Each column is subjected to binarization operation with a threshold set automatically according to the Otsu formula [25], i.e.: where: p r (n,i) -binarization threshold set automatically according to the Otsu formula for m∈(1,M), v r -automatic correction factor chosen arbitrarily at 0.5. According to formula (1) binarization threshold is determined on the basis of two values: p r (n,i) which is the binarization threshold according to the Otsu and v r which was chosen arbitrarily. Factor v r is the correction factor resulting from the image specificity L M (m,n,i). There is a necessity of the correction of the threshold value p r (n,i) by 50% (0.5) due to high contrast of the corneal edges occurring in some cases. The resulting binary image L B (m,n,i) is analysed to enable determination of the outer edge of the cornea. One of the simplest methods is the calculation of the XOR value between two images. L B (m,n,i) is one image and the other one is the result of its morphological erosion L E (m,n,i) = min(L B (m,n,i)) with a structural element SE sized M SE × N SE × I SE = 3 × 3 × 3 pixels. It results in the matrix L K (m,n,i): Since only the outer limit is determined, its contour L P (1) (n,i) can be defined as follows: The results are, to a large extent, influenced by the brightness of the iris cross-section which, in turn, affects the overvaluation of the threshold value p r . As a result, the proper edge of the cornea is not properly recognized. For the analysed 13'400 2D images, the above method of edge detection works well in about 80% of cases. The value of 80% of cases was defined as the number of correctly detected corneal contours in all analyzed cases (13'400 2D images). The correctly detected corneal contour was described as not having any discontinuity points and the difference of the position of each contour point was below 10 pixels (maximum displacement of the contour line) between following series images.
The other approach involves using the well-known Canny edge detection method. The obtained result L C (m,n,i) of edge detection, for the arbitrarily selected parameters of the threshold (0.1) and standard deviation of the Gaussian filter (0.99), requires further correction. This adjustment only involves morphological closing, i.e.: Next, the binary image L W (m,n,i) is subjected to determination of the upper edge analogously to (3), (4). The resulting contour of the outer edge of the cornea L P (2) (n,i) is correct in about 90% of cases for the analysed 13'400 2D images. The definition of the correctly detected corneal contours were analogous to the previously described method.
Incorrectly detected cornea outer edges coincide for both methods, (L P (1) (n,i)) and L P (2) (n,i)), and for this reason they cannot be mutually complementary. Due to better edge detection in the case of the second method (Canny), it was used further (L P (n,i) for simplification). The results are shown in Figure 3.

Processing
The outer edge of the cornea L P (n,i), determined in image pre-processing, forms the basis for further analysis. Automatic designation of two important parameters, namely the correction of the cornea deformation curve and the rapidly changing cornea deformation in dynamic states when applying force, is described below.

Cornea deformation correction
Under the influence of a pressure impulse from the Corvis tonometer, the cornea is deformed and the entire eyeball moves in the eye socket. For proper evaluation and interpretation of the results, displacement of the entire eyeball should be separated from the deflection curve of the cornea. For this purpose, in the first step, the shape of the cornea visible in time t = 0 (for i = 1) is removed, i.e.: The resulting image L T (n,i) and the graph L T (N/2,i) are shown in Figure 4. The graph L T (N/2,i) consists of two parts of the response: corneal deflection L TR (N/2,i) and the eyeball reaction L TO (N/2,i) to an air-puff stimulus, in general: The waveform L TO (n,i) at the point N/2 defined with the formula (8) can be subtracted from L T (N/2,i) (formula (7)). The correctness of the algorithm was pre-verified by determining the correlation between the maximum amplitude of deviation calculated by the Corvis software and the same quantity determined from the equation (7) max(L TR ). The results obtained are shown in Figure 5 -the error calculated as the mean standard deviation is less than 1%. Next, image analysis results are shown in Figure 6. These results ( Figure 6) are extremely interesting from the point of view of biomechanics of the eye. The movement of the eyeball in response to a pressure impulse can be observed. As a result, it is possible to correct the waveform of the cornea deformation itself. Also the symmetry in the waveforms L TR (N/2,i) and L TO (N/2,i) is quite clearly visible - Figure 7. Separation of the eyeball displacement L TO (n,i) from the waveform of the cornea deformation L T (n,i) allows to measure another new parameter -rapidly changing cornea deformations visible during force application.

Rapidly changing cornea deformations during force application
Rapidly changing cornea deformations (with a frequency above 150 Hz) are directly visible in Figure 4. These deformations are not measured by the Corvis tonometer and they contain significant information about the biomechanical parameters of the eye (in particular the cornea). Their measurement is possible owing to the profiled image analysis algorithm.
The image L TR (n,i), resulting from the previous processing stage, undergoes an operation which involves the removal of a slowly changing deformation. For this purpose, the image resulting from morphological opening is subtracted from the image L TR (n,i), i.e.: L Q n; i ð Þ ¼ L TR n; i ð Þ− max SE2 min SE2 L TR n; i ð Þ ð Þ ð9Þ Figure 4 Sample graphs L T (n,i) and L T (N/2,i). The results obtained concern a relative cornea deformation resulting from an air-puff stimulus. The cross-section along the main axis N/2 is shown on the left a). This graph consists of two overlapping waveforms. One refers to the cornea deformation and the other one to the entire eyeball response to an air puff. The object of the analysis is, inter alia, separation of these two waveforms -b).
for a structural element SE2 sized N SE2 × I SE2 = 33 × 33 pixels. The size of the structural element SE2 is matched to the low-frequency range, lasting more than a few dozen milliseconds. The resulting image L Q (n,i) is shown in Figure 8a (i) were obtained from harmonics marked in Figure 8c with areas I and II. These are respectively: I-constant and slowly changing component (below 150 Hz) and II-component responsible for periodic Figure 5 Correlation between the maximum amplitude of the cornea deformation calculated with the Corvis software and the author's algorithm. These results confirm the correctness of the adopted algorithm of analysis and processing of image sequences and their reading from the format *.cst. On this basis, other parameters described in the paper were further determined and calculated.   rapidly changing (after the above 200 Hz) positions of the cornea on the left side relative to the axis of symmetry (L QL (i)). The analysis of the graphs F QR and F QL (Figure 8c) provides very valuable information regarding the frequency range of the cornea changes under the influence of an air puff. It should be noted that the frequency range for the analysed cases of 96 eyes is 150 to 500 Hz. A block diagram of the proposed algorithm is shown in Figure 9.

Results
On the basis of the proposed algorithm, a pattern L V (n,i) of the response of the cornea and the entire eyeball to an air puff was created: The obtained images of 96 eyes were compared to the pattern created for each eye separately: A sample result, the error L Δ (n,i), is shown in Figure 10. L Δ (n,i) is the basis for further analysis. Based on L Δ (n,i), the following new features are calculated: Figure 9 Block diagram of the proposed image analysis algorithm. These are, in order, image acquisition from the Corvis tonometer, median filtration, determination of the cornea outer edge, separation of the relative cornea deformation, separation of reactions, eyeball displacement and rapidly changing cornea deformation. As a result, new indicators describing the biomechanical parameters of the eye were obtained during force application with a pressure impulse. w(1) -the maximum amplitude of the eyeball deflection during force application: w(2) -the maximum amplitude in the spectra F QR and F QL for the frequency range from 150 to 500 Hz. With a fixed number of frames per second during registration with the Corvis tonometer (one frame every 23 μs), it is: w(3) -the duration of the rapidly changing (from 150 to 500 Hz) cornea deformations calculated as a time range with a threshold of 80% of the maximum amplitude - Figure 10.
The choice of features and their number is closely related to the previously described automatically calculated parameters such as the correction of the cornea deformation or the entire eyeball reaction to a stimulus. The obtained results, the values of w(1), w(2) and w (3), measured with the Corvis tonometer for the first 10 eyes are shown in Table 1. The correlation results of the features w(1) and w(2) with w(3) are shown in Figure 11. The results presented in Figure 11 were divided into 4 different classes. The adopted number of classes is the result of the analysis of the distribution of individual clusters forming sets that do not have common features. On this basis, a decision tree was created which is shown in Figure 12. In all cases, a non-parametrical algorithm CART (Classification and Regression Trees) creating binary trees is used as the method for their induction. An increase in the nodes purity has been used as the criterion for assessing the quality of CART divisions. The Gini index has been used as the measure of nodes impurity. Because of a small number of cases, the tree creation was not limited by a minimum number of vectors in a node. In its current form (Figure 12), there   is no need to trim the treethere is no concern that the decision tree will be overfit to the data. For the analysed number of cases, all classification results are correct. The number of the eyes classified as Class 1 is 6 ( Figure 11 -red), Class 3-2 ( Figure 11 blue), Class 4 -2 ( Figure 11 magenta), and the rest belongs to Class 2 ( Figure 11 green). An increase in the number of the analysed eyes (patients) may, in the future, give more reliable results and provide additional evidence as to the number of classes. The individual classes, depending on the range of feature values (w(1), w(2) and w (3)) which is characteristic for them, define various types of reactions of the entire eyeball to a stimulus and different amplitudes of rapidly changing deformations of the cornea itself. A sample graph of the dependence of the maximum displacement of the entire eyeball w(1) on the cornea deformation max(L TR ) for each class is shown in Figure 13.
Both the features as well as the division into classes can be further used not only as an additional diagnostic tool but they can also enrich the existing software of the Corvis tonometer.

Critical summary
The paper presents the possibility of using the profiled algorithm of image analysis, proposed by the authors, to measure additional parameters of the corneal deformation. This algorithm enables separation of the eyeball reaction from the corneal deflection during force application. In addition, measurements of rapidly changing cornea deformations (from 150 to 500 Hz) were proposed. These measurements allowed to propose a pattern individually for each eye. The pattern, in turn, helped to determine deviations individually for each eye. On the basis of these deviations, Figure 12 Decision tree created on the basis of the features w(1), w(2), w(3) acquired from the image and their division into four classes. Due to good separation of individual classes and a small number of nodes, the created decision tree requires no pruning (there is no concern that the decision tree will be overfit to the data).
the features w(1), w(2) and w(3) as well as the classifier (decision tree) allowing a division into four classes were calculated. The paper only presents a tool for automatic measurement of additional new parameters when using the Corvis tonometer. A detailed clinical examination based on this method will be presented in subsequent papers. In addition, the accuracy and repeatability of measurements for one patient and a phantom, the correlation of features (primarily w(1)) for a single patient and the right and left eyes and the effect of the patient position during the test on the result, both IOP as well as the other features, will be verified. The issue of analysis of symmetry, or its lack, in the dynamic response of the entire eyeball to an air-puff stimulus, is quite interesting and remains open. The analysis of corneal deflection parallel in other cross-sections is an open topic here. From a practical point of view, higher harmonics (above 2 kHz), that for the current sampling frequency (every 23 μs) are difficult to analyse, are also interesting. Moreover, the presented methods of image analysis and processing may be carried out in another way, in another dedicated algorithm, for example, with the use of random methods of analysis of layers in tomographic eye images [26,27], the methods of mathematical morphology [28] or others [29]. The final choice of the analysis method is highly dependent on the results obtained in practice, which is an important element given the large interindividual variability. Indeed, biomechanical parameters of the cornea, measured in dynamic states, still represent a new area of interdisciplinary research.