Retinal optical coherence tomography image enhancement via shrinkage denoising using double-density dual-tree complex wavelet transform

Shahab Chitchian; Markus A. Mayer; Adam Boretsky; Frederik J. van Kuijk; Massoud Motamedi

doi:10.1117/1.JBO.17.11.116009

1 November 2012 Retinal optical coherence tomography image enhancement via shrinkage denoising using double-density dual-tree complex wavelet transform

Shahab Chitchian; Markus A. Mayer; Adam Boretsky; Frederik J. van Kuijk; Massoud Motamedi

Author Affiliations +

J. of Biomedical Optics, 17(11), 116009 (2012). doi:10.1117/1.JBO.17.11.116009

Abstract

Image enhancement of retinal structures, in optical coherence tomography (OCT) scans through denoising, has the potential to aid in the diagnosis of several eye diseases. In this paper, a locally adaptive denoising algorithm using double-density dual-tree complex wavelet transform, a combination of the double-density wavelet transform and the dual-tree complex wavelet transform, is applied to reduce speckle noise in OCT images of the retina. The algorithm overcomes the limitations of commonly used multiple frame averaging technique, namely the limited number of frames that can be recorded due to eye movements, by providing a comparable image quality in significantly less acquisition time equal to an order of magnitude less time compared to the averaging method. In addition, improvements of image quality metrics and 5 dB increase in the signal-to-noise ratio are attained.

Chitchian, Mayer, Boretsky, van Kuijk, and Motamedi: Retinal optical coherence tomography image enhancement via shrinkage denoising using double-density dual-tree complex wavelet transform

1. Introduction

The unique optical properties and physical characteristics of the human eye make it possible to noninvasively visualize the structure and morphology of the retina and diagnose retinal pathologies such as diabetic retinopathy, macular degeneration, and glaucoma using various optical imaging modalities.¹

Applications of optical coherence tomography (OCT), as a noninvasive optical imaging technique to perform high-resolution imaging of retinal microstructure, have well been demonstrated and commercialized in ophthalmology.² Automatic analysis of OCT images, such as retinal layer segmentation, provides objective quantification of retinal layers from the nerve fiber layer to the pigment epithelium layer. Measurements derived from such image analysis methods may improve the clinical diagnosis during glaucoma progression and age related macular degeneration.

OCT images are subject to various distortions that introduce artifacts in the OCT cross-sections such as floaters and saccades. Speckle noise, as a dominant source of artifacts, also occurs in the OCT image because particles composed in the underlying tissue structures are smaller than the coherence length of the light source.³ In other words, there is a limited spatial-frequency bandwidth of OCT interference signals.

Recently, wavelet techniques have been employed successfully in speckle noise reduction for OCT images.⁴5.6.^–⁷ Adler et al., Wagner et al., and Mayer et al., applied the discrete wavelet transform which provides the most compact representation, however, it has two main limitations which include the fact that it is not shift-invariant and not oriented in two dimensions.⁴^,⁶7.^–⁸ The dual-tree complex wavelet transform overcomes these limitations. The complex wavelet transform was used for wavelet denoising in OCT images where Chitchian et al. applied a denoising algorithm to reduce speckle noise using the dual-tree complex wavelet transform.⁷^,⁹

Advanced retinal tracking systems and registration algorithms are generally used to average multiple image frames to improve overall image quality and presentation. However, these methods have practical limitations due to the inability of patients to maintain fixation during examinations. By averaging multiple frames, these distortions may be suppressed, but information that could improve diagnosis or allow for more sophisticated image analysis is also lost. Using the same clinical system and wavelet denoising of eight frames, Wagner et al. and Mayer et al. recently reached an SNR which is comparable to an averaging of 29 frames.⁶^,⁷

More recently, Leigh et al. demonstrated that wavelet noise estimation technique is the overall preferred noise estimation technique for most image and video processing applications.¹⁰ Therefore, building on the earlier results, a denoising algorithm using double-density dual-tree complex wavelet transform is applied to a two-frame OCT image, the original image from OCT clinical system without preprocessing, to improve the image quality and overcome the limitations of commonly used multiple-frame averaging technique.⁹ Image quality metrics improvements and SNR increase that are better than those of 50-frame averaged images, as processed by the system, are achieved.

2. Wavelet Shrinkage Denoising

The wavelet shrinkage denoising algorithm requires the following four-step procedure,⁹

(1)

Y = W (X) λ = d (Y) Z = D (Y, λ) S = W^{- 1} (Z),

where

W (\cdot)

relates to the wavelet transform,

d (\cdot)

selects a data-adaptive threshold,

D (\cdot, λ)

denotes the denoising with threshold

λ

, and

W^{- 1} (\cdot)

performs the inverse wavelet transform. More detailed descriptions of each step are presented below.

2.1.

Double-Density Dual-Tree Complex Wavelet Transform

The dual-tree complex wavelet transform (CDWT) calculates the complex transform of a signal using two separate discrete wavelet transform (DWT) decompositions. By utilizing a Hilbert transform pair, it is possible for one DWT to produce the real coefficients and the other the imaginary coefficients.⁹

The double-density wavelet transform (DD-DWT) is defined by recursively applying the three-channel analysis filter bank to the signal. Selesnick previously presented approximate Hilbert transform pairs of wavelet frames that have the advantages of both types of wavelet frames described above.¹¹ The DWT based on these wavelet frames is called the double-density dual-tree complex wavelet transform (DD-CDWT). The Double-Density Wavelet Software by Selesnick et al. was used for implementation of the wavelet transform.¹²

2.2.

Shrinkage Denoising

Operator $d (\cdot)$ , in Eq. (1), selects a data-adaptive threshold and $D (\cdot, λ)$ denotes the denoising operator with threshold $λ$ . Bivariate shrinkage with local variance estimation algorithm is applied for shrinkage denoising.¹³

After estimating the signal components of the noisy coefficients in the wavelet domain, the inverse wavelet transform, $W^{- 1} (\cdot)$ , is taken to reconstruct the noise-free image.

3. Results

The top of Fig. 1 shows the scanning laser ophthalmoscope (SLO) image captured by the Spectralis OCT (Heidelberg Engineering, Heidelberg, Germany) during the eye scanning process. The B-Scans positions are marked, as well as an example of OCT B-Scans, is also shown. The scans were centered on the fovea in order to image the macular region of the retina. The B-Scans consisted of 768 A-Scans. Each A-Scan consists of 496 pixels. The B-Scans were cropped into $288 \times 288 pixels$ for analysis. The axial resolution of the Spectralis is 7 µm in tissue. The pixel length in the lateral direction is 5.73 µm.

Fig. 1

SLO image captured by the Spectralis OCT and a corresponding OCT B-Scan (top). The selected image: (a) original; (b) DD-CDWT denoised. Eight ROIs are shown in the original image. Combined application of the proposed denoising algorithm and the averaging method available in clinical systems: (c) 10-frame averaged image; and (d) DD-CDWT denoised of the averaged image.

For the present work, 7 B-Scans were acquired from a subject. 7 B-Scans are two-frame averaged images, the original images from the Spectralis OCT system without preprocessing.

The algorithm was executed on a Intel Core i5-2400, 3.10 GHz desktop personal computer. The total time for the denoising process was less than half a second using Matlab 7.14 (MathWorks, Natick, MA).

The unprocessed two-frame OCT image of the retina is shown in Fig. 1(a). Figure 1(b) shows the image after DD-CDWT denoising.

The performance of a new wavelet denoising procedure can be compared to the averaging approach which is now being used in latest clinical OCT systems including the Spectralis to remove speckle noise and improve SNR. Because of the eye movement and the time the eye can be kept open, especially in elder patients, there are limitations for clinicians not to average more than 20 frames. For quantitative comparison, the 50-frame averaged images were used because they have more improved results compared to the 20-frame averaged images. All averaged images were acquired using retinal tracker capability of the Spectralis. Also, the images were from a young and healthy subject with perfect fixation. It is clear that, in any cases of abnormalities, the results of denoising are more effevtive.

Image quality metrics were used to assess the performance of the denoising technique by measuring the contrast-to-noise ratio (CNR), the equivalent number of looks (ENL), the correlation $ρ$ , the structure similarity (SSIM), and full-width half-maximum (FWHM) over regions of interest (ROIs). Eight ROIs were defined in the original images, as shown in Fig. 1(a). The ROIs were chosen to cover different retinal layers and edges, and were used for different measurements. The metrics were calculated as the average over the ROIs used. In addition, the global SNR is calculated.

The CNR measures the contrast between image features [red and blue ROIs in Fig. 1(a)] and an area of background noise, while the ENL measures smoothness in areas that should have a homogeneous appearance [blue ROIs in Fig. 1(a)].⁹

Performance of edge preservation or sharpness is evaluated based on correlation. A parameter, $ρ$ , is calculated as:¹⁴

(2)

ρ = \frac{Γ (x - μ_{x}, \hat{x} - μ_{\hat{x}})}{\sqrt{Γ (x - μ_{x}, x - μ_{x}) \cdot Γ (\hat{x} - μ_{\hat{x}}, \hat{x} - μ_{\hat{x}})}},

where

x

is the unprocessed image and

\hat{x}

is the enhanced image in the location of the ROI, and

(3)

Γ (x_{1}, x_{2}) = \sum_{(i, j) \in ROI} x_{1} (i, j) \cdot x_{2} (i, j) .

Furthermore, SSIM index was also used to compare the structure of the unprocessed image, $x$ , and the enhanced image, $\hat{x}$ , in the ROI,

(4)

SSIM (x, \hat{x}) = \frac{(2 μ_{x} μ_{\hat{x}} + C_{1}) (2 σ_{x \hat{x}} + C_{2})}{(μ_{x}^{2} + μ_{\hat{x}}^{2} + C_{1}) (σ_{x}^{2} + σ_{\hat{x}}^{2} + C_{2})},

where

C_{1}

and

C_{2}

are constants used to avoid instability. SSIM compares local patterns of pixel intensities that have been normalized for luminance and contrast.¹⁵

ρ

and SSIM are measured in the red and blue ROIs in Fig. 1(a).

The sharpness of the edges is measured by FWHM at different edge ROIs in the image [green ROIs in Fig. 1(a)]. The examined edge in the ROI is deformed to a straight line. The values in the direction of the edge are summed up to reduce the influence of the noise. A logistic sigmoid function is fitted to the resulting values by using a nonlinear regression method. The FWHM is the width of the derivation of this fitted sigmoid function at its half maximum. The smaller this value is, the sharper the edge. We measured the FWHM at two sharp edges as well as one edge with low contrast. The sigmoid function is defined as follows

(5)

Θ (x) = q_{1} + \frac{q_{2}}{1 + \exp (- \frac{x + q_{3}}{q_{4}})},

where

q_{1}

to

q_{4}

are parameters that are optimized by the nonlinear regression using the Matlab nlinfit method.

The FWHM measures the width of $Θ^{'} (x)$ at the half of the maximum ( $Θ^{'} (- q_{3}) / 2$ ),

(6)

FWHM = | q_{4} \ln (\frac{- 2 d + 1 - \sqrt{1 - 4 d}}{- 2 d + 1 + \sqrt{1 - 4 d}}) |,

with

(7)

d = \frac{Θ^{'} (- q_{3}) q_{4}}{q_{2}} .

Table 1 shows image quality values for the original, 50-frame averaged, CDWT, and DD-CDWT wavelet denoised images. The values of CNR, ENL, SNR, $ρ$ , and SSIM for seven sample images show significant improvement. A paired $t$ -test was performed for comparison of CDWT and DD-CDWT wavelet denoising algorithms with statistical significance given by values of $P < 0.05$ . Comparison of CNR, ENL, SNR, $ρ$ , and SSIM results in $P$ values of 0.0083, 0.1633, 0.3788, 0.0001, and 0.0001, respectively. The values for ENL and SNR are considered not to be statistically significant. However, the values for CNR, $ρ$ , and SSIM are considered to be statistically significant.

Table 1

Image quality values for comparison of original, 50-frame averaged, CDWT, and DD-CDWT wavelet denoised images of normal retina and for comparison of original, CDWT, and DD-CDWT wavelet denoised images of a 20-frame averaged AMD case.

n=7	Image	CNR (dB)	ENL	SNR (dB)	FWHM	ρ (%)	SSIM (%)	Time (s)
	Original	6.4	352	24.85	3.46	100	100	0
Mean	Averaged (50-frame)	7.8	1975	25.96	3.61	53.70	33.45	4.4
	Denoised (CDWT)	12.2	555	29.67	3.72	96.07	79.86	0.4
	Denoised (DD-CDWT)	10.6	545	29.86	3.80	97.71	86.71	0.4
P value	Paired t-test (DD-CDWT vs. 50-frame averaged)	.0053	.0012	.0303	.7375	.0001	.0001	.0001
P value	Paired t-test (DD-CDWT vs. CDWT)	.0083	.1633	.3788	.6614	.0001	.0001	.0001
AMD	Averaged (20-frame)	8.4	1059	26.16	3.77	100	100	6.0
(n=1)	Denoised (CDWT) of Averaged (20-frame)	11.4	1279	28.69	3.82	98.05	88.79	6.4a
	Denoised (DD-CDWT) of Averaged (20-frame)	14.5	1410	31.07	3.85	97.54	87.88	6.4a

a

The averaging time was added to the denoising time.

DD-CDWT has better functioning to retain both the structure and sharpness of the unprocessed image in comparison to CDWT, the differences are extremely statistically significant. CDWT has superior CNR while ENL and SNR statistically remain the same for both approaches.

The averaging approach is now being used in latest clinical OCT systems including the Spectralis to remove speckle noise and improve SNR. Comparison of the averaged and DD-CDWT denoised images results in $P$ values for CNR, ENL, SNR, $ρ$ , and SSIM of 0.0053, 0.0012, 0.0303, 0.0001, and 0.0001, respectively. All values are considered to be statistically significant. Comparison of FWHMs shows sharpness in the edges are not reduced significantly using these methods.

While the averaging technique has a better ENL, almost four times the value for DD-CDWT, its performance is reduced compared to DD-CDWT denoising algorithm in terms of CNR and SNR. The main disadvantage is deterioration of the structure and sharpness of the unprocessed image, 33.45 percent and 53.70 percent, respectively. However, denoising, using DD-CDWT technique, preserves the structure and sharpness by factors of 86.71 percent and 97.71 percent, respectively. In other words, the enhancement approaches including averaging and the proposed technique, can have a false inducing effect, an artifact in the original image. By measuring sharpness and SSIM, these artifacts were quantified.

Feature and sharpness preservation will play an important role in accurately measuring cellular disruption to determine retinal pathology. In summary, the proposed denoising algorithm is a highly promising preprocessing approach which yields enhanced retinal OCT images for further computational analysis such as segmentation of retinal layers.

Last, it should be noted that enhancement based on averaged frames takes more time than the proposed wavelet shrinkage denoising algorithms and this is illustrated in Table 1. This imposes limitations for clinicians not to average more than 20-frame because of the eye movement and the time the eye can be kept open, especially in elder patients. An important advantage of the proposed denoising technique is not only to overcome this limitation with less than half a second enhancing time which is equal to an order of magnitude less time compared to the averaging method, but also to provide better results in terms of image quality metrics.

Using the same OCT system and similar type of images, Wagner et al. and Mayer et al. could improve SNR to the level of 29-frame averaged performance by applying wavelet denoising to the 8-frame averaged images.⁶^,⁷ However, the introduced approach went steps further to enhance SNR better than the 50-frame averaged level by using only two-frame images which are the original images from the Spectralis OCT system without preprocessing.

Finally, in Fig. 1(c) and 1(d), we show the combined application of the proposed denoising algorithm and the averaging method available in clinical systems. Figure 1(c) shows the commonly used ten-frame averaged image by clinicians from the Spectralis. The denoising result of Fig. 1(c) using DD-CDWT is shown in Fig. 1(d). Comparison of Fig. 1(c) and 1(d) demonstrates that the performance of the averaging method can also be improved by applying the proposed denoising technique in series.

As a demonstration of clinical applicability, Fig. 2(a) and 2(b) gives an example of an intermediate age-related macular degeneration (AMD) case before and after denoising, respectively. This case represents significant disruption in retinal structure, particularly around the retinal pigment epithelium and photoreceptors. The image before denoising, Fig. 2(a), is 20-frame averaged image acquired in the clinic. Table 1 quantifies the image quality metrics before and after denoising, demonstrating that the denoising algorithm substantially improves the image quality in terms of CNR, ENL, and SNR while preserving structure and sharpness. The proposed DD-CDWT denoising algorithm could attain 5 dB increase in SNR along with 6 dB increase in CNR with less than 3 percent and 13 percent reduction in $ρ$ and SSIM, respectively.

Fig. 2

Denoising of an intermediate AMD retinal image: (a) before (20-frame averaged); (b) after DD-CDWT denoising.

4. Conclusion

The proposed denosing algorithm provides significant improvements in image quality metrics while preserving subtle features of the retinal layers. This technique may be implemented in ophthalmic OCT systems currently being used to diagnose and monitor patients with diabetes, glaucoma, age related macular degeneration, and other retinal diseases.

Acknowledgments

The authors thank Dr. Ivan Selesnick of Polytechnic University, Brooklyn, NY for giving us advice to perform the software and providing Double-Density Wavelet Software. This research was supported by Seymour Fisher Scholarship from Ophthalmology Department at UTMB Health and funding for Research to Prevent Blindness (RPB) to UTMB Health and the University of Minnesota. The authors gratefully acknowledge funding of the Erlangen Graduate School in Advanced Optical Technologies (SAOT) by the German Research Foundation (DFG) in the framework of the German excellence initiative.

References

1.

M. AbramoffM. GarvinM. Sonka, “Retinal imaging and image analysis,” IEEE Rev. Biomed. Eng. 3, 169–208 (2010).1937-3333http://dx.doi.org/10.1109/RBME.2010.2084567Google Scholar

2.

D. Huanget al., “Optical coherence tomography,” Science 254(5035), 1178–1181 (1991).SCIEAS0036-8075http://dx.doi.org/10.1126/science.1957169Google Scholar

3.

S. XiangL. ZhouJ. Schmitt, “Speckle noise reduction for optical coherence tomography,” Proc. SPIE 3196, 79–88 (1998).PSISDG0277-786Xhttp://dx.doi.org/10.1117/12.297921Google Scholar

4.

D. AdlerT. KoJ. Fujimoto, “Speckle reduction in optical coherence tomography images by use of a spatially adaptive wavelet filter,” Opt. Lett. 29(24), 2878–2880 (2004).OPLEDP0146-9592http://dx.doi.org/10.1364/OL.29.002878Google Scholar

5.

A. Pizuricaet al., “Multiresolution denoising for optical coherence tomography: a review and evaluation,” Curr. Med. Imag. Rev. 4(4), 270–284 (2008).CMIRCV1573-4056http://dx.doi.org/10.2174/157340508786404044Google Scholar

6.

M. Wagneret al., “Wavelet based approach to multiple-frame denoising of OCT images,” in 2009 Russian Bavarian Conference on Bio-Medical Engineering (RBC), pp. 67–69, Russian Bavarian Conference (RBC), Munchen, Germany (2009).Google Scholar

7.

M. Mayeret al., “Wavelet denoising of multiframe optical coherence tomography data,” Biomed. Opt. Express 3(3), 572–589 (2012).BOEICL2156-7085http://dx.doi.org/10.1364/BOE.3.000572Google Scholar

8.

N. Kingsbury, “Complex wavelets for shift invariant analysis and filtering of signals,” Appl. Comput. Harmonic Anal. 10(3), 234–253 (2001).ACOHE91063-5203http://dx.doi.org/10.1006/acha.2000.0343Google Scholar

9.

S. ChitchianM. FiddyN. Fried, “Denoising during optical coherence tomography of the prostate nerves via wavelet shrinkage using dual-tree complex wavelet transform,” J. Biomed. Opt. 14(1), 014031 (2009).JBOPFO1083-3668http://dx.doi.org/10.1117/1.3081543Google Scholar

10.

A. Leighet al., “Comprehensive analysis on the effects of noise estimation strategies on image noise artifact suppression performance,” in Proc. IEEE Sym. Multimedia, pp. 97–104, IEEE, Dana Point, CA (2011).Google Scholar

11.

I. Selesnick, “The double-density dual-tree DWT,” IEEE Trans. Signal Process. 52(5), 1304–1314 (2004).ITPRED1053-587Xhttp://dx.doi.org/10.1109/TSP.2004.826174Google Scholar

12.

C. WagnerI. Selesnick, “Double-density wavelet software,” http://taco.poly.edu/selesi/DoubleSoftware/.Google Scholar

13.

L. SendurI. Selesnick, “Bivariate shrinkage functions for wavelet-based denoising exploiting interscale dependency,” IEEE Trans. Signal Process. 50(11), 2744–2756 (2002).ITPRED1053-587Xhttp://dx.doi.org/10.1109/TSP.2002.804091Google Scholar

14.

F. Sattaret al., “Image enhancement based on a nonlinear multiscale method,” IEEE Trans. Image Process. 6(6), 888–895 (1997).IIPRE41057-7149http://dx.doi.org/10.1109/83.585239Google Scholar

15.

Z. Wanget al., “Image quality assessment: from error visibility to structural similarity,” IEEE Trans. Image Process. 13(4), 600–612 (2004).IIPRE41057-7149http://dx.doi.org/10.1109/TIP.2003.819861Google Scholar

Citation Download Citation

Shahab Chitchian, Markus A. Mayer, Adam Boretsky, Frederik J. van Kuijk, Massoud Motamedi, "Retinal optical coherence tomography image enhancement via shrinkage denoising using double-density dual-tree complex wavelet transform," Journal of Biomedical Optics 17(11), 116009 (1 November 2012). http://dx.doi.org/10.1117/1.JBO.17.11.116009

Access the abstract

JOURNAL ARTICLE
4 PAGES

DOWNLOAD PDF + SAVE TO MY LIBRARY

SHARE

GET CITATION

< Previous Article

|

Next Article >

KEYWORDS

Denoising

Optical coherence tomography

Image enhancement

Wavelet transforms

Signal to noise ratio

Image quality

Wavelets

1.

Introduction

2.

Wavelet Shrinkage Denoising

(1)

2.1.

Double-Density Dual-Tree Complex Wavelet Transform

2.2.

Shrinkage Denoising

3.

Results

Fig. 1

(2)

(3)

(4)

(5)

(6)

(7)

Table 1

a

Fig. 2

4.

Conclusion

Acknowledgments

References

Show All Keywords

KEYWORDS/PHRASES

SEARCH IN:

PUBLICATION YEARS