2.1.

2-D-DWT

2-D-DWT was used for image decomposition. In detail, the speckled image is decomposed, and the sub-band components are obtained, as shown in Fig. 1(a), where $C_{\mathrm{A}}$ is the approximation component and $C_{\mathrm{H}}$, $C_{\mathrm{V}}$, and $C_{\mathrm{D}}$ are the detail components of the horizontal, vertical, and diagonal directions, respectively. $C_{\mathrm{A}}i$, $C_{\mathrm{H}}i$, $C_{\mathrm{V}}i$, and $C_{\mathrm{D}}i$ are the corresponding coefficients at the decomposition level $i$. The orthogonal fourth-order Daubechies wavelet was utilized as it has longer filter lags, which can help to reduce the pixilation effect by bringing smoothing effects.¹⁶ The purpose of doing wavelet decomposition is to decompose the image speckle noise into different levels. Wavelet coefficients corresponding to true edge features tend to cluster spatially across multiple decomposition levels, whereas coefficients resulting from the speckle noise tend to be approximately decorrelated.^14,21 Because the multilevel DWT is a linear transform, speckle noise, which is additive in logarithmically processed OCT images, is still additive in the wavelet domain. The sub-band images, therefore, can be compounded across different decomposition levels to achieve speckle reduction by averaging.⁵

Fig. 1

(a) Ordering of the approximation and detailed coefficients of a three-level two-dimensional discrete wavelet transform (2-D-DWT) and (b) step-by-step flowchart of the proposed method (three-level DWT). Solid arrows represent 2-D-DWT. Dashed arrows represent the procedure for compounding process and 2-D-IDWT.

2.2.

Compounding

Figure 1(b) shows the proposed compounding strategy flowchart for a multilevel (three-level in the figure) 2-D wavelet analysis. For convenience of illustration, the detailed coefficients at the decomposition level of $i$, $C_{\mathrm{H}}i$, $C_{\mathrm{V}}i$, and $C_{\mathrm{D}}i$, are represented by $C_{\mathrm{S}}i$. However, $C_{\mathrm{H}}i$, $C_{\mathrm{V}}i$, and $C_{\mathrm{D}}i$ are handled separately, and any operation on the $C_{\mathrm{S}}i$ is equivalent to the operation on $C_{\mathrm{H}}i$, $C_{\mathrm{V}}i$, and $C_{\mathrm{D}}i$, respectively. Assuming that three-level 2-D-DWT is performed, the compounding strategy $C(·)$ in Eq. (2) has the following detailed procedures [denoted by dashed arrows in Fig. 1(b)].

After the above stage of compounding, speckle noise will be effectively reduced in the wavelet domain. This stage achieves the same goal as the thresholding stage for the conventional wavelet-based thresholding techniques, which reduce the speckle noise by thresholding the transformed coefficients.^14,16–18

2.3.

2-D-IDWT

$C_{\mathrm{A}}{1}^{\prime }$, $C_{\mathrm{H}}{1}^{\prime }$, $C_{\mathrm{V}}{1}^{\prime }$, and $C_{\mathrm{D}}{1}^{\prime }$ correspond to a sub-band component with speckle noise significantly reduced. At last, single-level 2-D-IDWT is performed directly on $C_{\mathrm{A}}{1}^{\prime }$, $C_{\mathrm{H}}{1}^{\prime }$, $C_{\mathrm{V}}{1}^{\prime }$, and $C_{\mathrm{D}}{1}^{\prime }$ to reconstruct a despeckled image.

3.1.

Feasibility of WDC

To demonstrate the effectiveness of the proposed WDC algorithm, we imaged a human fingertip in vivo using a swept-source OCT system operating at the center wavelength of 1550 nm.²⁴ The obtained image, as shown in Fig. 2, has an axial resolution of $\sim 12\mu \text{m}$.

Fig. 2

(a) Fingertip image without despeckling and (b) despeckled fingertip image. Arrows point to sweat glands.

Figure 2(a) shows the original fingertip image, and Fig. 2(b) demonstrates the despeckled fingertip image processed by the proposed WDC method. Results show substantial smoothing of the grainy appearance due to the existence of speckle noise and enhanced visualization of structures such as sweat glands indicated by white arrows. It should be noted that strong speckle reduction is applied in order to greatly enhance the discernibility of the sweat glands from the surrounding speckle noise, and one can observe some blurring effects consequently. The blurring effect can be reduced by applying slight speckle reduction. However, slight speckle reduction will in turn reduce the visibility of the sweat glands. The trade-off between speckle reduction and edge preservation will be discussed in detail in the subsequent sections.

To further explain the effects of each stage of the algorithm on the output, we show the sub-band images obtained after the compounding step (bottom row in Fig. 3). For comparison, we also show the corresponding sub-band images by applying single-level 2-D-DWT to the original image (top row in Fig. 3). As shown in Fig. 3, compared with the original sub-band images, the newly obtained sub-band images show clear speckle reduction. Moreover, the edge information is also preserved at the same time. It indicates that the compounding step can effectively reduce speckle noise in the wavelet domain with simultaneous edge preservation.

Fig. 3

Sub-band images of the original fingertip image after wavelet decomposition (top row). Obtained sub-band images after the compounding process (bottom row). A: approximation sub-band; H: horizontal sub-band; V: vertical sub-band; and D: diagonal sub-band.

3.2.

WDC for IVOCT

Subsequently, in order to demonstrate the applicability of the proposed WDC algorithm to reduce speckle noise in the context of intravascular imaging where the visibility of image morphological structures is severely deteriorated by speckle noise, in vivo intravascular OCT imaging of porcine carotid artery image was obtained using a commercial Lightlab C7-XR Fourier Domain OCT system. The details of animal imaging protocol have been previously described by Cheng et al.¹³ The original rotary scanning image ($1937\times 1937$) and the despeckled image based on the proposed WDC method are demonstrated in Figs. 4(a) and 4(c), respectively. The image processed by wavelet-based soft global thresholding²⁰ is also showed in Fig. 4(b) for comparison. All analyses were executed on an Intel Core i5-750, 2.67-GHz desktop personal computer using MATLAB 7.11.0 (MathWorks, Natick, Massachusetts).

Fig. 4

(a) Original intravascular optical coherence tomography (IVOCT) image of porcine carotid artery. Yellow rectangle indicates the noise region used in the metrics calculation. Green rectangle indicates the zoomed region in Fig. 5. Red rectangle indicates the signal regions. Blue rectangles are used for equivalent number of looks (ENL) calculations; (b) despeckled IVOCT image of porcine carotid artery based on the soft global thresholding; and (c) despeckled IVOCT image of porcine carotid artery based on the proposed wavelet domain compounding (WDC). For direct comparison, the images are shown on the same gray scale.

Most of the features were well preserved in the despeckled image. The grainy appearance was mostly smoothed, which makes some hidden features in the original image more obvious. The distinctive ring structure encircling the vessel lumens are the arterial vessel walls. The external elastic lamina, as pointed by the white arrow in Fig. 4(c), is visualized with a much clearer discernibility after speckle reduction. The media-adventitia border is also obviously observed in the despeckled image. For comparison, we have also performed conventional wavelet-based soft global thresholding,²⁰ as shown in Fig. 4(b). The threshold is chosen to be 2.4 times the noise variance, which is obtained from the robust median estimator of the highest sub-band of the wavelet transform.^20,25 It should be noted that the number 2.4 is chosen to obtain the optimal trade-off between speckle reduction and edge preservation. The detailed procedures have been described in other reports.^13,14 In brief, by examining the trend of the speckle reduction performance [such as contrast-to-noise ratio (CNR)] and edge preservation parameter for various choices of threshold and choosing the optimal trade-off point, the number 2.4 is obtained. For the subsequent images, the optimal thresholds for different images are different, but the method finding the optimal threshold follows the above-discussed way.

To better appreciate the performance of the WDC on speckle suppression, Fig. 5 shows the enlarged view of the green region of interest (ROI) for the original image, image after speckle reduction based on soft global thresholding, and image after speckle reduction based on the proposed WDC, respectively. The despeckled image based on WDC is much clearer than the original one. Most importantly, this is achieved when the image edge sharpness is well preserved demonstrating its ability of the algorithm to preserve the signals, while attenuating speckle noise. Structures with clearer visualization are indicated by red arrows in Fig. 5(c). While Fig. 5(b) may also achieve comparable speckle reduction, it shows reduced edge preservation compared with Fig. 5(c). When the soft global thresholding technique is employed, speckle noise may be further reduced by setting a larger threshold, but at the cost of edge preservation.

Fig. 5

(a) Enlarged view of the green region of interest (ROI) in the original image; (b) despeckled image processed by wavelet-based soft global thresholding method; and (c) despeckled image processed by the proposed WDC algorithm. The red arrows indicate structures with clearer discernibility in the processed image.

The output of the proposed speckle reduction algorithm is dependent on three main parameters: (a) the wavelet type used for the wavelet decomposition, i.e., Daubechies (db) and symlets (sym), (b) the interpolation method for up-sampling, i.e., linear and cubic, and (c) the level of the wavelet decomposition. The effects of the level of the wavelet decomposition will be discussed in detail in the following section. For the first two parameters, we have tried to vary the wavelet bases and interpolation methods to study the possible effects on the final output. While both the wavelet type and the interpolation method would affect the final speckle reduction performance, results indicate that different wavelet bases or different interpolation methods did not have a noticeable change in the image quality. There was no observable difference in the image quality in terms of wavelet types or interpolation methods. This also agrees well with what has been observed for the conventional wavelet-based thresholding techniques, which is that the most important factor in wavelet denoising is the decomposition level rather than the wavelet or threshold type.²⁶

To quantitatively evaluate the performance of the proposed algorithm, image quality metrics were calculated over ROIs including high signal regions (media), low signal regions (adventitia), and a noise background, which have been depicted in Fig. 4(a). The CNR was defined as $\mathrm{CNR}=10\log ({\mu }_{\mathrm{r}}-{\mu }_{\mathrm{b}})/\sqrt{{\sigma }_{\mathrm{r}}^2+{\sigma }_{\mathrm{b}}^2}$, and the equivalent number of looks (ENL) was defined as $\mathrm{ENL}={\mu }_{\mathrm{r}}^2/{\sigma }_{\mathrm{r}}^2$, where ${\mu }_{\mathrm{r}}$ and ${\mu }_{\mathrm{b}}$ are the mean of the intensity from a set of ROIs and a background noise region, respectively. ${\sigma }_{\mathrm{r}}$ and ${\sigma }_{\mathrm{b}}$ are the standard deviation of the intensity over the ROIs and the background noise region, respectively. The CNR measures the contrast between image features [red and blue ROIs in Fig. 4(a)] and an area of background noise (yellow ROI), while the ENL measures smoothness in areas, which should have a homogeneous appearance (blue ROIs).¹⁴ The metrics were calculated as the average over the ROIs used. In addition, the global signal-to-noise ratio (SNR) was calculated as $\mathrm{SNR}=10\log [\max {(X_{\mathrm{lin}})}^2/{\sigma }_{\mathrm{lin}}^2]$, where $X_{\mathrm{lin}}$ is the 2-D matrix of pixel values in the OCT image and ${\sigma }_{\mathrm{lin}}^2$ is the noise variance, both on linear intensity scales. An edge-preserving parameter, $\rho$, was calculated as¹⁵

The value of $\rho$ ranges between 0 and 1. The larger the parameter, the more edges are preserved. $\rho$ should be close to unity for an optimal effect of edge preservation. For the values of SNR, CNR, and ENL, there are no upper bounds. In the situation of OCT speckle reduction, a larger value of SNR, CNR, or ENL indicates a better image quality and better speckle reduction performance.

Table 1 shows image quality metrics for the original and despeckled images with $i$ varying from 2 to 8. In all cases, the SNR, CNR, and ENL showed improvement compared with the original image. The level of decomposition $i$ governs the performance between speckle reduction and edge preservation parameter.

Table 1

Image quality metrics.

The WDC algorithm showed effective speckle reduction on OCT imaging of human fingertip and carotid arteries. As the implementation of the algorithm WDC relies on DWT in MATLAB, there is a maximum wavelet decomposition level for the processed image. Whereas a wavelet decomposition level can be set to any positive integer value in the algorithm, it is still necessary to obtain the maximum level of wavelet decomposition for a specific image as it can help to avoid unreasonable levels of decomposition beyond the maximum one. The maximum level of decomposition is determined by the size of the original image and the wavelet we choose. The rule for computing the maximum decomposition level is the last level for which at least one coefficient is correct.²⁷ The mathematical rationale is that if we want to ensure at least one coefficient is correct, we have to satisfy $(l_w-1)\times 2^{\mathrm{lev}}<l_x$, where $l_w$ is the length of the wavelet filter, lev is the maximum decomposition level, and $l_x$ is the size of the input signal. For the demonstrated image size of IVOCT image of the porcine carotid artery and the chosen wavelet, a maximum level of wavelet decomposition $i_{\max }$ is 8, where $i$ denotes the level of 2-D wavelet decomposition. The level $i$ may vary in order to obtain the best speckle reduction performance with minimum edge preservation degradation. Figures 4(a) and 4(c) show the original carotid artery image and the image after speckle reduction processed by the proposed WDM method at the decomposition level of $i=7$, respectively.

Edge preservation may be compromised in the process of speckle reduction with a small 8.7% loss of edge sharpness as shown in Table 1. More levels of wavelet decomposition do not necessarily produce better results as it is expected that it may unavoidably cause degradation of the image edge sharpness as more levels of interpolation and image fusion process (multiplication) are also involved. There is, therefore, a trade-off between speckle reduction and edge sharpness preservation. The level of wavelet decomposition $i$ acts as an adjustable parameter that controls the trade-off between speckle reduction and edge preservation. Such adjustment is needed, when one wants to see different features such as in the case of evaluating atherosclerotic plaques or neovascularization.¹³ As shown in Table 1, the edge preservation parameter increases at the beginning as $i$ increases, and then drops with further increase of $i$, while the CNR and ENL increase (improve) as $i$ increases, and then head toward a “saturate” state beyond that level. This behavior agrees well with the above discussions. Therefore, a wavelet decomposition level of 7 is believed to be optimal to achieve the best performance between speckle reduction and edge preservation for the demonstrated carotid artery image.

For comparison, the performance of the wavelet-based soft global thresholding speckle reduction technique is also shown in the Table 1. Comparing with wavelet-based soft global thresholding method, the proposed WDC method is superior as it can reduce more speckle noise without significant edge preservation compromise, whereas the speckle reduction based on the soft global thresholding method results in more edge distortion. It is also interesting to point out that optimizing the threshold of soft global thresholding techniques may require highly knowledgeable personnel to perform, as indicated by the previous discussions, while the implementation of the WDC method requires almost no human intervention as the performance controlling parameter is only the decomposition level. It, from another perspective, shows the advantage of the WDC method over conventional wavelet-based methods for its relatively simple comfortable range of trade-off control. More complex wavelet-based despeckling schemes, such as spatially adaptive wavelet filter, may preserve the edge well with a comparable speckle reduction performance; it, however, can be computationally inefficient, taking up to $\sim 7\min$.¹⁴ In contrast, the proposed WDC method provides comparable image quality with significantly less processing time, taking $\sim 2\mathrm{s}$. The computation time increases slightly with the increasing decomposition levels, and the processing time for the maximum decomposition levels is only 2.1 s. Without resorting to complex and computationally inefficient wavelet-based data-adaptive schemes, the presented WDC methods can be expected to be used in many applications. Besides, future implementation in C++ and/or GPU parallel processing may be useful to perform real-time speckle reduction for IVOCT application, in which case, soft global thresholding, however, may not be possible as the speckle pattern fluctuates significantly between B-mode frames due to motion artifacts¹³; and thus the optimal threshold chosen for the soft global thresholding will differ in between B-mode frames.

3.3.

WDC for Clinical Applicability of IVOCT

As a demonstration of clinical applicability, Figs. 6(a) and 6(c) gives OCT images of a stented porcine carotid artery before and after denoising by the proposed WDC method, respectively. Figure 6(b) shows the corresponding despeckled image by wavelet-based soft global thresholding technique for comparison. Parameters (threshold value) were adjusted, so that the ENL metric of Figs. 6(b) and 6(c) were the same. The structure of the vessel wall and the struts of the stent are more observable in Figs. 6(c) than in 6(b). As pointed by arrows in Fig. 6(c), these structures show a clearer visualization after denoising. The despeckled image clearly depicts stent apposition. The degree of apposition of the stent to the vessel wall can also be clearly discerned. This further justifies that the WDC method is robust in the application of intravascular OCT systems.

Fig. 6

(a) Original stented image; (b) despeckled image processed by wavelet-based soft global thresholding method; (c) despeckled image processed by the proposed WDC algorithm. Images are shown on the same gray scale. The red arrows point to structures with clearer discernibility in the processed image. The scale bar represents 1 mm.

3.4.

WDC for Ophthalmic OCT

We also applied the proposed WDC technique to other OCT images and obtained similarly good results. The example is illustrated in Fig. 7 for an in vivo image of human retinal structure obtained by a spectral domain OCT (SD-OCT) system. The SD-OCT system (Vivolight LTD, China) uses a superluminescent diode of 840 nm and allows the capture of 30,000 $\mathrm{A}\text{-}\text{scans}/\mathrm{s}$. The axial and lateral resolutions are 8 and 15 μm, respectively. The sensitivity is 102 dB with 800-μW light incident on the sample. The despeckled image in Fig. 7(c) shows the ability of the proposed WDC method to reduce speckle noise in the ophthalmic OCT images. For comparison, the despeckled image by wavelet-based soft global thresholding technique is also shown in Fig. 7(b). For the same speckle reduction performance [the ENL metric of Figs. 7(b) and 7(c) is the same], we can observe that the proposed WDC method outperforms the wavelet-based soft global thresholding in the edge preservation. As indicated by red arrows, the retinal layer structures in Fig. 7(c) are more discernible than those in Fig. 7(b). This suggests that the WDC is applicable in the ophthalmic OCT systems as well, which justifies the feasibility of the WDC method in the OCT image speckle reduction.

Fig. 7

(a) Original retina image; (b) despeckled image processed by the wavelet-based soft global thresholding method; (c) despeckled image processed by the proposed WDC algorithm. Images are shown on the same gray scale. The red arrows point to structures with clearer visibility in the processed image. The scale bar represents 0.5 mm.