Synthetic OCT data for image processing performance testing

2011 18th IEEE International Conference on Image Processing SYNTHETIC OCT DATA FOR IMAGE PROCESSING PERFORMANCE TESTING P. Serranho C. Maduro, T. Santos, J. Cunha-Vaz, R. Bernardes IBILI, Faculty of Medicine, University of Coimbra, Azinhaga de Santa Comba, Celas, 3000-548 Coimbra, Portugal IBILI, Faculty of Medicine, University of Coimbra, and AIBILI, Coimbra, Portugal ABSTRACT The use of synthetic images is needed for testing the performance of image processing methods in order to establish a ground truth to test performance metrics. However, these synthetic images do not represent real applications. The aim of this paper is to build a mathematical model to obtain a synthetic noise-free image mimicking a real Optical Coherence Tomography (OCT) B-scan or volume from the human retina, in order to establish a ground truth for filtering performance metrics in this context. Moreover we also suggest a method to add speckle noise to this image based on the speckle noise of the given OCT volume. In this way we establish a replicable method to obtain a ground truth for image processing performance metrics that actually mimics a real case. Index Terms— Optical Coherence Tomography, Retina, Image Processing, Synthetic Image, Synthetic Noise. 1. INTRODUCTION Current image processing techniques (eg. despeckling filtering methods) with application to optical coherence tomography (OCT) rely on the respective qualitative evaluation of its results. Qualitative evaluation is usually subjective, since it depends strongly on the expert evaluating the results. On the other hand, quantitative approaches are reduced to using synthetic images which consists of an homogeneous background and a set of abstract geometric objects as concentric rings, squares or triangles with smooth or abrupt edges [1, 2] or are accepted by the image processing community as Lena or the photographer [3, 4, 5]. These solutions are therefore far from optimal, since these synthetic images do not represent the features of real medical image data. Approaches for optical testing using phantoms mimicking tissue have already been developed [6], however to the best of our knowledge no approach has been made in order to obtain a noise-free synthetic image of the retinal tissue as a ground truth to test the performance of image processing tools nor a noise model for OCT speckle. This work was partially funded by PTDC/SAU-BEB/103151/2008 and program COMPETE (FCOMP-01-0124 FEDER-010930). 978-1-4577-1302-6/11/$26.00 ©2011 IEEE 409 The purpose of the paper at hand is to establish a noisefree synthetic image that mimics real OCT data. In this way, we suggest a method to establish a ground truth for the computation of performance metrics in the context of OCT image processing tools. We stress that it only makes sense that a feature as fundamental as the mean square error is computed if a noise-free ground truth is known. Moreover we suggest a method to mimic the speckle noise present in the real OCT data (as in most optic based mechanisms [7, 8]), so that noise can be added to the created synthetic image. In this way we have a complete model for OCT image processing performance metrics, since we establish a realistic ground truth and noise pattern model. The computation of the synthetic image relies on features of the real OCT data. Each of the considered eye scans was processed in order to extract required parameters. For healthy subjects, only the segmentation of the inner limiting membrane (ILM), the outer segment layer (OSL) and the retinal pigment epithelium (RPE) is required. We note that for pathologic eyes, the segmentation of additional structures (cyst, epiretinal membrane, macular hole, etc...) to be preserved might also be needed. In each segmented region, OCT data is surface fitted using appropriate basis functions. As for the computation of the noise distribution, we do not aim at mimicking the physical process of speckle formation but to mimic the speckle noise grainy aspect of the given OCT scan. To achieve this, we study the intensity distribution in the vitreous. In this area, all the variability is due to noise, and therefore the model for the noise distribution is fitted accordingly, having in mind the grainy appearance of speckle noise. The paper is organized as follows. In section 2 we present in detail the method to obtain the synthetic retina and the synthetic noise. In section 3, we illustrate the results obtained, as well as results for the performed statistical tests to support the good performance of the proposed methods. Finally, in section 4 we summarize the results with some conclusions. 2011 18th IEEE International Conference on Image Processing 2. METHODS an approximation of the form Eye scans of 10 healthy volunteers were used following Cirrus OCT (Carl Zeiss Meditec, Dublin, CA, USA) scans using both the 200x200x1024 and the 512x128x1024 Macular Cube Protocols. 2.1. Synthetic image For each B-scan we consider the segmentation of the ILM, the interface of the upper OSL and the interface between the RPE and the choroid, which in our case was performed manually. In this way, the image is divided in four different areas as illustrated in figure 1 (red contour in the left image). Each area is then mapped to the unitary square [0, 1]2 such that the left, right, upper and lower boundary of each segmented region are mapped to the correspondent subsets within the unitary square of x = 0, x = 1, y = 0 and y = 1, respectively. According to the variability of the signal in each region, a set of basis functions is chosen in the unitary square in order to fit the given data. For the stability of the solution, we considered Tikhonov regularization. For instance, as the vitreous is supposed to be almost constant and this region should vary only in depth, we considered a linear approximation fvitreous (x, y) = a0 + a1 x + a2 y, x, y ∈ [0, 1], where x is the horizontal direction, y is the depth (vertical) one and a0 , a1 , a2 are the coefficients that best fit the given OCT data in a least squares sense. The following two regions (upper retina and OSL/RPE) present a much higher variability, so one needs a higher degree basis. In order to cope with high oscillations we considered a trigonometric series in each of the two regions considered of the form fretina (x, y) = K X N X ak,n cos(kπx) cos(nπy) K X N X bk,n cos(kπx) sin(nπy) k=0 n=1 + K X N X ck,n sin(kπx) cos(nπy) k=1 n=0 + K X N X 2 ak,n e−σx (x−xk ) −σy (2n −1)y , (2) k=0 n=0 for x, y ∈ [0, 1], with appropriate choices of K, N, σx and σy , and xk = k/K, for k = 0, 1, . . . , K, with ak,n being the coefficients that best fit the data for this region in a least square sense. (a) (b) (c) Fig. 1. (a) Original B-scan with given segmented areas; (b) obtained synthetic image; (c) Comparison over the retina of both A-scan profiles (green dashed line on (a)). 2.2. Noise model For each B-scan we consider a 100×100 pixel matrix Avitreous within the vitreous, in which we consider the distribution of intensities. We note that in the vitreous no oscillation in the signal intensities is expected, being therefore all the variability due to speckle noise. We consider Anoise = Avitreous − µAvitreous , k=0 n=0 + fchoroid (x, y) = K X N X dk,n sin(kπx) sin(nπy), (1) k=1 n=1 for x, y ∈ [0, 1] with appropriate choices of K and N , where ak,n , bk,n , ck,n , dk,n are again the coefficients that best fit to the given OCT data in each of these regions in a least squares sense. Finally, for the last region (choroid), it is clear that the signal is decreasing in depth y (due to absorption) and may vary in the horizontal coordinate x. Hence we chose 410 where µ holds for the mean. It is clear that µAnoise = 0 and that the standard deviations coincide (σAnoise = σAvitreous ). We stress that the statistical properties of the OCT speckle noise are represented in Anoise and one can easily determine its probability density function. As a starting point we create a random matrix B (0) with the same dimensions of the original image and the same density probability function of Anoise . In order to get the usual grainy appearance of OCT speckle noise in the synthetic noise matrix, we considered three different scales, namely the original matrix B (0) and the matrices obtained by taking the mean in m1 × n1 (denoted by B (1) ) and m2 × n2 blocks (denoted by B (2) ). We also considered an outlier matrix in each scale, given by  (k) (k)  Bi,j , Bi,j > 3σB (k) (k) Oi,j = , k = 0, 1, 2,  0 otherwise, 2011 18th IEEE International Conference on Image Processing where σB (k) is the standard deviation of the entries of B (k) . For each of these scale matrices we considered two low-pass Gaussian filters h1 and h2 with different standard variation σ1 and σ2 , respectively. Finally, we computed a weighted combination of these matrices     C (k) = α1 O(k) + α2 B (k) ∗ h1 + α3 B (k) ∗ h2 , (3) respectively. We also considered two low-pass Gaussian filters h1 and h2 such that the first is anisotropic with standard variation σ1 = [0.15, 1.5] and the second is isotropic with standard variation σ2 = [0.5, 0.5]. Finally, in equations (3) and (4) we considered for k = 0, 1, 2, where ∗ holds for the discrete convolution, consider In figure 3 we present the sum of the synthetic B-scan and the synthetic noise in order to compare with the original scan. The main visible difference between the synthetic and original data is at the beginning of the choroid, which is not a problem for the intended application, since the region of interest is the retina. M = β1 C (0) + β2 C (1) + β3 C (2) . α1 = 0.2, β1 = 0.6, (4) We note that the sum of both the triplets of the coefficients αi and βi , i = 1, 2, 3, should be equal to one. We consider M̃ = M − µM and make Asynthetic−noise = α2 = 0.2, α3 = 0.6, β2 = 0.3, β3 = 0.1. σAnoise M̃ , σM̃ in order to the synthetic noise matrix Asynthetic−noise to have the same global statistics as the real noise matrix Anoise from the OCT data. The results are illustrated in figure 2, showing that both the statistics and the visual appearance of the original and synthetic noise are similar. (a) (b) (c) (d) (e) (f) (g) (h) (i) (a) (b) (c) Fig. 3. (a) Original OCT B-scan; (b) Synthetic OCT B-scan with synthetic noise added; (c) Comparison of the probability density function of the intensities of the original (red dashed) and synthetic (blue) B-scan. Fig. 2. Original noise (top), synthetic noise (middle) and comparison of both probability density functions (bottom) with several signal-to-noise ratios, namely SNR=3 dB (left), SNR=7 dB (center) and SNR=10 dB (right). 3. RESULTS A total of 10 B-scans from 10 healthy eyes were processed resulting in the synthetic data representing the major characteristics of the respective real OCT scans. Having in mind the characteristic of the OCT signal in each region, we considered K = N = 8 for both regions referring to fretina in (1) and K = 8, N = 3, σx = 8, σy = 36 for fchoroid in (2). Moreover, for the construction of the matrices B (1) and B (2) we considered block of size 1 × 2 and 1 × 5, 411 We note that the standard variation for the difference between the original and the noise-free synthetic image is 14.29 ± 0.34 (mean ± std. dev.) while the noise level in the vitreous is 10.88 ± 0.27. This shows that the error of our approximation is of the same order of magnitude as the noise level. We have also compared the probability density functions of the intensities of both the original and noisy synthetic images, which are illustrated in figure 3(c). A Rank-sum test for the comparison of both probability functions for the sample of dimension 10 gave p-values 0.37±0.28, which shows that the distribution of intensities is similar in the synthetic and original images. A comparison within each segmented layer is also presented in figure 4. As expected, the main differences are presented within the retina, since we force the synthetic retina to be smoother than the original signal in order to have a smooth ground truth. In order to validate the synthetic noise model alone, we considered 100 × 100 pixel regions of the vitreous of the original and the correspondent synthetic data. A Rank-sum Test for the comparison of the two distribution gave p-values of 0.83 ± 0.12 in our 10-dimensional sample, showing that the speckle synthetic model mimics the real one with high accuracy. Moreover, the same 20 images were shown to 3 OCT 2011 18th IEEE International Conference on Image Processing given that the segmentation of additional structures (eg. cysts, epiretinal membranes, macular holes, etc...) is provided. The possibility of using less smooth approximation spaces in each segmented region will also be considered. (a) Acknowledgments (b) The authors would like to thank Dr. Melissa Horne and Carl Zeiss Meditec (Dublin, CA, USA) for their support on getting access to OCT data, and AIBILI Clinical Trial Center technicians for their support in managing data, working with patients and performing scans. (c) (d) 5. REFERENCES Fig. 4. Comparison of the probability density function of the intensities of the original (red dashed) and synthetic (blue) Bscan in each region ((a) Vitreous, (b) Upper Retina, (c) OSL and RPE, (d) choroid) for the example shown in figure 3. technicians in order for them to classify them into original or synthetic noise. Results demonstrate the difficulty in correct classification, with 52 misclassifications and only 8 correct classifications out of 60, that is, 86.6% and 13.3%, respectively. Additionally, out of the 30 classifications of synthetic noise images, 22 were classified as original (73.3%) and only 8 were classified as synthetic (26.7%). No original noise image was classified as such, illustrating that the synthetic noise images were very close to the original ones. 4. DISCUSSIONS AND CONCLUSIONS We propose a method to generate synthetic OCT data that mimics real one. The approach is two-folded. We obtain a noise-free synthetic B-scan that mimics the main features of a given original one with an appropriate segmentation, in order to establish a ground truth for processing methods. The method is designed so that for healthy eyes one only needs the segmentation of 4 areas (that is, 3 interfaces), which can be easily accomplished. Moreover, we showed the statistic validation of an OCT speckle noise model, based on the characteristics of the given real OCT data. If one adds the synthetic noise and the synthetic noise-free B-scan generated by our methods, one has an adequate synthetic image in the context of OCT processing methods. In this way, this process makes it possible to have a quantitative evaluation of the performance of any image processing procedure (eg. filtering) by providing adequate synthetic data as ground truth. These results have been used for testing the performance of an improved complex diffusion despeckling method [1] subsequently to its publication. We also intend to apply this method to pathologic eyes. This extension is straighforward, View publication stats 412 [1] R. Bernardes, C. Maduro, P. Serranho, A. Araújo, S. Barbeiro, and J. Cunha-Vaz, “Improved adaptive complex diffusion despeckling filter,” Opt. Express, vol. 18, no. 23, pp. 24048–24059, 2010. [2] H. Salinas and D. Fernández, “Comparison of pde-based nonlinear diffusion approaches for image enhancement and denoising in optical coherence tomography,” IEEE Trans Med Imaging, vol. 26 (6), pp. 761–771, 2007. [3] G. Gilboa, N. Sochen, and Y. Zeeni, “Image enhancement and denoising by complex diffusion processes,” IEEE Trans Pattern Anal Mach Intell, vol. 26, no. 8, pp. 1020 –1036, 2004. [4] P. Perona and J. Malik, “Scale-space and edge detection using anisotropic diffusion,” Pattern Analysis and Machine Intelligence, IEEE Transactions on, vol. 12, no. 7, pp. 629 –639, July 1990. [5] J. Rajan and M.R. Kaimal, “Image denoising using wavelet embedded anisotropic diffusion (wead),” IET Conference Publications, vol. 2006, no. CP522, pp. 589– 593, 2006. [6] D.C. Fernandez and H.M. Salinas, “A tissue phantom for investigating volume quantification on retinal images obtained with the stratus oct system,” in Engineering in Medicine and Biology Society, 2004. IEMBS ’04. 26th Annual International Conference of the IEEE, 2004, vol. 1, pp. 1225 –1228. [7] J. Dainty, A. Ennos, M. Françon, J. Goodman, T. McKechnie, G. Parry, and J. Goodman, “Statistical properties of laser speckle patterns,” in Laser Speckle and Related Phenomena, vol. 9 of Topics in Applied Physics, pp. 9–75. Springer Berlin / Heidelberg, 1975, 10.1007/BFb0111436. [8] J. W. Goodman, “Some fundamental properties of speckle,” J. Opt. Soc. Am., vol. 66, no. 11, pp. 1145– 1150, Nov 1976.