Real-time speckle reduction in optical coherence tomography using the dual window method.

Speckle is an intrinsic noise of interferometric signals which reduces contrast and degrades the quality of optical coherence tomography (OCT) images. Here, we present a frequency compounding speckle reduction technique using the dual window (DW) method. Using the DW method, speckle noise is reduced without the need to acquire multiple frames. A ~25% improvement in the contrast-to-noise ratio (CNR) was achieved using the DW speckle reduction method with only minimal loss (~17%) in axial resolution. We also demonstrate that real-time speckle reduction can be achieved at a B-scan rate of ~21 frames per second using a graphic processing unit (GPU). The DW speckle reduction technique can work on any existing OCT instrument without further system modification or extra components. This makes it applicable both in real-time imaging systems and during post-processing.


Introduction
Optical coherence tomography (OCT) has been proven to be a successful tool for imaging through highly scattering media, such as biological tissues, with high resolution. However, due to the nature of coherent imaging, OCT suffers from speckle noise which degrades image contrast when imaging through turbid samples [1]. Speckle suppression techniques generally depend on compounding several images with uncorrelated speckle noise. These includes spatial [2], polarization [3,4], angular [5] and frequency [1,6,7] compounding techniques. Most speckle reductions techniques require multiple acquisitions at the same location so that uncorrelated speckle patterns can be recorded and eliminated by averaging. Frequency compounding techniques are powerful tools for speckle suppression, but result in images with reduced axial resolution since the total spectral bandwidth is divided into narrow bands to create uncorrelated speckle patterns. The Dual Window (DW) method was originally conceived to perform time-frequency analysis such that the spectral and temporal resolutions are determined independently using two orthogonal spectral windows [8]. Therefore, the DW method can be regarded as a frequency division technique that preserves the original image resolution. In this manuscript, we extend the DW method to perform speckle reduction. Frequency compounding is achieved by processing each A-scan with multiple narrow spectral windows each with independent speckle realizations which are then combined, Significantly, there is minimal loss of axial resolution by using the DW approach where each low-resolution windowed image is multiplied by the full-bandwidth image.

DW speckle reduction method
Spectral windows of two different sizes were used to process each A-scan in the DW speckle reduction (DW-SR) method. Narrow spectral windows (Δλ 1 ) evenly spaced in k-space were used to create multiple images with uncorrelated speckle patterns. These images were averaged to suppress speckle noise. The averaged image was then multiplied by another processed image which used a wide spectral window (Δλ 2 ) so that high axial resolution was retained. The flow chart in Fig. 1(a) is an illustration of the data processing workflow. There are three main concerns when discussing the performance of the DW-SR method. They are preservation of resolution, the choice of spectral window size and number, and the computational performance of the implementation. We will discuss each of the factors in detail in the following three sections.

Preservation of resolution
Since spectral bandwidth is manipulated in the DW-SR process, we carefully characterized the degradation of axial resolution using a commercial spectral-domain OCT system (Spark DRC, Wasatch Photonics Inc., Durham, NC). The OCT system has a center wavelength of 830 nm, bandwidth of 155 nm, imaging depth range of 2 mm, and axial and lateral resolutions in air of 2.1 μm and 8 μm, respectively. The degradation of axial resolution was characterized under two different conditions. First, the scatter-free axial resolution of the system was measured using the full width at half maximum (FWHM) of a mirror reflector. Then, we characterized the resolution degradation in turbid media by imaging the same mirror reflector under ~250 µm of scattering phantom.
The FWHM of the mirror peak was measured to be 2.9 µm when using the full bandwidth of the OCT system. The DW method was used to process the same interferometric data with 30 nm windows, and we found out that the FWHM of the mirror peak was degraded to 4.4 µm. This outperformed a Short Time Fourier Transform (STFT) with the same 30 nm spectral window, which yields a FWHM of 9.5 µm. The axial resolution imaging through the highly scattering sample was characterized by measuring the FWHM of the mirror peak. The scattering phantom was fabricated using polystyrene microspheres (d = 1.09 µm; Thermo Scientific, Duke Scientific, Microgenics Corporation, Fremont, CA) in polydimethylsiloxane (PDMS) (Dow Corning, Midland, MI). The phantom had a reduced scattering coefficient (µ s ') of approximately 2.08 mm −1 (calculated using µ s = 26 mm −1 measured in a transmission experiment and g = 0.92 from Mieplot Software (v4.5)), which is similar to human dermal tissue [9]. The axial resolution was calculated under three different conditions. First, we used the full 155 nm bandwidth of the OCT system, which provided us with the axial resolution of the system. Then, a single narrow window (Δλ 1 = 5 nm, 15 nm, 30 nm or 50 nm) was used to process the interferometric data, which provided us with the axial resolution achieved when a STFT is used. Finally, we used the DW method, in which five narrow spectral windows of different spectral widths (Δλ 1 = 5 nm, 15 nm, 30 nm or 50 nm) were used to produce images with decorrelated speckle noise. To implement the DW processing, the average of the image stack was then multiplied by another image recovered from a single Δλ 2 = 150 nm spectral window. The measured axial resolutions by all three methods are shown in Fig. 1(b). A FWHM of 7.5 µm was measured using the full bandwidth of the system. As expected, for both the DW and STFT methods, a wider spectral window generally provides better axial resolution. However, good axial resolution (< 30 µm) is retained by the DW method even using spectral windows as narrow as 5 nm.

Choice of spectral window size and number
Two competing factors govern our choice of window size. A narrower spectral window causes a greater loss in spectral bandwidth, resulting in worse spatial resolution. In contrast, narrower spectral windows can provide better decorrelation between speckle patterns since more separate windows can be included for a given bandwidth. We therefore must determine the size and number of spectral windows that can provide us with the best balance between speckle reduction capability, spatial resolution and computational complexity.
A scattering Spectralon phantom (Labsphere, North Sutton, NH) was imaged to quantitatively demonstrate the capability of the DW-SR method in speckle reduction. Spectralon samples offer high diffuse reflectance and are highly Lambertian over a broad wavelength range, and therefore can be used as an ideal diffuse sample that exhibits minimal intrinsic variation in signal. The distribution of pixel intensities from a conventional OCT image of the Spectralon sample is shown in blue in Fig. 2(a) and are observed to follow a Rayleigh distribution, which is an intrinsic property of OCT speckle noise [10]. The normalized standard deviation C of the image is calculated using: where σ and <I> are the standard deviation and average intensity of a region of interest in the acquired image, respectively. Note that the calculated C of an OCT image results from the combined contribution of both speckle noise and the intrinsic variation in signal. It has been reported that the intrinsic speckle noise of an OCT signal yields a value of C spk = 0.52 [6]. Therefore, the C OCT = 0.56 measured from the Spectralon sample suggests that there is very limited intrinsic variation in signal, primarily due to speckle noise.
To quantitatively analyze speckle suppression using DW-SR, we first captured multiple OCT images of the Spectralon sample from the same location. A single Spectralon image was processed with different window sizes (5 nm, 15 nm, 30 nm and 50 nm) and window numbers . The spectral windows are evenly spaced in k-space and can be overlapping. Multiple frames were also acquired, and averaged from 1 to 25 times without applying a window. The calculated C OCT and C DW-SR for the Spectralon sample are shown in Fig. 2(b). As expected, simple averaging of multiple OCT images did not provide appreciable speckle noise reduction, yielding only a 4.5% reduction in C OCT with N = 25 averages. The DW-SR method, on the other hand, provides better speckle reduction capability. We achieved ~17% reduction in C DW-SR with five windows for all the window sizes used in the experiment. In general, both narrower windows and a larger window number help in reducing speckle noise. However, narrow windows reduce axial resolution and a large number of windows imposes greater computational demand. For all window sizes, optimal speckle reduction performance was achieved by using ~10 windows. Further increasing the window number did not contribute to the suppression of speckles since not enough spectral separation was achieved, which resulted in greater correlation in the speckle patterns of different images.
The distribution of pixel intensities of the Spectralon image after DW-SR with five 30 nm windows is shown in Fig. 2(a) (in yellow). The distribution has shifted from a Rayleigh distribution towards a Poisson distribution, which is the expected intensity distribution for scatterers randomly distributed in a phantom [10]. The calculated C DW-SR of the DW-SR image also reduced from 0.56 to 0.46, well below the theoretical C spk of 0.52. Fig. 2. (a) Pixel intensity distribution of Spectralon OCT images before (blue) and after (yellow) DW-SR (five 30 nm windows). The pixel intensity exhibits a Rayleigh distribution before DW-SR due to speckle noise, and shifted towards a Poisson distribution after DW-SR due to speckle reduction. (b) Normalized standard deviation C before and after DW-SR, plotted against a theoretical speckle reduction curve (dotted black line) if the speckle noise in the compounding images are completely uncorrelated Finally, we compared the speckle reduction capability of the DW-SR method with the theoretical limit. The normalized standard deviation C of an OCT image can be viewed as the combined contribution from the intrinsic variation in signal (C int ) and the speckle noise (C spk ), which is represented by, The intrinsic variation of signal (C int ) in the Spectralon sample image was calculated to be 0.21 based on the measured C OCT = 0.56 and the theoretical C spk = 0.52. The reduction of speckle noise C spk is inversely proportional to (N)^1/2, where N is the number of compounding images that have completely uncorrelated speckle noise patterns [10]. The theoretical C theo = (C int 2 + C spk 2 )^1 /2 curve is plotted in Fig. 2(b) (black dotted line). Comparing C DW-SR to C theo , we found that increasing N reduces speckle effectively when N is small, particularly for smaller window sizes. However, the gains diminish at larger N, causing the performance of the DW-SR method C DW-SR to diverge from the theoretical decrease of C theo . This is primarily caused by the correlation in speckle between overlapping windows; little benefit is derived from amalgamating highly correlated speckle patterns.

Real-time speckle reduction
Most OCT speckle reduction techniques, although effective, require compounding of several images with uncorrelated speckle patterns during multiple acquisitions. This not only increases the instrumentational complexity of the system to physically generate uncorrelated speckle fields, but also requires longer acquisition times for multiple scans. Here, we demonstrate that the DW-SR method can be used for speckle reduction in real-time OCT imaging. Since additional computational power is needed to implement DW-SR, we utilized the parallel computing power of a graphical processing unit (GPU) to implement DW-SR. The PC system used for performance testing was equipped with an i7 quad-core CPU and a nVidia GeForce GTX Titan X graphics card. All performance analysis was conducted using Matlab (R2017b, Mathworks, Natick, MA). Wavelength interpolation, fast Fourier Transforms (FFTs) and image rendering were performed on the GPU. Each raw interferogram consists of 2048 data points and 512 A-scans were processed for each B-scan.
We first used various window numbers with different spectral widths to test the performance of DW-SR. As shown in Table 1, the processing time increased with the number of windows, but remained relatively constant for different spectral widths. The processing time of DW-SR using five 30 nm windows was approximately 25 ms per frame, which converts to approximately 40 frames per second (fps). This suggests that the DW-SR method can be used for real-time speckle reduction at video rate. The actual frame rate achieved using the above configuration was ~21 fps. By profiling, we found out that our current frame rate is limited by the input/output (I/O) overhead between the camera and the frame grabber rather than the processing time of images. The theoretical frame rate (~40 fps) under current API (Matlab 2017b) may be achieved if a more efficient data transfer scheme is employed. Based on the above analysis, to balance speckle reduction and computation performance, we chose to demonstrate the DW-SR method using five 30 nm windows.

Quantification of contrast enhancement
The improvement in image contrast in both lateral and axial directions by the DW-SR method was quantified by imaging two custom-made phantoms. The first phantom was made of titanium oxide (TiO 2 ) powder (Alfa Aesar, Haverhill, MA) dispersed in PDMS. A small gap (~20 µm) was created in the phantom to form a high-contrast target in the lateral dimension. The second was similar TiO 2 /PDMS phantom, additionally embedded with 25 µm polystyrene microspheres (Thermo Scientific, Duke Scientific, Microgenics Corporation, Fremont, CA), which was used as high-contrast targets in the axial dimension. The contrast enhancement was measured by calculating the CNR represented by: where µ r and σ r are the mean and standard deviation of all region of interest (R), and µ b and σ b are the mean and standard deviation of background noise. All parameters above are derived from logarithmic OCT images (hence the computation of the ratio by subtraction). Figure 3(a) shows images of the ~20 µm gap before and after DW-SR. The lateral crosssections of the images were (averaged along the axial dimension) plotted on the right. A 2.98 dB increase in CNR (7.11 dB and 10.09 dB before and after DW-SR) along the lateral dimension was observed by DW-SR. Improvement in contrast was also observed along the axial dimension using the phantom embedded with 25 µm microspheres. To visually observe the reduction in speckle noise, the intensities of the microsphere images were plotted in 3D on the same linear scale in Fig. 3(b). We observed a significant reduction in the noise floor after DW-SR, which provided a 2.08 dB (9.97 dB and 12.05 dB before and after DW-SR) improvement in CNR along the axial dimension.

Speckle reduction in biological tissue
As a demonstration of speckle reduction in biological tissue, we applied DW-SR to images acquired from healthy mouse skin. Data were acquired in vivo from female C57BL/6 mice (Jackson Laboratories, Bar Harbor, ME). After anesthetization, the mice were shaved and interrogated by the aforementioned commercial OCT system. Representative images before and after DW-SR are shown in Fig. 4. We used five 30 nm spectral windows for DW-SR. Both images are displayed on the same logarithmic scale. Qualitatively, the structural details of mouse skin, especially those from the subcutaneous tissue, appear to be better resolved after DW-SR. The noise floor after DW-SR also appears decreased, since the contribution of shot noise from each spectral window is statistically independent, like speckle, and is similarly suppressed. To quantitatively compare the DW-SR method with existing speckle reduction techniques, we calculated some well-established speckle-reduction performance metrics using the mouse skin images. The signal-to-noise ratio (SNR) was calculated using the equation below: where S lin and σ lin are the magnitude of signal and standard deviation of background noise of the image on a linear scale. The SNR of the DW-SR image is 30.42 dB, compared to 27.94 dB before DW-SR, a 2.48 dB improvement. Another metric used is the contrast-to-noise ratio (CNR) defined in the previous section. The CNRs were calculated to be 4.27 dB and 6.24 dB before and after DW-SR, respectively. This means that a ~25% (1.97 dB) increase in CNR was achieved by using the DW-SR method.

Conclusion
In conclusion, we have presented a new frequency compounding speckle reduction technique using the DW method. This novel DW-SR technique can work on any existing OCT instrument without further system modification or extra components. We explored the choice of spectral window size and number to find the best balance between speckle reduction capability, the loss of resolution and computational speed. With five 30 nm spectral windows, a ~25% improvement in CNR is achieved with only minimal (~17%, from 7.5 µm to 8.8 µm) reduction in resolution. More importantly, we have demonstrated OCT with real-time speckle reduction at video rate, the first such implementation to our knowledge. This was accomplished by exploiting the parallel computing power of a Titan X GPU to assist in both signal processing and image rendering. For comparison, a recent study using speckle-modulating OCT achieved almost 73% reduction in normalized standard deviation C when averaging 100 images with uncorrelated speckle patterns [10]. Another angular compounding technique also achieved 58% (3.99 dB) improvement in CNR by averaging just 2 images with angular separation [5]. Though our improvement of 13-21% reduction in normalized standard deviation C and ~25% (1.97 dB) increase in CNR (using five 30 nm spectral windows) is lower compared to these other studies, these gains are achieved with no physical modification of the OCT system or need for multiple acquisitions. DW-SR can be utilized to reduce speckle both in real-time imaging systems and during post-processing.