- Research
- Open access
- Published:
Automatic method of analysis of OCT images in assessing the severity degree of glaucoma and the visual field loss
BioMedical Engineering OnLine volume 13, Article number: 16 (2014)
Abstract
Introduction
In many practical aspects of ophthalmology, it is necessary to assess the severity degree of glaucoma in cases where, for various reasons, it is impossible to perform a visual field test - static perimetry. These are cases in which the visual field test result is not reliable, e.g. advanced AMD (Age-related Macular Degeneration). In these cases, there is a need to determine the severity of glaucoma, mainly on the basis of optic nerve head (ONH) and retinal nerve fibre layer (RNFL) structure. OCT is one of the diagnostic methods capable of analysing changes in both, ONH and RNFL in glaucoma.
Material and method
OCT images of the eye fundus of 55 patients (110 eyes) were obtained from the SOCT Copernicus (Optopol Tech. SA, Zawiercie, Poland). The authors proposed a new method for automatic determination of the RNFL (retinal nerve fibre layer) and other parameters using: mathematical morphology and profiled segmentation based on morphometric information of the eye fundus. A quantitative ratio of the quality of the optic disk and RNFL – BGA (biomorphological glaucoma advancement) was also proposed. The obtained results were compared with the results obtained from a static perimeter.
Results
Correlations between the known parameters of the optic disk as well as those suggested by the authors and the results obtained from static perimetry were calculated. The result of correlation with the static perimetry was 0.78 for the existing methods of image analysis and 0.86 for the proposed method. Practical usefulness of the proposed ratio BGA and the impact of the three most important features on the result were assessed. The following results of correlation for the three proposed classes were obtained: cup/disk diameter 0.84, disk diameter 0.97 and the RNFL 1.0. Thus, analysis of the supposed visual field result in the case of glaucoma is possible based only on OCT images of the eye fundus.
Conclusions
The calculations and analyses performed with the proposed algorithm and BGA ratio confirm that it is possible to calculate supposed mean defect (MD) of the visual field test based on OCT images of the eye fundus.
Introduction
Analysis of OCT images of the eye fundus is a technique known for many years which enables to acquire large amounts of information useful for medical diagnosis. Today, there are many known profiled algorithms for the analysis of both successive layers of the eye fundus [1–20] as well as the choroid layer [21]. Depending on the approach and practical implementation, these algorithms operate in the area of texture analysis [1, 2], algorithms based on segmentation with Canny method [4], SVM [5], mathematical morphology [6], active contour [7, 8], fuzzy algorithms [9], Markov model [10], robust segmentation [11], random contour analysis [12] and others [13–20]. Implementation of these algorithms most commonly occurs in the C, C# programming environment if they are used in commercially available applications, or, for example, in Matlab if they are at the stage of testing and modification [6]. The analysis time of a single OCT image of the fundus is several tens of milliseconds to a few seconds with the most computationally intensive applications [1, 17]. In practice, the biggest difficulty is large inter-individual variability and a possible high degree of retinal pathology. This translates into problems with the designation of successive layers of the eye fundus (RNFL–retinal nerve fibre layer, RPE–retinal pigment epithelium, and others) in the following cases (Figure 1):
-
The layer (u 1 , u 3 ) is not continuous - Figure 1a,d,
-
The layer (u 2 , u 3 ) merges with another layer (u 1 ) - Figure 1b,
-
One layer ends (u 1 or u 2 ) and another one starts (u 3 ) in the same place - Figure 1c.
The mentioned problems can be handled with in a limited way using artificial intelligence (knowledge gained from other images) or anthropometric evidence. In any case, image analysis and the proposed approach to solving such cases are difficult.
From a medical viewpoint, the RNFL thickness is essential in the case of the optic disk analysis. On the basis of its average thickness and further measurements, other parameters such as a cup/disk area ratio or a cup/disk vertical diameter are obtained. These parameters have an influence on the diagnosis of an ophthalmologist as well as the applied methods of treatment.
In the case of glaucoma, it is essential to test the visual field. Its measurement is usually carried out using a static perimeter. A reliable parameter obtained from perimetric data is the average depth of a defect. In many cases, for example, in a number of different additional diseases, it is impossible to obtain this type of measurement or the results obtained are not reliable. Then, it is necessary to search for other methods which will enable to determine the supposed visual field results. One such method is the analysis of OCT images of the eye fundus proposed by the authors. It allows not only to determine the severity degree of glaucoma but also to determine the supposed MD value of visual field test. It is further shown which selected parameters in OCT images correlate, to the greatest extent, with the results obtained from perimetry. Moreover, the new algorithm proposed by the authors and the quantitative ratio of the optic disk condition are described.
Material
OCT images of the eye fundus of 55 glaucoma patients (110 eyes) were obtained from the SOCT Copernicus (Optopol Tech. SA, Zawiercie, Poland). Patients ranged in age from 37 to 88 years. The obtained images (B-scans) had a resolution of M × N = 1000 × 600 pixels (where M – number of rows, N – number of columns) at a colour resolution of 8 bits/pixel. For each patient I = 99 B-scans were obtained allowing for full 3D reconstruction of the eye fundus image. At the same time, 110 results were acquired from the static perimeter OCTOPUS 300Series V 6.07c. All data were analysed in a source format. The presented analysis was carried out retrospectively on the basis of data obtained during typical medical examinations. The data were anonymised and the patients were examined in accordance with the Declaration of Helsinki.
Methods
The method of analysis and processing of OCT images proposed by the authors was preceded by the analysis of the results obtained from conventional methods. This analysis indicates features which should be taken into account when constructing the target image analysis algorithm.
Selection of the features of the image
As stated in the previous section, a visual field test was performed on 55 glaucoma patients. The features of the eye fundus were measured using a Copernicus tomograph. The results obtained in the software available for these devices are shown in Figure 2 and Figure 3. The features obtained on this basis are presented in Table 1. Examples of values for each feature for the first 10 patients are shown in Table 2. Each of these characteristics will be further, for simplicity, denoted by the symbol “w”. In this case, the features from w(1) to w(28) were acquired from the tomograph, and w(29) and w(30) from the perimeter. The parameters from w(1) to w(30) were further subjected to exploratory factor analysis in order to detect their cross correlation and to remove those features which are statistically insignificant. For this purpose, the following items were designated in succession: correlation of features w, eigenvalues of the correlation matrix, the number of classes using the Kaiser’s and Cattell’s criteria, rotation of the coordinate system, correlation of individual features with the newly created classes. The results (for the whole group of patients) are shown in Table 3 whereas the correlation chart of various features with three classes are shown in Figure 4. The analysis of Table 3 and Figure 4 suggests several conclusions useful for further analysis:
Class 1 significantly (above 0.5) correlates with the features w(16), w(13), w(6), w(5), w(8), w(3), w(15) and w(17) and is in a negative correlation with the features w(19), w(7) and w(4). It means that the cup/disk diameter and the cup diameter are a representative of the class 1. The ratio of the cup vertical/horizontal and the area and volume of the rim are in the inverse correlation.
Class 2 significantly correlates with the features w(22), w(20), w(24), w(25), w(23). These features are responsible for the RNFL parameters.
Class 3 correlates with the features w(10), w(2), w(12), w(11) and w(14) which are responsible for the parameters of the disk diameter in the vertical and horizontal axes.
Analysis of this division into classes, their correlation with the results obtained from perimetry and interpretation of individual features allow to specify the features of OCT images necessary for further analysis. The particularly important features of OCT image are: disk diameter and its relation to cup as well as the RNFL parameters in the area. This information is the basis for the design of a profiled image analysis algorithm, correlation analysis with the results of perimetry and for suggesting a new quantitative ratio of the optic disk condition.
Pre-processing
Pre-processing of images is related to automatic unpacking of the source file with the extension *.oct and thus acquiring images in the *.jpg format and calibration data. These data mainly include information recorded by the tomograph in the info.ini file, for example, the position of subsequent B-scans in space, the numbers of A-scans and B-scans, date of test, the patient data and others. The sequence of images L GRAY (m,n,i) obtained on this basis is further subjected to filtration with a median filter whose mask h 1 is 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. The filtered image L M (m,n,i) (Figure 5) is subjected to further preliminary transformations. These include determination of the RNFL.
Processing – RPE, ILM, RNFL, cup, disk
In order to determine the RNFL, the RPE layer (retinal pigment epithelium) is pre-determined automatically. The RPE layer is the most characteristic and the brightest layer in an OCT image, therefore, it is used to determine the target RNFL. Additionally, the RPE is used to determine the disk boundary. For this purpose, every n-th column of the image L M , each A-scan, is analysed. Analogously to the method described in [21, 22], the position of the maximum brightness for each column is pre-determined, i.e.:
where:
m,n – coordinates of the rows and columns of the matrix m∈(1,M) and n∈(1,N).
The formula (1) can be directly applied only when for all the analysed rows of each column there is only one maximum brightness value. As shown in [6], in practice it occurs in about 80% of the cases analysed at a resolution of 8 bits per pixel. On this basis, two essential elements are calculated: disk and ILM layer. The disk location, the image L disk , is here calculated as the area:
where: mean – mean value,
p m – threshold equal 2/3.
The value of 2/3 (p m ) was chosen arbitrarily and is associated with the maximum acceptable shift of the RPE layer in the area of the optic disk. The other element, namely the ILM layer area, was calculated as the area from the first pixel of each row to L RPE (n,i). It is very easy to perform binarization in this area using Otsu method [23], thus setting the binarization threshold p r . The resulting binary image L B (m,n,i) is:
Based on the binary image L B (m,n,i), the edge of the ILM layer was defined as:
The waveform of L ILM (n,i) and the RPE boundary, i.e. L RPE (n,i), are the basis for further analysis of the RNFL layer and cup area. The area between the layers ILM and RPE is further subjected to an operation which involves determination of the edge using two orthogonal Prewitt masks h h and h v sized 3 × 3 pixels. Filtration results in the image L G (m,n,i), i.e.:
where: L Gh (m,n,i) and L Gv (m,n,i) are the result of convolution of the image L M (m,n,i) with the masks h h and h v.
The RNFL is characterized by a contour visible in each column of the image L M (m,n,i) as a transition from a light area (below the ILM) to a dark one (the other layers). It is therefore the first edge below the ILM for L G (m,n,i) > 0, i.e. (analogously to (4), (5)):
The protection (m > L ILM (n,i))∧(m < L RPE (n,i)) is additional and prevents detection of the RNFL on the edge of the analysed object with the background. The resulting edge has a number of disadvantages. The biggest one is the lack of continuity caused by both occurring shadows as well as insufficient contrast in the place of its detection. For such cases, the method of active edge correction was suggested. This method involves adjustment of the position of the layer L RPE (n,i). For each point of the layer L RPE (n,i), the mean brightness value above and below L RPE (n,i) in the areas L u and L d is determined - Figure 6. The new corrected position of the layer L RPE *(n,i) is determined as the maximum difference between the average values in the areas L u (Δud,n,i) and L d (Δud,n,i). These areas are shifted in the column axis in the range Δud = ±med(|L ILM (n,i)-L RPE (n,i)|). The adjusted position of the RNFL is therefore calculated as:
The resulting new position of the layer L RNFL *(n,i) is the basis for further analysis. On the basis of the positions of the ILM and RPE layers, i.e.: L ILM (n,i), L RPE (n,i), the cup area L cup (n,i) is determined. This area is calculated as the intersection of the boundary of L disk (n,i) at the offset height equal to 150 μm with the layer L ILM (n,i).
As evident from the exploratory factor analysis presented in the previous section, apart from the RNFL thickness, also limitation of the analysis area - ROI (region of interest) is vital. The ROI ranges from 2∙r disk to 3∙r disk in an OCT image and in order to determine it, it is necessary to determine the position of coordinates of the centre (n d0 ,i d0 ), i.e.:
The value of r disk is the average radius calculated from the area of a circle:
The calculation of the centre of gravity is followed by the conversion of the ROI to the polar coordinate system, i.e.:
The function atan2 returns θ the same size as n and i containing the element-by-element, four-quadrant inverse tangent (arctangent) of the real parts of n and i. The obtained results, the ROI in the polar coordinate system L RNFL *(θ,ρ), are shown in Figure 7. The proposed algorithm is shown globally in a block diagram form in Figure 8. On the basis of the presented ROI, the following results (features) were obtained: the disk diameter and its relation to the cup and RNFL parameters. Their thorough analysis and correlation with the results from perimetry are presented in the next section.
Results
The obtained result L RNFL *(θ,ρ), the ROI in the polar coordinate system (Figure 7), is the basis for further analysis. Based on the presented exploratory factor analysis and the proposed method of image analysis and processing, the features representing different classes (Figure 4, Table 3), associated with L RNFL *(θ,ρ), were formulated in the following way:
w g (1) – the relative difference in the positions of centres of gravity calculated as:
where:
n d0 , i d0 – coordinates of the centre of gravity of the disk,
n c0 , i c0 – coordinates of the centre of gravity of the cup,
w g (2) – relative, inverse, minimum distance between the boundaries of the cup and disk i.e.:
where: r min – minimum distance between the boundaries of the cup and disk.
w g (3) – RNFL thickness deviation in the ROI (L RNFL ROI(θ)) for particular areas: temporal (TM), superior (SU), nasal (NA), interior (IN) from the standard thickness changes:
where: L RNFL WZ(θ) – the reference value of the RNFL thickness (the graph highlighted in green in Figure 2).
The mentioned features are components of the proposed biomorphological glaucoma advancement ratio (BGA) which is defined as:
After value substitution and simplification:
Figure 9 presents a list of various types of cup arrangements relative to the disk and the impact of their sizes on the values of the features w g (1), w g (2) and w g (3) as well as w BGA . In addition, for comparison, the values of the DDLS (Disk Damage Likelihood Scale) were also shown [24]. A change in the position of the centre of gravity of the cup relative to the disk is expressed by the feature w g (1). The closer the centres of gravity of the cup are to each other, the closer the feature w g (1) is to zero. If the centres of gravity of the cup and disk are far apart, the extreme position of the cup is on the edge of the disk, then the value of w g (1) is close to one. In turn, the feature w g (2) is a measure of centricity of the distribution of cup edges relative to the disk. When the edge of the cup is substantially spaced from the edge of the disk, the value of w g (2) is close to zero. When the edges of the cup and the disk overlap in any area, then the value of w g (2) is close to one. The situation is similar in the case of advanced glaucoma where r min is close to zero and w g (2) is close to one. A change in the value of w g (3) is related to a deviation with respect to the standard waveform of the RNFL thickness (Figure 2). The greater this difference is, the larger the value of w g (3) is, and vice versa. The features w g (1), w g (2) and w g (3) are summed up and thus determine the ratio w BGA. The ratio w BGA reaches values close to one if all three features w g indicate pathology. In the case of the right features w g , the ratio w BGA is equal to zero.
Based on the proposed ratio w BGA , its correlation with the results of perimetry was determined. In addition, it was compared with the results obtained only on the basis of known software (features w(1) to w(29)). The following results were obtained:
-
Correlation of classes 1, 2, 3 (obtained from features w(1) to w(29)) with the results obtained from perimetry is 0.78,
-
Correlation of the proposed ratio w BGA referred to the described method of image analysis and processing with the results of perimetry is 0.86 – Figure 10. The function describing the dependence of w(29) (denoting the MD value from perimetry) on w BGA takes the following form:
The difference in the correlation results stems from the profiled method of image analysis and processing and the measurement methodology of the described issue.
Summary
The paper presents a new methodology for image analysis and processing profiled for the optic disk analysis with a particular emphasis on its correlation with static perimetry. The calculations and analyses performed with the proposed algorithm confirm the possibility of prediction of the visual field result based on OCT images of the eye fundus. The analysis time of a single OCT image of the eye fundus with the use of the proposed algorithm on the computer Core i5 CPU M460 @ 2.5GHz 4GB RAM does not exceed 0.2 s. Additionally, a new ratio reliable with respect to the quantitative assessment of the optic disk condition was proposed. On its basis, a comparison of the results of the correlation between the known methods and the methodology described in this paper was made.
The present paper does not exhaust this interesting subject. The proposed method of image analysis and processing is one of the possible approaches to solving this type of task. Certainly, other methods of image analysis and processing can be used here [25, 26] such as texture analysis [27, 28] or other methods applied to solving similar issues [29–31].
In further studies, the authors intend to test practical usefulness of the proposed algorithm and the ratio w BGA . The results obtained for large inter-individual variability (diversity of shapes and topology of the optic disk), repeatability of the data and the impact of technological diversity and device parameters (tomograph and perimeter) on the result are particularly vital here.
Abbreviations
- AMD:
-
Age-related macular degeneration
- ROI:
-
Region of interest
- RNFL:
-
Retinal nerve fiber layer
- TE:
-
Temporal
- SU:
-
Superior
- NA:
-
Nasal
- IN:
-
Inferior
- BGA:
-
Biomorphological glaucoma advancement
- MD:
-
Mean defect
- LV:
-
Loss variance.
References
Quellec G, Lee K, Dolejsi M, Garvin MK, Abràmoff MD, Sonka M: Three-dimensional analysis of retinal layer texture: identification of fluid-filled regions in SD-OCT of the macula. IEEE Trans Med Imaging 2010,29(6):1321–1330.
Kajić V, Povazay B, Hermann B, Hofer B, Marshall D, Rosin PL, Drexler W: Robust segmentation of intraretinal layers in the normal human fovea using a novel statistical model based on texture and shape analysis. Opt Express 2010,18(14):14730–14744. 10.1364/OE.18.014730
Park SY, Kim SM, Song YM, Sung J, Ham DI: Retinal thickness and volume measured with enhanced depth imaging optical coherence tomography. Am J Ophthalmol 2013,156(3):557–566. 10.1016/j.ajo.2013.04.027
Koprowski R, Wrobel Z: Identification of layers in a OCT image of an eye based on the canny edge detection. Inf Technol Biomed Adv Intell Soft Comput 2008, 47: 232–239. 10.1007/978-3-540-68168-7_26
Wang S, Zhu WY, Liang ZP: ShapeDeformation: SVM regression and application to medical image segmentation. Eight Int Conf Comput Vis 2001, 2: 209–216.
Koprowski R, Wróbel Z: Image Processing in Optical Coherence Tomography: Using Matlab. Katowice, Poland: University of Silesia; 2011. http://www.ncbi.nlm.nih.gov/books/NBK97169/
Yazdanpanah A, Hamarneh G, Smith B, Sarunic M: Intra-retinal layer segmentation in optical coherence tomography using an active contour approach. Med Image Comput Comput Assist Interv 2009,12(2):649–656.
Yezzi A, Tsai A, Willsky A Proceedings of the seventh IEEE international conference on computer vision. In A Statistical Approach to Snakes for Bimodal and Trimodal Imagery. Washington, DC: IEEE Computer Society; 1999:898–903.
Mayer MA, Tornow RP, Bock R, Hornegger J, Kruse FE: Automatic nerve fiber layer segmentation and geometry correction on spectral domain OCT images using fuzzy C-means clustering. Invest Ophthalmol Vis Sci 1880, 2008: 49.
Koozekanani D, Boyer KL, Roberts C: Retinal thickness measurements in optical coherence tomography using a Markov boundary model. Med Imaging IEEE Trans 2001,20(9):900–916. 10.1109/42.952728
Gregori G, Knighton RW: A robust algorithm for retinal thickness measurements using optical coherence tomography (stratus OCT). Invest Ophthalmol Vis Sci 2004, 45: 3007.
Koprowski R, Wróbel Z: Layers recognition in OCT eye image based on random contour analysis. Computer recognition systems 3. Adv Intell Soft Comput 2009, 57: 471–478. 10.1007/978-3-540-93905-4_56
Shahidi M, Wang Z, Zelkha R: Quantitative thickness measurement of retinal layers imaged by optical coherence tomography. Am J Ophthalmol 2005,139(6):1056–1061. 10.1016/j.ajo.2005.01.012
Mujat M, Chan R, Cense B, Park B, Joo C, Akkin T, Chen T, Boer J: Retinal nerve fiber layer thickness map determined from optical coherence tomography images. Opt Express 2005,13(23):9480–9491. 10.1364/OPEX.13.009480
Mishra A, Wong A, Bizheva K, Clausi DA: Intra-retinal layer segmentation in optical coherence tomography images. Opt Express 2009,17(26):23719–23728. 10.1364/OE.17.023719
Lu S, Cheung C, Liu J, Lim S, Leung C, Wong T: Automated layer segmentation of optical coherence tomography images. IEEE Trans Biomed Eng 2010,57(10):2605–2608.
DeBuc C, Somfai GM: Early detection of retinal thickness changes in diabetes using optical coherence tomography. Med Sci Monit 2010,16(3):MT15-MT21.
Haeker M, Abràmoff M, Kardon R, Sonka M: Segmentation of the surfaces of the retinal layer from OCT images. Med Image Comput Comput Assist Interv 2006,9(1):800–807.
Ishikawa I, Stein DM, Wollstein G, Beaton S, Fujimoto JG, Schuman JS: Macular segmentation with optical coherence tomography. Invest Ophthalmol Vis Sci 2005,46(6):2012–2017. 10.1167/iovs.04-0335
Hu Z, Abramoff MD, Kwon YH, Lee K, Garvin M: Automated segmentation of neural canal opening and optic Cup in 3-D spectral optical coherence tomography volumes of the optic nerve head. Invest Ophthalmol Vis Sci 2010,51(11):5708–5717. 10.1167/iovs.09-4838
Koprowski R, Teper S, Wrobel Z, Wylegala E: Automatic analysis of selected choroidal diseases in OCT images of the eye fundus. Biomed Eng Online 2013, 12: 117. 10.1186/1475-925X-12-117
Koprowski R, Teper S, Weglarz B, Wylęgała E, Krejca M, Wróbel Z: Fully automatic algorithm for the analysis of vessels in the angiographic image of the eye fundus. Biomed Eng Online 2012, 11: 35. 10.1186/1475-925X-11-35
Otsu N: A threshold selection method from gray-level histograms. IEEE Trans Sys Man Cyber 1979,9(1):62–66.
Wei H, Zangalli C, Spaeth GL: Valid estimation of the appearance of the optic nerve: the case against using Cup/disk ratios. J Curr Glaucoma Pract 2010,4(3):133–136.
Sonka M, Michael Fitzpatrick J Handbook of Medical Imaging. In Medical Image Processing and Analysis. Belligham: SPIE; 2000.
Tadeusiewicz R, Ogiela MR: Automatic understanding of medical images new achievements in syntactic analysis of selected medical images. Biocybern Biomed Eng 2002,22(4):17–29.
Korzynska A, Zychowicz M: A method of estimation of the cell doubling time on basis of the cell culture monitoring data. Biocybern Biomed Eng 2008,28(4):75–82.
Korzynska A, Iwanowski M: Multistage morphological segmentation of bright-field and fluorescent microscopy images. Opt-Electron Rev 2012,20(2):87–99.
Siedlecki D, Nowak J, Zajac M: Placement of a crystalline lens and intraocular lens: retinal image quality. J Biomed Opt 2006,11(5):054012. 10.1117/1.2358959
Porwik P: Efficient spectral method of identification of linear Boolean function. Control Cybern 2004,33(4):663–678.
Korus P, Dziech A: Efficient method for content reconstruction with self-embedding. IEEE Trans Image Process 2013,22(3):1134–1147.
Acknowledgements
No outside funding was received for this study.
Author information
Authors and Affiliations
Corresponding author
Additional information
Competing interests
The authors declare that they have no competing interests.
Authors’ contributions
RK suggested the algorithm for image analysis and processing, implemented it and analyzed the images. MR, ZW performed the acquisition of the tomography images and consulted the obtained results. All authors have read and approved the final manuscript.
Authors’ original submitted files for images
Below are the links to the authors’ original submitted files for images.
Rights and permissions
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Koprowski, R., Rzendkowski, M. & Wróbel, Z. Automatic method of analysis of OCT images in assessing the severity degree of glaucoma and the visual field loss. BioMed Eng OnLine 13, 16 (2014). https://doi.org/10.1186/1475-925X-13-16
Received:
Accepted:
Published:
DOI: https://doi.org/10.1186/1475-925X-13-16