Superresolving, artifact-free optical coherence tomography with deconvolution-random phase modulation

Deconvolution is widely used to improve resolution in microscopy. Unfortunately, deconvolution in optical coherence tomography (OCT) is sensitive to noise which inherently results in additive artifacts. These artifacts severely impair the �delity of deconvolved images and limit its use in OCT imaging. Here we propose a framework that encodes numerical random phase masks in the Fourier space of a deconvolved OCT image, to produce a sub-resolution image with the artifacts removed. The optimized joint operation of an iterative Richardson-Lucy deconvolution and a numerical synthesis of random phase masks (RPM), termed as Deconv-RPM, enables a 2.7-fold enhancement in both transverse and axial resolutions. We demonstrate that the Deconv-RPM method allows imaging previously unresolved cellular-level details in the human labial mucosa in vivo with signi�cantly improved clarity.


Introduction:
Optical coherence tomography (OCT) is a powerful technique that provides non-invasive images with micrometer resolution and millimeter depth.To assess the tissue microstructures at the cellular and subcellular level, advancements [1][2][3][4][5] have been made to improve spatial resolutions and image contrast.However, high-resolution OCT images often suffer from broaden or distorted point spread function (PSF) due to inherent trade-off between resolution and depth-of-focus, artifacts resulting from strong scattering, de ciencies of optical components, speckle noise caused by a coherent light source, and imaging system's noise.Without introducing more optical components, many software-based approaches have been harnessed to recover the degraded OCT images by deterministic models [6,7] or a class of evaluation metrics [8][9][10].Even that, breaking the resolution limit in high-resolution OCT images remains a challenge.
One promising super-resolution algorithm is deconvolution, which has been extensively used in the uorescence microscopy.Many attempts have also been made in the OCT images for the mitigation of transverse and axial resolution blurring with respect to depth [11][12][13][14][15][16][17].For example, deconvolution has shown the potential for super-resolution in the sparse scatterers [13].However, when amplitude-based data is used to invert a coherent process, the image quality will be inevitably deteriorated as the operation creates sidelobes from the noise in the image around the scatterers.This drawback impedes the application of deconvolution in OCT imaging [18].Here, we describe a post-processing algorithm Deconv-RPM, in which the Fourier or pupil plane of a deconvolved frame is separately multiplied with multiple random phase marks before the products are inverse Fourier transformed to the image plane and summed together.Random phase masks (RPM) are employed to remove the artifacts while preserving the superresolved features with deconvolution.With an analogous idea, numerical random masks methods have also been applied to digital holography for image denoising [19][20][21].
Transverse and axial resolutions are separable in OCT systems, where the axial resolution is governed by the central wavelength and the bandwidth of the detected light while the transverse resolution is limited by the effective numerical aperture.With a separable one-dimensional (1-D) scheme, each depth slice or A-line of an OCT frame is deconvolved by the unregularized Richardson-Lucy algorithm.While some of the regularization constraints such as Tikhonov-Miller regularization [22] and Total Variation regularization [23] have been used to suppress noise ampli cation, they typically call for some prior knowledge on the measured data and have a drawback of either smooth edges or a contrast loss in the image restoration.Our approach was validated experimentally where the transverse and the axial resolution are improved 2.7-fold and 2.7-fold in the confocal region, respectively.In a turbid sample such as human labial mucosa, the unprecedented performance by Deconv-RPM was demonstrated in crosssectional or volumetric imaging compared to the unprocessed images.Deconv-RPM was used to extract scattering signals from nonkeratinized epithelial cells and blood cells in vivo using normal and highspeed line scan cameras, respectively.Finally, we presented an intuitive interpretation of the effects of Deconv-RPM through simulation tests.This robust and reliable method is expected to break new ground for OCT in biomedical applications.

Method (Principle Of Deconv-rpm):
In this section, we elaborate on the Deconv-RPM method that is schematically illustrated in Fig. 1.For coherent imaging, image formation of one B-scan in OCT can be deimages are averaged,scribed as a convolution of the electric eld of the original undistorted sample object u obj (x, z) and the 2-D PSF h(x, z) followed by adding noise n, The backscattered magnitude of the conventional OCT frame u img (x, z) is obtained by the Fourier transform of the raw interferogram in k z -space.As demonstrated, the convolution kernels can be approximated by absolute values instead of the complex values [13].The combined PSF h(x, z)is the product of the individual PSF h(r), can be approximated by, where the pixel index r represents the coordinate (x, z) or y for 3D z-stacks, w r is related to the full width at half-maximum (FWHM) Δrof h(r)according to Δr = 2 √ ln2w r .Shown in the owchart, the Deconv- RPM method mainly comprises two processing steps, where physically breaking the resolution limit with deconvolution and improving the image quality with random phase masks.In step 1, Lucy-Richardson deconvolution algorithm is implemented on one B-scan U img (r) , performed by the built-in MATLAB where o t (r) is the estimated object image after t-th iteration., − u m (r) is the mean value of u m (r), ℑ denotes the 1-D Fourier transform, φ is normal random numbers with a mean of 0. RPM algorithm has two tuning parameters: the number of phase masks m, and the standard deviation of normal random numbersσ. Results:

Demonstration of Deconv-RPM method on sub-resolution particles
To demonstrate the superresolving performance of the Deconv-RPM method, we imaged suspended subresolution particles in water.Detailed system construction of the spectral domain OCT (SD-OCT) has been reported in our earlier work [24,25].In the proof-of-concept study, we rst switched 2 Superluminescent Diodes (SLD1 and SLD3) off (Superlum Broadlighters cBLMD-T-850-HP, Ireland) to provide a light source with a center wavelength of 850 nm and a FWHM bandwidth of 55 nm, resulting in an axial resolution of 5.80 µm in air.The use of a self-made low numerical aperture objective lens allows a large confocal region of 177.84 µm at the expense of the transverse resolution was 9.81 µm in air.The polystyrene calibration spheres (CLB9, Sigma-Aldrich, Diameter 0.9 µm) are ideal as sub-resolution scatterers for characterizing the transverse and axial PSF of the OCT system.
The standard SD-OCT provides low-resolution images of the particles with both the coherent noise and the incoherent noise visible in the background (Fig. 2a).Using xed Gaussian kernels derived from the theoretical PSF, Lucy-Richardson deconvolution improves both the transverse and the axial resolution but also gives rise to ringing artifacts (Fig. 2b) as expected.These artifacts are effectively suppressed by the random phase masks, while the improvement in resolutions by deconvolution is preserved (Fig. 2c).
FWHM and Peak signal-to-noise ratio (PSNR) de ned as PSNR = 10log 10 (MAX 2 (u img )/Var(n)) [26] are used to compare the resolution and signal-to-noise level between the original image and the reconstructed image.In the confocal region, Deconv-RPM achieves 2.71 ± 0.46 times ner lateral resolution and 2.71 ± 0.51 times ner axial resolution than those of the standard SD-OCT (Figs. 2d and   2e).As shown in Fig. 2f, owing to an improvement of 13.4 ± 2.0 dB in PSNR, the image contrast is signi cantly enhanced in Deconv-RPM image.

Super-resolution imaging of nonkeratinized epithelial cells in vivo
To demonstrate the applicability of Deconv-RPM to static scatterers, we recorded the B-scans of human labial mucosa in vivo using µ-OCT (the system used to produce high-resolution images in Supplementary Information, A 30 KHz A-line rate allows for a B-scan rate of 1024 lines at 29.3 frames per second.).The IRB at the Nanyang Technological University (IRB 2016-1-015) approved the human study in vivo, in which written informed consent was obtained from participating subjects.Speckle noise is inherent in coherent imaging [27], which reduces contrast and makes it di cult to resolve boundaries between highly scattered structures in the tissue.Standard B-scan averaging (Spatial compounding) [28] was applied to both unprocessed OCT images and Deconv-RPM images for speckle reduction that the reduction of speckle noise is inversely proportional to √ N (N = 50).However, there is an inherent trade-off between the speckle reduction and the resolution improvement: averaging B-scan comes at the expense of image resolution.It is obvious that the Deconv-RPM image (Fig. 3b) is much sharper than the conventional OCT image (Fig. 3a) and the deconvolved image (Supplementary Fig. 2) throughout the depths.Speci cally, nonkeratinized epithelial cells seen as nucleocytoplasmic boundaries can be unambiguously discerned (Figs.3c and 3d, pink arrowheads in green dashed-line box).Moreover, inside a blurred, noisy vessel in the deep tissue ne structures (pink arrowheads in orange dashed-line box), presumably red blood cells, can be clearly visualized by Deconv-RPM, illustrating the expected resolutions gain and the reduction of noise.
Figures 3g and 3h show the volume-rendered epithelium with both the unprocessed and the reconstructed data.10 frames averaging of the en-face images shows the high scattering signals from the nucleocytoplasmic boundaries.Depths relative to the epithelium in vivo are shown in Supplementary Movie 1.For reference, histology from swine oor of mouth indicates details of anatomical structures such as cytoplasm and nuclei [29].

Resolution-enhanced imaging of moving blood cells imaging in vivo
To demonstrate the potential of Deconv-RPM in moving scatterers, we performed labial mucosa imaging of a human volunteer using µ-OCT (A 126 KHz A-line rate allows for a B-scan rate of 512 lines at 246.1 frames per second.).We scanned the same location on the mucosa multiple times consecutively that moving structures such as red blood cells could be found at each B-scan.Figure 4 presents capillary blood ow containing a mass of blood cells, where the ow direction is marked by red arrow.Figure 4c was created by manually segmenting the capillary region in Fig. 4a and then applying Deconv-RPM to Fig. 4b.For better comparison between images with and without Deconv-RPM, 3X close-up view (Figs.4b' and 4c') on the regions marked in Figs.4b and 4c showing the real-time feature of blood cells.A singlestep discrete Fourier transform (DFT) algorithm [30] was used for alignment of the B-scans that the capillary region was registered.Details of the moving structures for super-resolution imaging in human labial mucosa in vivo are shown in Supplementary Movie 2. The results presented above show that our proposed method is capable of improving the resolution for moving structures using µ-OCT in a noninvasive and real-time manner.

Discussions:
Current deconvolution algorithms applied struggles in the presence of noise such as that in OCT images.
To quantitatively evaluate the performance of Deconv-RPM, we implement 1-D simulations following the processing steps in Fig. 1.The second row in Fig. 5a shows the noisy image u img , generated by the convolution of ideal points (top row of Fig. 5a) and a complex Gaussian PSF (transverse) followed by adding Gaussian white noise n.We measured the FWHM and the PSNR of intensity pro les.The variance of normally distributed noise Var(n) was determined by a reasonable PSNR value of 25 dB, speci ed by the PSNR of the conventional OCT image (Fig. 2a).Sub-optimal w for the initial PSF was set, which amounts to a 1.1-fold increase in optimal w opt .Note, w opt is the FWHM of the absolute values of the PSF.
Tuning parameters m and σwere 10000 and 0.7 respectively.Apparent ringing artifacts can be seen in the unregularized deconvolved image (third row of Fig. 5a) while Deconv_RPM provides artifact-free images with a 2.5-fold improvement in resolution (bottom row of Fig. 5a, Fig. 5b).Image restoration by Deconv-RPM typically took from roughly 70 s for a 1000 pixels × 100 pixels matrix, performing 10 iterations in step 1 and 10000 resamplings in step 2, on a 3.4 GHz Intel core i7-3770 processor with 8 GB memory.
Each resampling is independent that a GPU-based parallel computing design could signi cantly accelerate the computational speed.
We also evaluated the resilience of Deconv-RPM method, clarifying the contributions of the parameters w, σand m to the resolution improvement and the PSNR enhancement.Figure 5c illustrates a slight decrease in resolution improvement when w is away from the optimal w opt .Put in a different way, the PSF estimation has little impact on the performance.By varying the standard deviationσ, we nd that there exists an optimal σ for both resolution improvement and the PSNR enhancement, once other parameters are set (Fig. 5d). Figure 5e shows the ceiling of the algorithm performance when m increases.Under different noise levels of the unprocessed image, shown in Fig. 5f, Deconv-RPM is rather robust to improve resolution, with a PSNR gain of more than 10 dB in the low PSNR scenarios (unprocessed images, PSNR below 32.5 dB).For the experimental results, the FWHM w was set as the estimate of spatial resolution and the number of iterations t for Lucy-Richardson deconvolution was chosen as 4-10; σ and m were chosen empirically to ensure the suppression of artifacts.The pixel size in both spatial and frequency domains should ful ll the Nyquist sampling criterion.
Although our approach is successful in recovering cellular and sub-cellular features of human labial mucosa in vivo, there is room for future development.Some of the weak signals along with the artifacts generated by deconvolution could be also suppressed in step 2 of our approach.Usually, OCT systems encounter limited depth-of-focus when employing a high numerical aperture in the sample arm.Deconv-RPM is compatible with numerical refocusing methods such as interferometric synthetic aperture microscopy [7], aperture synthesis [31] and subaperture correlation [32].In practice, scatterers outside the confocal region could be refocused before the Deconv-RPM processing.Strictly speaking, Deconv-RPM should rely on the phase stability.Since the unstable phase fronts among axial scans do not greatly modify the backscattered intensity, Deconv-RPM is also robust to these errors as Lucy-Richardson deconvolution [13].

Conclusion:
The standard deconvolution, in concert with the numerical random phase masks, permitted artifact-free imaging with an approximate 2.7-fold resolutions improvement, which signi cantly exceeds the performance of these two methods separately.In conclusion, this method could be directly applied to the conventional OCT image without any system modi cation, opens up new perspectives to perform superresolution imaging with potentially drastic noise reduction for nding out the speci c sub-cellular structures in the turbid tissues.More broadly speaking, Deconv-RPM method can work on any other coherent or incoherent imaging modality suffered from severe noise artifacts to render a super resolution image.
See  d-f) was performed on the linear-scale images.Also, see Supplementary Figure 1 for further quanti cation of the performance of the Deconv-RPM on the high-resolution image acquired with micro-OCT (μ-OCT), and its comparison to the deconvolved-only image.
is its transposed matrix.Typically, ∧ h (r) is set as space-invariant in Lucy-Richardson deconvolution.The main challenging tasks for deconvolution methods to superresolve images in OCT are the absence of accurate complex PSF and the inevitable ringing artifacts during the deconvolution.In step 2, the deconvolved image o(r) is further modulated by the phase masks of exp(iφ m ) in its respective Fourier domainO(k r ).Through the forward and backward Fourier operation, a set of modulated images are averaged,

image above for gure legend Figure 2
Figure 2

Figure 3 In
Figure 3

Figure 4 In
Figure 4