Cluster-based filtering framework for speckle reduction in OCT images

M. Hossein Eybposh; Zahra Turani; Darius Mehregan; Mohammadreza Nasiriavanaki

doi:10.1364/BOE.9.006359

1. Introduction

Optical coherence tomography (OCT) is an optical imaging modality comparable to ultrasound imaging, except that OCT uses light while ultrasound uses sound waves [1,2]. OCT is used for performing high-resolution cross sectional imaging and works based on low-coherence interferometry [3]. The interferometry relies on the temporal and spatial coherence of optical waves that are backscattered from the tissue [4]. If the central wavelength of the light source is equal to or larger than the scattering compartments within the sample under investigation, the interference of the reflected light with different amplitudes and phases generates a grainy texture in the image called speckle. Speckle degrades the quality of OCT images and conceals the clinically important features [5]. By suppressing the speckle, the quality of the images is improved, and the diagnostically relevant features become more visible.

Methods for speckle reduction are divided into two main categories; hardware based methods, and software based methods [6]. The main hardware-based speckle reduction methods are compounding techniques [7,8]. It has been proven that the averaging successfully reduces the noise by the factor of $\sqrt{N}$ where N is the number of B-scan images to be averaged if the images are sufficiently un-correlated [9]. In 2012, Szkulmowski et al. proposed a shifting beam method that has been utilized for speckle reduction of synthetic aperture radar (SAR), ultrasound, and OCT images [6]. In this method scan beams are shifted orthogonal to both light beam propagation and lateral scanning directions. The images are then averaged. In another method, introduced by Wang et al. in 2013 [10], the probe beam is decentered from the pivot of the scanning mirror to create multiple images that are finally averaged to obtain a single enhanced image [11–13].

Software-based approaches (also called digital filters), process the acquired images offline [14–17]. Total variation (TV) [18], and block matching and 3D filter (BM3D) [19] are two popular de-speckling methods. TV estimates new pixel values by minimizing the amount of variation in the image, ignoring the small scale anatomical structures in the image. Expanding the idea of TV, Wu et al. in 2015 [20], estimated the despeckled image based on local statistics of the speckle. In 2004, Allende et al. proposed a despeckling method that works by detecting outliers in local patches and then cluster analysis within each patch [21]. Each pixel in the patch is then assessed as normal pixel or outlier. In the end, outlier pixels are eliminated while the normal pixels are left with minimal changes. While this method yields adequate images, it does not guarantee edge preservation. Lee [22] developed a local linear minimum mean square error filter, also known as Lee filter (LF), that is a locally adaptive estimation of the Wiener filter (WF). With the assumption of additive noise model, the filter works pixel-wise and estimates the new pixel values based on local statistics. Although the filter successfully degrades noise, it suffers from edge smoothing effects [23]. To avoid edge smoothing, the filter needs to estimate the local statistics in an edge-aware fashion. Jin et al. [23] proposed an adaptive Wiener filter (AWF) that estimates the despeckled pixel values in such a way that avoids over-smoothing of the edges. Assuming that the image and the speckle noise are stationary Gaussian processes, they model the image and the noise with the nonstationary mean and nonstationary variance model [24,25]. In this method, each pixel value is estimated based on local mean and local variance [26]. Although the local mean and variance are determined adaptively, this method still suffers from edge smoothing as it is unable to detect the edges effectively. In 2007, Ozcan et al. discussed several digital filtering methods to decrease the speckle in OCT images [26]. The authors have already implemented six digital filtering methods including enhanced LF [27], hybrid median filter [28], Kuwahara filter, wavelet filtering [29], methods based on artificial neural network [16,17,30] and AWF [31]. From the comparison of the obtained results, they concluded that the enhanced LF and the WF improve the signal-to-noise ratio (SNR) and quality of the OCT images.

The major challenge of current speckle reduction methods is the deterioration of small structures and edge smoothing in the image. The goal of a speckle reduction algorithm is to deconvolve the noise from the original image [4]. Although some algorithms find the optimum solution to the deconvolution problem, the original development of these algorithms require the power spectrum of the noisy image as well as the gold standard (GS) image, which are not available in practice. Practical versions of WF (e.g. LF and AWF) are developed under the assumption that the noise model is additive. This assumption does not apply to multiplicative speckle noise in OCT images. Moreover, LF and AWF assume that the noise and image to be locally stationary. This assumption is only valid in homogenous regions, and not valid around the edges. To restore the edges in the despeckled image, the filtering algorithm needs to effectively detect the edges and avoid applying the filter on the edges. Most current algorithms have not been successful in this regard.

In this study, we developed a cluster-based speckle reduction framework (CSRF) to prevent edge smoothing and small structure deterioration by despeckling in a cluster-wise fashion. For this purpose, we first detect the edges by clustering pixels in the OCT image (using conventional clustering techniques). A despeckling method (i.e. adaptive filtering) is then applied to pixels from individual clusters to preserve the edges. Using this methodology, the borders of the clusters are enhanced, and a pattern similar to gray level quantization may be created in the image that can easily be removed by a mean filtering. This is because the filtering only smooths the clusters and not the cluster edges. This framework is generalizable to any combination of clustering and despeckling methods. This study presents a general framework for OCT despeckling. Therefore, the experiments are focused on demonstrating the effectiveness of this new paradigm.

2. Methodology

2.1 OCT system imaging

The OCT system used in this study is a multi-beam, Fourier-domain, swept-source OCT (Vivosight, Michelson Diagnostic TM Inc., Kent, UK) with a central wavelength of 1305 ± 15 nm. The lateral and axial resolution of our system is 7.5 µm and 10 µm, respectively. The 10-kHz sweep rate determines the time to generate one reflectivity profile. The penetration depth of the system was measured as 1.5 mm in healthy human skin [32].

The OCT is based on multi-beam technology, similar to the technology used in dynamic focus OCT [33], in which four 0.25 mm wide consecutive confocal gates are combined to provide a total confocal gate of 1 mm. Utilizing the multi-beam technology, the images obtained from the four channels are averaged. In OCT, the reflectivity profile is termed as an axial scan (A-scan or A-line). By grouping together several A-lines for different transversal positions of the incident beam on the sample, a cross section image or a B-scan is generated [34]. The images obtained with this OCT system are B-Scan images with a size of 6 mm × 2 mm and software inferred C-scan images with a size of 6 mm × 6 mm. The lateral and axial resolutions of the OCT system are measured as 7.5 and 10 microns, respectively. A schematic of our OCT system is shown in Fig. 1.

Fig. 1 Schematic diagram of the Michelson multi-beam swept-source OCT (SS-OCT); M: mirror. C: optical coupler, M: mirror, PD: photo detector, OA: Optical attenuator.

Download Full Size | PDF

2.2 Cluster-based speckle reduction framework

CSRF is based on clustering pixels into clusters/regions with similar optical properties, then applying the despeckling method on each cluster individually. The benefits of this approach are three-fold: 1) the multiplicative speckle noise in each cluster could be converted to an additive noise, therefore, filtering methods that were originally developed for additive noise could successfully be employed, 2) the image and noise could safely be considered as locally stationary random fields, 3) the edges, which are the borders of two clusters, are well preserved. Eventually, the enhanced images will contain a gray level quantization pattern which can easily be removed by a mean filtering. Figure 2(a) shows how CSRF works.

Fig. 2 Cluster-based Speckle Reduction Framework (CSRF) block diagram. (a) A generalized flowchart of the algorithm, (b) clustering flowchart, (c) despeckling flowchart.

Download Full Size | PDF

2.2.1 Speckle model in CSRF

Assuming the result of the clustering algorithm is regions with homogenous optical properties, the multiplicative speckle noise in individual clusters can be modelled to additive noise. The explanation is as follows. Let us consider $I (i, j)$ as an OCT image pixel value at the location $i, j$ and $S (i, j)$ as the multiplicative speckle noise.

I (i, j) = X (i, j) \times S (i, j)

Taking the logarithm from Eq. (1), the multiplicative noise is converted to an additive noise:

Z (i, j) = \ln (X (i, j)) + \ln (S (i, j))

The Taylor expansion of $I (i, j)$ in an arbitrary point $(i_{0}, j_{0})$ in cluster $k$ is:

\begin{array}{l} I (i, j) = X_{0} \times S_{0} + (i - i_{0}) (S_{0} \times \frac{\partial X (i, j)}{\partial i} |_{i = i_{0}, j = j_{0}} + X_{0} \times \frac{\partial S (i, j)}{\partial i} |_{i = i_{0}, j = j_{0}}) \\ + (j - j_{0}) (S_{0} \times \frac{\partial X (i, j)}{\partial j} |_{i = i_{0}, j = j_{0}} + X_{0} \times \frac{\partial S (i, j)}{\partial j} |_{i = i_{0}, j = j_{0}}) \end{array}

where $X_{0}$ and $S_{0}$ refer to $X (i_{0}, j_{0})$ and $S (i_{0}, j_{0})$ respectively. Pixel intensities in a cluster slightly differ from each other. This variation could be effectively estimated by the Taylor expansion up to the first derivative. For Eq. (2) the Taylor expansion is:

\begin{array}{l} Z (i, j) = \ln (X_{0}) + \ln (S_{0}) + (i - i_{0}) (\frac{\frac{\partial X (i, j)}{\partial i} |_{i = i_{0}, j = j_{0}}}{X_{0}} + \frac{\frac{\partial S (i, j)}{\partial i} |_{i = i_{0}, j = j_{0}}}{S_{0}}) \\ + (j - j_{0}) (\frac{\frac{\partial X (i, j)}{\partial j} |_{i = i_{0}, j = j_{0}}}{X_{0}} + \frac{\frac{\partial S (i, j)}{\partial j} |_{i = i_{0}, j = j_{0}}}{S_{0}}) \end{array}

Integrating Eq. (3) with Eq. (4), we obtained:

Z (i, j) = \ln (X (i, j)) + \ln (S (i, j)) = \frac{1}{X_{0} \times S_{0}} (X_{0} \times S_{0} \times \ln (\frac{X_{0} \times S_{0}}{e}) + I (i, j))

This results in:

I (i, j) = a \times \ln (X (i, j)) + a \times \ln (S (i, j)) - b

where $a = X_{0} \times S_{0}$ and $b = a \times \ln (\frac{a}{e})$ . Assuming that the image is normalized to [-1, 1], Eq. (7) is obtained:

I (i, j) = a \times X (i, j) + a \times S (i, j) - b

If the cluster is a homogenous region with similar pixel values, then $a$ and $b$ are considered equal for all $(i_{0}, j_{0})$ in the cluster. Intuitively, multiplicative noise would change relatively as the pixel intensities change. In other words, higher pixel intensities correspond to higher intensity noise. In a homogenous region where the intensity variation is negligible, the multiplicative noise could be considered uncorrelated with the pixel intensities/location. From Eq. (7) we can see that the multiplicative form of the noisy image (Eq. (1) could be estimated by an additive model. Therefore, all the despeckling methods developed based on additive noise model (e.g., LF, and AWF) could effectively be used in CSRF.

Locally adaptive filtering methods (e.g. LF and AWF), assume the image and noise to be locally stationary. This assumption is only valid in regions with homogenous optical properties [35]. Stationarity is not preserved around the edges since both the image and the multiplicative speckle noise rapidly change near the edges [22,23,35]. In CSRF, on the other hand, each cluster is considered as an individual image with homogenous optical properties. Considering that there are no edges in individual clusters, the image and noise are guaranteed to be locally stationary in each cluster.

2.2.2. Clustering

We used unsupervised clustering methods to cluster pixels in the OCT images. Figure 2(b) demonstrates the clustering algorithm. In the clustering algorithm, each pixel is considered as a data point with two features: pixel intensity (PI) and attenuation coefficient (AC). The AC is estimated for each pixel in the feature extraction phase. The data points are then clustered by k-means (KM) algorithm. The clustering refinement includes filtering the clustering results to eliminate small clusters. The details of the KM clustering algorithm are beyond the scope of this study and we only provide a brief description of the methods. We refer the readers to references [36,37] for a more detailed description of the algorithm.

2.2.2.1. Feature extraction

Each pixel is assigned a set of features, i.e., PI and AC, to describe its optical properties. AC was estimated using the Vermeer et al. [38] approach. Equation (8) estimates the corresponding attenuation coefficient for each pixel.

A C [i, j] = \frac{1}{2 Δ} \ln (1 + \frac{I [i, j]}{\sum_{z = i + 1}^{\infty} I [z, j]})

where

A C [i, j]

is the attenuation coefficient of pixel

i

, at the A-line number

j

,

Δ

is the pixel spacing, and

I [i, j]

is the pixel intensity in the OCT image.

By including pixel locations in the feature set, the clustering method is forced to be sensitive to the position of pixels during clustering. This is in contrast to the goal of clustering in CSRF, which is to detect the edges of homogenous regions. Moreover, column positions in OCT images are of no significant value for clustering, since regions with similar optical properties are normally stretched along the imaging surface. Moreover, the row position for thin and curved layers does not have sufficient discriminability among pixels from different layers/regions.

2.2.2.2. K-means clustering

The K-Means method works by assigning a cluster center to each cluster $i$ , $m_{i}^{(k)}$ , which is the mean of the data points corresponding to that cluster. The algorithm starts ( $k = 0$ ) with $M$ randomly initialized cluster centers, $m_{i}^{(0)}$ . Next, it iterates between the assignment step and the update step. In the assignment step, a data point is assigned one of the clusters according to its distance to cluster center $m_{i}^{(k)}$ . The distance $d (m_{i}^{(k)}, x_{n})$ between point $x_{n}$ and cluster center $m_{i}^{(k)}$ is defined as:

d (m_{i}^{(k)}, x_{n}) = \sqrt{\sum_{m} {(m_{i, m}^{(k)} - x_{n, m})}^{2}}

In the update step, the cluster centers are calculated according to the new clustering results. These two steps are iterated until the cluster centers do not change significantly, i.e for all clusters $i$ :

m_{i}^{(k)} \approx m_{i}^{(k + 1)}

2.2.2.3. Clustering refinement

Clustering algorithms usually do not yield spatially connected clusters. We therefore eliminate small clusters by filtering. The filtering method depends on the selected clustering method. For KM clustering, the clustering results are cluster labels, thus a median filter is suitable.

2.2.3. Despeckling

We used two locally adaptive despeckling algorithms. When processing a pixel from an arbitrary cluster, the CSRF forces the filtering method to disable the effect of pixels from other clusters. In other words, for a filtering method to be incorporated into CSRF, it is necessary that the filtering method allow disabling the effect of pixels from other clusters. This property enables edge preservation, however, limits the type of filtering methods that are suitable for CSRF.

Our experiments focus on demonstrating the effect of integrating conventional filtering methods into CSRF. In order to have a fair evaluation of the CSRF, we qualitatively and quantitatively compared the despeckled images obtained from the CSRF with those obtained from the original despeckling methods. We ignored the comparison between such results, with those obtained from other despeckling methods [18–21]. Figure 2(c) demonstrates the despeckling algorithm.

2.2.3.1 Lee filter

Local linear minimum mean square error filter developed by Lee [22], also known as the Lee filter, was employed in this study. LF is a simplified version of the LATDF [25] that estimates the local mean and local variance of a sliding window, $H (i, j)$ , of size $n$ by $m$ pixels around each pixel located at the $i^{t h}$ row and $j^{t h}$ column of the image, and estimates the despeckled pixel value of $\hat{I} (i, j)$ using Eq. (11):

\hat{I} (i, j) = μ_{x} + \frac{σ_{x}^{2}}{σ_{x}^{2} + σ_{n}^{2}} (I_{(i, j)} - μ_{x})

where $μ_{x}$ and $σ_{x}^{2}$ are the local mean and variance in a window of $n$ by $m$ pixels surrounding the central pixel $(i, j)$ , and $σ_{n}^{2}$ is the noise variance in the image. The weights are then normalized so that they sum to 1. To transform the multiplicative speckle noise into additive noise when applying the original LF and AWF, we take the logarithm of the OCT image if the OCT image is not already logarithmic. This step is omitted when these filtering methods are used in CSRF.

To incorporate the Lee filter into CSRF, the pixels in $H_{i, j}$ are weighted according to Eq. (12) so that only the pixels from the same cluster as that of the pixel in $(i, j)$ are considered when estimating the local mean and variance (i.e., $μ_{x}$ and $σ_{x}^{2}$ ).

ω (i, j, p, q) = {\begin{cases} 1, L (p, q) = L (i, j) \\ 0, L (p, q) \neq L (i, j) or (p, q) = (i, j) \end{cases}

where

ω (i, j, p, q)

corresponds to the weight of the pixel

(p, q)

in the sliding window around the pixel in

(i, j)

, and

L (i, j)

corresponds to the cluster label of the pixel at location

(i, j)

.

The noise variance, $σ_{n}^{2}$ in Eq. (11), is computed for each cluster individually. Estimation of a different noise variance for each cluster accounts for the correlation between noise and pixel intensity, which is due to the multiplicative nature of speckle noise. In the remainder of the manuscript, CSRF-LF refers to LF incorporated into CSRF.

2.2.3.2 Adaptive Wiener filter

Jin et al. proposed a modified version of LF, and called it AWF [23]. In this approach, different weights are assigned to pixels in $H_{i, j}$ according to Eq. (13). These weights are determined based on each pixel’s similarity to the central pixel. The weights were then considered in the calculation of local statistics, i.e., $μ_{x}$ and $σ_{x}^{2}$ :

ω (i, j, p, q) = \frac{K (i, j)}{1 + a (\max [\in^{2}, {(I (i, j) - H_{i, j} (p, q))}^{2}])}

where

a

and

\in^{2}

are parameters of the filter, and

K (i, j)

is the normalization factor such that all weights sum to 1. Moreover, the weight of the central pixel is set to zero. To incorporate AWF into CSRF, we assign zero weight to pixels with a different cluster label than that with the central pixel. In the remainder of this paper, CSRF-AWF refers to AWF incorporated into CSRF. Algorithm 1 provides a detailed description of the algorithm used in CSRF.

Algorithm 1. Cluster-based Speckle Reduction Framework (CSRF)

View Table

2.2.4. Mean filtering

In the final step of CSRF, a mean filtering with 3 by 3 window size is applied on the despeckled images to remove a pattern that we call quantization pattern. This pattern is generated when all the cluster edges are left completely intact. Although some of the cluster edges correspond to the edges of tissue layers, there are cluster edges that do not represent clinically significant edges. The less significant edges are responsible for the layered pattern in the output images (i.e. quantization pattern). These edges are smoothed by mean filtering. We chose a small window size for the mean filter to avoid over-smoothing the more significant edges.

2.3. Denoising assessment

2.3.1 Gold standard denoised image

In order to evaluate the performance of CSRF, we need to compare its denoising performance with the state of the art denoising technique. One of the most straightforward approaches for OCT despeckling is B-scan averaging. In this study the gold standard images are generated by averaging 170 successive scans from the same site. Due to movement artefacts, it is necessary to register all images to a reference image before averaging. In this study, the reference image is chosen arbitrary from the 170 B-scans (see Fig. 4(b) for an example), and other images are registered to it using the enhanced correlation coefficient registration algorithm [39]. The registered B-scans are then averaged (see Fig. 4(a) for an example of a gold standard B-scan). It is important to note that the averaging approach for denoising, while yielding appealing results, is computationally expensive and requires a longer duration of sampling (with minimal movements) and therefore is not a practical solution.

2.3.2 Denoising assessment metrics

In order to compare the results from the CSRF with the GS denoising approach, we used two quantitative assessment metrics including signal-to-noise ratio (SNR), and structural similarity index (SSIM). The equations explaining these quality metric measures are provided in Eq. (14) and Eq. (15), respectively.

PSNR compares the signal of the OCT image to its background noise [40]. SSIM score measures the image quality based on structural similarity between the GS and despeckled images [41].

P S N R = 10 \log_{10} (\frac{\max (I^{2})}{σ_{n}^{2}})

S S I M (I_{G S}, \hat{I}) = \frac{(2 μ_{I_{G S}} μ_{\hat{I}} + C_{1}) (2 σ_{I_{G S}} μ_{\hat{I}} + C_{2})}{(μ_{I_{G S}}^{2} + μ_{\hat{I}}^{2} + C_{1}) (σ_{I_{G S}}^{2} + σ_{\hat{I}}^{2} + C_{2})}

where

I_{G S}

,

I

, and

\hat{I}

are GS, noisy, and the estimated (despeckled) images, respectively.

σ_{n}^{2}

is the speckle noise variance. In Eq. (15),

μ_{I_{G S}}

and

μ_{\hat{I}}

are the mean of the GS and despeckled images respectively.

σ_{I_{G S}}

is the covariance between the GS and the despeckled images.

σ_{I_{G S}}^{2}

and

σ_{\hat{I}}^{2}

are variance of GS and despeckled images respectively,

C_{1}

and

C_{2}

are constants;

C_{1} = 6.5025

and

C_{2} = 58.5225

[41].

3. Results and discussion

To evaluate the performance of the CSRF, OCT skin images were acquired from fourteen different body sites of a healthy, 25-year-old, male volunteer. The OCT machine is FDA approved. The institutional review board at Wayne State University (Independent Investigational Review Board, Detroit, MI) approved the study protocol, and informed consent was obtained from the patient before enrollment in the study. The body sites included anterior neck, buccal region, calves, chest, dorsum of foot, dorsum of hand, ears, forearm, forehead, lips, nose, orbit, palm, and upper back. For each of the fourteen body sites, 170 B-scans were acquired from 4 to 5 different locations, generating 56 data sets for evaluation. The proposed de-speckling method was applied on only one of the B-scans in each data set, forming 56 B-scans in total.

The initial experiments showed that the number of clusters significantly affected performance of the framework. In order to determine the optimum number of clusters, we investigated the effect of the number of clusters in a small subset of the OCT data sets (10 data sets). In the initial experiments, the quantitative and qualitative results indicated that with small numbers of clusters (less than 10) the major edges were preserved; but the small structures were deteriorated. We hypothesized that the number of clusters should be large enough so that the specific patterns of the small structures could be clustered separately. Our results showed that the optimum performance was achieved for 20 clusters per image. However, we believe that the optimum number of clusters might vary between imaging device since each imaging device has a different level of noise. Also, we believe that the optimum number of clusters is subject-independent, due to the anatomical similarities between subjects. We recommend testing different numbers of clusters when images from other body sites are used.

Our results from 56 data sets showed that on average, integration of AWF with CSRF, improves the SNR, and SSIM metrics by 13.63 dB, and 0.04, respectively (see Fig. 3(a) and (b)). The results also show an average improvement of 13.88 dB, and 0.05 in SNR, and SSIM for integration of LF with CSRF (see Fig. 3(c) and (d)).

Fig. 3 Results of SNR, and SSIM on 61 OCT image data sets when CSRF used, in comparison with when the original despeckling algorithms used. AWF refers to Jin et al. (23) adaptive Wiener filter, and LF stands for Lee filter. (a) PSNR calculated for AWF and CSRF-AWF, (b) SSIM calculated for AWF and CSRF-AWF, (c) PSNR calculated for LF and CSRF-LF, (d) SSIM calculated for LF and CSRF-LF.

Download Full Size | PDF

When testing different window sizes (i.e., 3 × 3 to 13 × 13), we observed that a window size of 9 by 9 pixels for the despeckling methods yields optimum qualitative and quantitative results. Window sizes smaller than 9 by 9 pixels did not effectively improve the quality of the images. Similarly, window sizes greater than 9 by 9 pixels smoothed the edges and deteriorated small structures. It is worth mentioning that the window size 9 by 9 pixels may not be appropriate for OCT images of other organs, e.g. retina. For the mean filter in the third and final step of CSRF, a 3 × 3 window was used to smooth the borders of the clusters and alleviate the problem of quantization pattern. Please note that a window size of 3 × 3 has a negligible effect on major edges, i.e., the ones that are diagnostically important.

As shown in Fig. 4(c), the results of the AWF filtering show a noisy pattern around the edges. A checkerboard pattern is also observed throughout the image, which might be due to the inability of AWF’s pixel weighting approach (see Eq. (13) to differentiate between edges and speckle. Moreover, the small structures, shown in the red box in Fig. 4(c) is deteriorated. Notably, the CSRF-AWF seems to solve the problem of checkerboard pattern. We can say that CSRF-AWF outperforms the original AWF method, both in edge preservation and smoothing the homogenous regions (see Fig. 4(d)). Figure 4(e) shows that CSRF-AWF significantly outperforms AWF in eliminating speckle. The major edges are deteriorated after the application of AWF, while CSRF-AWF leaves significant edges almost intact. As shown in Fig. 4(f), the original LF does not effectively reduce speckle around edges or small structures and shows that the despeckled image contains small black dots. In the results of CSRF-LF however (Fig. 4(g)), such patterns have been removed. One can say that the images despeckled by CSRF-LF are smoother than those obtained from the original LF. This improvement is not however as significant as that seen in CSRF-AWF compared to AWF. The results from Fig. 4(h) indicate that LF fails to eliminate the speckle pattern (manifested in the form of rapid changes of pixel intensities in the A-line) and deteriorates the edges, while CSRF-LF completely eliminates the speckle and preserves the major edges from deterioration. In Fig. 5, despeckling results of several OCT images of skin taken from different body sites of a 25-year-old healthy male are shown.

Fig. 4 Despeckling results of an OCT skin image taken from the palm of a 25-year-old healthy male. The window size was set to (9 × 9) pixels for all the despeckling methods. (a) GS image, (b) an original (noisy) OCT image, despeckled using (c) AWF, (d) CSRF-AWF, (f) LF, (g) CSRF-LF. (e,h) comparison of a part of an A-line #100 profile in the original, CSRF-AWF, and CSRF-LF methods, respectively. Vertical yellow line on (a) corresponds to the A-Line #100.

Download Full Size | PDF

Fig. 5 Despeckling results of several OCT images of skin taken from different body sites of a 25-year-old healthy male. The window size was set to (9 × 9) pixels for all the despeckling methods. (a,b,c,d,) original (noisy) OCT images, (e,f,g,h) despeckled using CSRF-LF. The PNSR improvement in image a, b, c, and d, are 11.5 dB, 10.1 dB, 12.2 dB, and 9.2 dB, respectively.

Download Full Size | PDF

3.1. Computational performance of CSRF

The computational complexity of CSRF is determined as the sum of the complexity of the clustering algorithm and that of the despeckling algorithm. The complexity of KM clustering method is in the order of $0 (i k n d)$ . Where $n$ is the number of pixels, $k$ the number of clusters, $i$ the number of iterations until the clustering is converged, and $d$ the dimensionality of input features. Although CSRF increases the complexity of the despeckling algorithms, the clustering algorithms are highly parallelable and there are a number of multi-core implementations of these methods available [42–44]. Moreover, there are implementations of the KM method that use graphical processing units [42,43], with which the execution time of the algorithm is significantly reduced.

The computational complexity of the filters (AWF and LF) does not depend on the number of clusters or the structure of the image (i.e. the shape of the clusters). Regardless of these factors, each pixel is processed only once therefore the computational complexity of the filtering step relies merely on the filtering approach that is going to be used (in this case LF, and AWF). On average the despeckling of each image took 5.3 seconds on personal computer with Core-i7 6700K central processing unit and 32 giga bytes of memory. Processing each image using the averaging method (the approach that we used to get the GS images) took 42 minutes on average.

3.2. Final remarks and future work

Although speckle (noise) decreases the image quality, blurs the image and conceals the diagnostically relevant features in OCT images, it carries submicron structural information of the tissue being imaged. Therefore an intelligent despeckling algorithm is required to make the image more eligible while preserving the major features in the image. In this study, we developed a cluster-based filtering framework, called CSRF. The framework enhances the performance of despeckling algorithms by applying the despeckling algorithm on the regions with the same optical properties: clusters. Our results show that images despeckled with our proposed framework are qualitatively and quantitatively improved (see Fig. 3 and Fig. 4). Adding other statistical features of OCT images to the feature vector of clustering, and including first and higher order statistics, is something that should be explored in future works. Other optical properties such as the scattering coefficient, anisotropy factor, and geometrical properties such as shape and thickness can also be used in clustering algorithms to improve its performance. We observed that the integration of a common despeckling method, e.g. AWF, in CSRF significantly increases its despeckling capabilities (see Fig. 3 and Fig. 4). Visual inspection of the despeckled images validated the results from qualitative stand point. Another important finding was the consistent improvement of both quality metrics, PSNR and SSIM when the original filtering methods were used in CSRF.

However, the scope of this study is merely to introduce CSRF as a de-speckling framework and validate its effectiveness, we believe that the outperformance of CSRF-based methods compared to the original methods, indicates the great potential of this framework for despeckling. This framework is proposing a new viewpoint on OCT despeckling. Future studies could focus on developing new filtering methods that are tailored for CSRF.

4. Conclusion

In this paper, we proposed a cluster-based speckle reduction framework (CSRF) for OCT images to reduce the speckle. The method was successfully tested on 56 sets of OCT images of human skin, in vivo. The results showed an average improvement of 13.88 dB, and 0.05 for Lee filtering, and 13.63 dB, and 0.04 for adaptive Wiener filter in PSNR and SSIM, respectively. The proposed method was tested on OCT images of skin, however, it could easily be used for OCT images of other sites, e.g., retina. This framework, if a more sophisticated clustering algorithm is used, helps in further functionalizing speckle reduction algorithms to enhance the visibility of a specific characteristic in the OCT images.

Funding

Michelson Diagnostics; Wayne State University Startup fund.

Acknowledgments

This project has been partially funded by Michelson Diagnostics and partially by Wayne State University Startup fund.

Disclosures

The authors have no relevant financial interests in the manuscript and no other potential conflicts of interest to disclose.

References

1. J. G. Fujimoto, C. Pitris, S. A. Boppart, and M. E. Brezinski, “Optical coherence tomography: an emerging technology for biomedical imaging and optical biopsy,” Neoplasia 2(1-2), 9–25 (2000). [CrossRef] [PubMed]

2. S. Adabi, Z. Turani, E. Fatemizadeh, A. Clayton, and M. Nasiriavanaki, “Optical Coherence Tomography Technology and Quality Improvement Methods for Optical Coherence Tomography Images of Skin: A Short Review,” Biomed. Eng. Comput. Biol. 8(3), 1179597217713475 (2017). [CrossRef] [PubMed]

3. D. Huang, E. A. Swanson, C. P. Lin, J. S. Schuman, W. G. Stinson, W. Chang, M. R. Hee, T. Flotte, K. Gregory, C. A. Puliafito, and et, “Optical coherence tomography,” Science 254(5035), 1178–1181 (1991). [CrossRef] [PubMed]

4. J. M. Schmitt, S. H. Xiang, and K. M. Yung, “Speckle in optical coherence tomography,” J. Biomed. Opt. 4(1), 95–105 (1999). [CrossRef] [PubMed]

5. J. W. Goodman, Speckle Phenomena in Optics: Theory and Applications (Roberts and Company Publishers, 2007).

6. M. Szkulmowski, I. Gorczynska, D. Szlag, M. Sylwestrzak, A. Kowalczyk, and M. Wojtkowski, “Efficient reduction of speckle noise in Optical Coherence Tomography,” Opt. Express 20(2), 1337–1359 (2012). [CrossRef] [PubMed]

7. M. R. Avanaki, A. G. Podoleanu, M. C. Price, S. A. Corr, and S. A. Hojjatoleslami, “Two applications of solid phantoms in performance assessment of optical coherence tomography systems,” Appl. Opt. 52(29), 7054–7061 (2013). [CrossRef] [PubMed]

8. M. R. Avanaki, R. Cernat, P. J. Tadrous, T. Tatla, A. G. Podoleanu, and S. A. Hojjatoleslami, “Spatial compounding algorithm for speckle reduction of dynamic focus OCT images,” IEEE Photonics Technol. Lett. 25(15), 1439–1442 (2013). [CrossRef]

9. J. Cheng, L. Duan, D. W. K. Wong, D. Tao, M. Akiba, and J. Liu, “Speckle Reduction in Optical Coherence Tomography by Image Registration and Matrix Completion,” presented at the Medical Image Computing and Computer-Assisted Intervention, 2014.

10. J. H. Chen, C. S. Huang, K. C. Chien, E. Takada, W. K. Moon, J. H. Wu, N. Cho, Y. F. Wang, and R. F. Chang, “Breast density analysis for whole breast ultrasound images,” Med. Phys. 36(11), 4933–4943 (2009). [CrossRef] [PubMed]

11. N. Iftimia, B. E. Bouma, and G. J. Tearney, “Speckle reduction in optical coherence tomography by “path length encoded” angular compounding,” J. Biomed. Opt. 8(2), 260–263 (2003). [CrossRef] [PubMed]

12. H. Zhang, Z. Li, X. Wang, and X. Zhang, “Speckle reduction in optical coherence tomography by two-step image registration,” J. Biomed. Opt. 20(3), 036013 (2015). [CrossRef] [PubMed]

13. J. J. G. Valverde, et al., “Speckle reduction process based on digital filtering and wavelet compounding in optical coherence tomography for dermatology,” presented at the European Conferences on Biomedical Optics, 2015.

14. S. Adabi, H. Mohebbikarkhoran, D. Mehregan, S. Conforto, and M. Nasiriavanaki, “An Intelligent Despeckling method for swept source optical coherence tomography images of skin,” in SPIE Medical Imaging, 2017, vol. 3, pp. 101372–101376: International Society for Optics and Photonics.

15. M. R. Pishgoo, M. R. N. Avanaki, and R. Ebrahimpour, “The application of multi-layer artificial neural networks in speckle reduction (methodology),” Journal of Electrical and Computer Engineering Innovations 2(1), 37–42 (2014).

16. M. R. Avanaki, P. P. Laissue, T. J. Eom, A. G. Podoleanu, and A. Hojjatoleslami, “Speckle reduction using an artificial neural network algorithm,” Appl. Opt. 52(21), 5050–5057 (2013). [CrossRef] [PubMed]

17. M. Avanaki, P. P. Laissue, and A. Hojjatoleslami, “De-Noising Speckled Optical Coherence Tomography Images Using an Algorithm Based on Artificial Neural Network,” J. Neurosci. Neuroeng. 2(4), 347–352 (2013). [CrossRef]

18. L. I. Rudin, S. Osher, and E. Fatemi, “Nonlinear total variation based noise removal algorithms,” Physica D 60(1), 259–268 (1992). [CrossRef]

19. B. Chong and Y.-K. Zhu, “Speckle reduction in optical coherence tomography images of human finger skin by wavelet modified BM3D filter,” Opt. Commun. 291, 461–469 (2013). [CrossRef]

20. T. Wu, Y. Shi, Y. Liu, and C. He, “Speckle reduction in optical coherence tomography by adaptive total variation method,” J. Mod. Opt. 62(21), 1849–1855 (2015).

21. H. Allende and J. Galbiati, “A non-parametric filter for digital image restoration, using cluster analysis,” Pattern Recognit. Lett. 25(8), 841–847 (2004). [CrossRef]

22. J.-S. Lee, “Digital image enhancement and noise filtering by use of local statistics,” IEEE Trans. Pattern Anal. Mach. Intell. 2(2), 165–168 (1980). [CrossRef] [PubMed]

23. F. Jin, P. Fieguth, L. Winger, and E. Jernigan, “Adaptive Wiener filtering of noisy images and image sequences,” in International Conference on Image Processing, (IEEE, 2003), volume 2, pp. 349–352. [CrossRef]

24. D. T. Kuan, A. A. Sawchuk, T. C. Strand, and P. Chavel, “Adaptive noise smoothing filter for images with signal-dependent noise,” IEEE Trans. Pattern Anal. Mach. Intell. 7(2), 165–177 (1985). [CrossRef] [PubMed]

25. L. Yaroslavsky, K. Egiazarian, and J. Astola, “Transform domain approaches for image denoising,” J. Electron. Imaging 11(2), 149–156 (2002). [CrossRef]

26. A. Ozcan, A. Bilenca, A. E. Desjardins, B. E. Bouma, and G. J. Tearney, “Speckle reduction in optical coherence tomography images using digital filtering,” J. Opt. Soc. Am. A 24(7), 1901–1910 (2007). [CrossRef] [PubMed]

27. A. Lopes, R. Touzi, and E. Nezry, “Adaptive speckle filters and scene heterogeneity,” IEEE Trans. Geoscience and Remote Sensing 28(6), 992–1000 (1990). [CrossRef]

28. A. Nieminen, P. Heinonen, and Y. Neuvo, “A new class of detail-preserving filters for image processing,” Pattern Analysis and Machine Intelligence, IEEE Transactions on, no. 1, (IEEE, 1987), pp. 74–90.

29. J.-L. Starck, F. D. Murtagh, and A. Bijaoui, Image Processing and Data Analysis: The Multiscale Approach (Cambridge University Press, 1998).

30. S. Adabi, S. Conforto, A. Clayton, A. G. Podoleanu, A. Hojjat, and M. Avanaki, “An intelligent speckle reduction algorithm for optical coherence tomography images,” International Conference on Photonics, Optics and Laser Technology, PHOTOPTICS 2016, pp. 38–45, 2016. [CrossRef]

31. Y. Funama, Y. Noguchi, and M. Shimamura, “Reduction of artifacts in degraded CT image by adaptive Wiener filter,” Journal of Med. Electron. Biol. Eng. 40(1), 1–6 (2002).

32. J. Holmes and J. Welzel, “OCT in Dermatology,” in Optical Coherence Tomography: Technology and Applications (2015), 2189–2207.

33. M. R. N. Avanaki and A. Podoleanu, “En-face time-domain optical coherence tomography with dynamic focus for high-resolution imaging,” J. Biomed. Opt. 22(5), 56009 (2017). [CrossRef] [PubMed]

34. A. M. Schmitt, “Principles and application of optical coherent tomography in dermatology,” Dermatology (Basel) 217(1), 12–13 (2008). [CrossRef] [PubMed]

35. J. M. Schmitt, “Optical coherence tomography (OCT): a review,” IEEE J. Sel. Top. Quantum Electron. 5(4), 1205–1215 (1999). [CrossRef]

36. D. J. MacKay, Information Theory, Inference and Learning Algorithms (Cambridge University Press, 2003).

37. J. C. Bezdek, R. Ehrlich, and W. Full, “FCM: The fuzzy c-means clustering algorithm,” Comput. Geosci. 10(2–3), 191–203 (1984). [CrossRef]

38. K. A. Vermeer, J. Mo, J. J. Weda, H. G. Lemij, and J. F. de Boer, “Depth-resolved model-based reconstruction of attenuation coefficients in optical coherence tomography,” Biomed. Opt. Express 5(1), 322–337 (2014). [CrossRef] [PubMed]

39. G. D. Evangelidis and E. Z. Psarakis, “Parametric image alignment using enhanced correlation coefficient maximization,” IEEE Trans. Pattern Anal. Mach. Intell. 30(10), 1858–1865 (2008). [CrossRef] [PubMed]

40. G. Andria, F. Attivissimo, A. M. Lanzolla, and M. Savino, “A suitable threshold for speckle reduction in ultrasound images,” IEEE Trans. Instrum. Meas. 62(8), 2270–2279 (2013). [CrossRef]

41. Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, “Image quality assessment: from error visibility to structural similarity,” IEEE Trans. Image Process. 13(4), 600–612 (2004). [CrossRef] [PubMed]

42. M. Al-Ayyoub, A. M. Abu-Dalo, Y. Jararweh, M. Jarrah, and M. Al Sa’d, “A gpu-based implementations of the fuzzy c-means algorithms for medical image segmentation,” J. Supercomput. 71(8), 3149–3162 (2015). [CrossRef]

43. R. Farivar, D. Rebolledo, E. Chan, and R. H. Campbell, “A Parallel Implementation of K-Means Clustering on GPUs,” PDPTA 13(2), 212–312 (2008).

44. K. Kerdprasop and N. Kerdprasop, “Parallelization of k-means clustering on multi-core processors,” Proceedings of the 10th WSEAS International Conference on Applied Computer Science,ACS, vol. 10 (2010), pp. 472–477.

Initialization $I [i, j]$	intensity of a pixel at row $i$ and column $j$ in the original OCT image
$A C [i, j]$	The estimated attenuation coefficient of a pixel at row $i$ and column $j$ in the original OCT image
$Δ$	pixel size in mm
$R$	number of pixels in each A-line (number of rows)
$C$	number of A-lines (number of columns)
$N$	total number of pixe^{ls in the OCT image} $(R \times C)$
$F$	feature vector including ( $I [i, j]$ , $A C [i, j]$ )
$K$	desired number of clusters
$C L_K M$	KM clustering result
$H_{i, j}$	sliding window around pixel $[i, j]$ in the original image
$N 1 \times N 2$	size of sliding window
$W_{i, j}$	Weight matrix ^{for pixels in} $H_{i, j}$
$n_{l}$	Number of pixels in $H_{i, j}$ that belong to the same cluster as the central pixel
$1_{3 * 3}$	A matrix of ones with size 3*3
$A v e r a g e (d a t a, w e i g h t s) = 1 / N_{d} \sum d a t a \times w e i g h t$
1. Feature Extraction 1.1. For each pixel in position $[i, j]$ , calculate and assign a set of features. We only use the attenuation coefficient [38]: $A C [i, j] = \frac{1}{2 Δ} \ln (1 + \frac{I [i, j]}{\sum_{z = i + 1}^{\infty} I [z, j]})$ 1.2. For each pixel at position $[i, j]$ , assign the feature vector $F [i, j]$ : $F [i, j] = (I [i, j], A C [i, j])$ 2. Clustering 2.1. Apply KM clustering which is $C L_K M = K M (F, K)$ 2.2. Apply median filter on the clustering results: $C L_K M = M e d i a n_F i l t e r (C L_K M)$ 3. Apply Filtering: for a sliding window of size $N 1 \times N 2$ around pixel $[i, j]$ , $H_{i, j}$ 3.1. Estimate the noise variance for each individual cluster $σ_{n k}^{2}$ 3.2. Compute the weight matrix $W_{i, j}$ for all the pixels in $H_{i, j}$ . The elements of $W_{i, j}$ _, $ω (i, j, p, q)$ are determined as follows: 3.2.1. Lee filter (LF): all the pixels in the window are weighted the same $ω (i, j, p, q) = 1 / N 1 \times N 2$ 3.2.2. Jin’s adaptive Wiener filter (AWF) $ω (i, j, p, q) = \frac{K (i, j)}{1 + a (\max [\in^{2}, {(I (i, j) - H_{i, j} (p, q))}^{2}])}$ 3.2.3. Cluster-based LF (CSRF-LF) $ω (i, j, p, q) = {\begin{cases} 1, L (p, q) = L (i, j) \\ 0, L (p, q) \neq L (i, j) or (p, q) = (i, j) \end{cases}$ 3.2.4. Cluster-based adaptive Wiener filter (CSRF-AWF) $ω (i, j, p, q) = \frac{K (i, j)}{{\begin{cases} 1 + a (\max [\in^{2}, {(I (i, j) - H_{i, j} (p, q))}^{2}]), \| L (p, q) = L (i, j) \\ 0, \| L (p, q) \neq L (i, j) o r (p, q) = (i, j) \end{cases}}$ 3.3. For each window, calculate the local statistics: $μ_{i, j} = A v e r a g e (H_{i, j}, W_{i, j})$ $σ_{i, j}^{2} = A v e r a g e (H_{i, j} - μ {[i, j]}^{2}, W_{i, j})$ 3.4. Calculate the new pixel value: $\hat{I} (i, j) = μ_{x} + \frac{σ_{x}^{2}}{σ_{x}^{2} + σ_{n}^{2}} (I_{(i, j)} - μ_{x})$ 4. Smoothing the cluster edges with a 3 × 3 Mean filter: $\hat{I} (i, j) = A v e r a g e (H_{i, j}, 1_{3 * 3})$

Cluster-based filtering framework for speckle reduction in OCT images

Abstract

1. Introduction

2. Methodology

2.1 OCT system imaging

2.2 Cluster-based speckle reduction framework

2.2.1 Speckle model in CSRF

2.2.2. Clustering

2.2.2.1. Feature extraction

2.2.2.2. K-means clustering

2.2.2.3. Clustering refinement

2.2.3. Despeckling

2.2.3.1 Lee filter

2.2.3.2 Adaptive Wiener filter

2.2.4. Mean filtering

2.3. Denoising assessment

2.3.1 Gold standard denoised image

2.3.2 Denoising assessment metrics

3. Results and discussion

3.1. Computational performance of CSRF

3.2. Final remarks and future work

4. Conclusion

Funding

Acknowledgments

Disclosures

References

Cited By

Figures (5)

Tables (1)

Equations (15)

Biomedical Optics Express