2.1.

Optical Setup

To image thick tissues, a custom ODT setup was built [Figs. 1(a) and 1(b)]. The system was based on a Mach–Zehnder interferometer, which was equipped with a long-working-distance objective lens for imaging thick tissues, an automated sample stage for raster-scanning, and a digital micromirror device (DMD; DLPLCR6500EVM, Texas Instrument) for illumination beam control. A blue continuous-wave laser (central $\text{wavelength}=457\mathrm{nm}$$\text{wavelength}=457\mathrm{nm}$, Cobolt Twist, Cobolt) was selected to avoid the H&E-staining absorption peak. A plane wave illuminated the samples with a specific illumination angle, which was systematically controlled by projecting time-multiplexed hologram patterns on the DMD.²⁶ The intensity of the light projected on the sample by the condenser objective [LUCPLFLN40X, numerical aperture $\left(\mathrm{NA}\right)=0.6$$(\mathrm{NA})=0.6$, Olympus] was of $0.2{\mathrm{mW}/\mathrm{mm}}^2$$0.2{\mathrm{mW}/\mathrm{mm}}^2$ when using a camera exposure of 2 ms. The scattered field from the sample was collected by the objective lens (UPlanSAPO20X, $\mathrm{NA}=0.75$$\mathrm{NA}=0.75$, Olympus) and projected to the camera (LT425M-WOCG, Lumenera Inc.), where the scattered field from the sample interfered with a reference beam and generated a spatially modulated interferogram, from which both amplitude and phase images, as a function of illumination angles, were retrieved using Fourier analysis²⁷ [Fig. 1(c)]. These multiple 2D field images were used to reconstruct a 3D RI tomogram of the sample using the Fourier diffraction theorem.²⁸ The theoretical resolution of our imaging system, as defined by the maximum bandwidth of the optical transfer function, was 0.17 and $1.4\mu m$$1.4\mu \mathrm{m}$ in the lateral and axial directions, respectively.²⁹

There are two main technical challenges in the reconstruction of wide and thick biological tissues: (1) alignment of multiple overlapping image tiles into a single stitched image and (2) improvement of image contrast against sample-induced aberrations. We addressed these issues computationally using correlation-based subpixel adjustments and digital field refocusing during RI reconstruction, respectively.

2.2.

Evaluation of Optical Diffraction Tomography in Thick Samples

RI reconstruction using the Fourier diffraction theorem²⁸ is based on a weak scattering approximation and, as such, does not hold well as samples become thicker. Recently, different methods mitigating multiple scattering effects on image reconstruction have been developed.^30–32 However, these are computationally expansive and are still not suitable for a wide FoV imaging. However, it has also been shown that refocusing the field at every focal plane before reconstructing the RI can increase contrast in thick samples.^33,34 This method is less efficient than the previous ones, but it is computationally less expansive and has been used in our algorithm.

Furthermore, we only refocused the field at a few focal planes and then reconstructed the RI around those virtual focal planes. This decreases the computational cost while conserving the quality improvement from the original refocusing method. To estimate the efficiency of our algorithm, the same data were reconstructed with and without digital refocusing. Then, the imaging quality was assessed by computing the RI contrast using its standard deviation. Indeed, sample-induced aberrations cause the fields to add up destructively, reducing both contrast and resolution. From Figs. 2(b) and 2(c), one can see that when digital refocusing was not used, the contrast rapidly decreases from the optical focal plane. However, when refocusing is used, contrast again increases at every virtual focal plane. The refocusing interval is chosen as a compromise between the computation speed and contrast lost between the virtual focal planes. In this experiment, a refocusing interval of $9\mu m$$9\mu \mathrm{m}$ was used. Even when using digital refocusing, due to multiple light scattering, the quality of the retrieved tomograms degraded as the tissue thickness was increased. Due to this, the tissue thickness is limited to $\sim 100\mu m$$\sim 100\mu \mathrm{m}$ when using the current setup and sample preparation. The maximum tissue thickness might be lower for fresh tissues, but it could be increased using longer wavelength light³⁵ or ODT techniques using spatiotemporally incoherent illumination.³⁶ The developments of the advanced reconstruction algorithm considering multiple light scattering would also further improve the imaging thickness.^37,38

A last limiting factor is that the illumination and detection size of single tiles is spatially limited. This limits the maximal axial FoV. Indeed, the volume that receives the detection at full numerical aperture and is illuminated by every illumination angle is limited to a conical region above and below the focal plane. Out of this region, the image will have both lower resolution and ghost images due to the circular symmetry of the discrete Fourier transform used in the reconstruction algorithm. It is important to overlap the tiles sufficiently so that these erroneous regions of the tomograms are removed when stitching. In all the experiments, an overlapping distance corresponding to 20% of the camera FoV was used.

2.3.

Sample Preparation

It is crucial when imaging thick samples to use RI matching to limit multiple light scattering. Formalin-fixed paraffin-embedded tissue blocks from the normal small intestine, colon, and pancreas, and from the pancreatic neuroendocrine tumor and intraductal papillary neoplasm of the bile duct were selected. Tissue blocks were sliced with $100\mu m$$100\mu \mathrm{m}$ thickness, deparaffinized using three consecutive 10 min xylene bath, and finally mounted between two number zero coverslips using Permount mounting medium (Permount, $\mathrm{RI}=1.52$$\mathrm{RI}=1.52$, Fisher Chemical). Thin $3\text{–}\mu m$$3\text{–}\mu \mathrm{m}$-thick sections were mounted on 1-mm-thick slide glass with Acrymount mounting medium.

2.4.

Stitching Segmented Tomograms

Accurate image alignment between overlapping image tiles was achieved by correlation-based subpixel adjustment [Figs. 1(d)–1(g)]. During data acquisition, overlapping regions between adjacent 3D tiles were set to 20%, which guaranteed border artifact suppression during RI reconstruction [Fig. 1(d)]. The precise relative 3D position between overlapping regions was then determined with subpixel accuracy using a phase-correlation algorithm³⁹ [Fig. 1(e)]. This position-correcting process enabled artifact reduction due to image misalignment, as evidenced by a mean squared errors analysis [Fig. 1(f)]. Next, the global positions of multiple image tiles were further adjusted using least square minimization.⁴⁰ Minimization was weighted using the Pearson correlation coefficient of overlapping sections. This allowed quantifying similarity between overlapping regions and thus gave less importance to regions without tissue, where finding the relative position was more difficult. Weighting and the use of robust least square minimization avoided error propagation. Finally, images were placed at the obtained position and seamlessly blended together [Fig. 1(g)].

An important consideration when using stitching with ODT is that the illumination and detection size of single tiles is spatially limited. This limits the maximal axial FoV. Indeed, the volume that receives the detection at full numerical aperture and is illuminated by every illumination angle is limited to a conical region above and below the focal plane. Out of this region, the image will have both lower resolution and ghost images due to the circular symmetry of the discrete Fourier transform used in the reconstruction algorithm. It is, as such, important to overlap the tiles sufficiently so that these erroneous regions of the tomograms are removed when stitching. In all the experiments, an overlapping distance corresponding to 20% of the camera FoV was used.

2.5.

Focus Finding in Tissues

The use of coherent light in ODT allows the focal plane to be adjusted numerically, avoiding blurring artifacts due to focus mismatching,⁴¹ which occurs in bright-field microscopy and other incoherent imaging methods. In our reconstruction algorithm, when imaging thin samples, the sample position is first automatically detected, and then an image containing only the most in-focus part of the tissue is generated. The optimal focus plane was found using the Tamura of the gradient focus criterion⁴² on the tomogram, which is a function that reaches its maximum when an image is well focused. The focal position was researched in this way every $40\mu m$$40\mu \mathrm{m}$ and interpolated in regions without samples. This method was also used in thick tissue to quickly find high contrast well-focused lateral cross-sections. Alternatively, the lateral shift between the fields taken from different illumination angles could also be used to find the focus position in a thin sample.⁴³

2.6.

Reconstruction of Refractive Index Distributions

Initially, reconstructed tomograms of thick tissues suffered from weak image contrast and distortion due to sample-induced aberrations. Considering that samples nearest the optical focal plane display sharper image contrast when reconstructed with the Rytov approximation, we exploited holographic refocusing to reconstruct 3D RI maps with higher image contrast^33,34 [Figs. 2(a)–2(c)]. We digitally refocused the obtained light field images to virtual focal planes using a diffraction kernel based on Green’s function [Fig. 2(a)]. The refocusing interval is chosen as a compromise between the computation speed and contrast lost between the virtual focal planes. Here, we set the refocusing interval to $9\mu m$$9\mu \mathrm{m}$, reconstructed volumetric RI tomograms $9\mu m$$9\mu \mathrm{m}$ thick with a voxel size of $170\mathrm{nm}\times 170\mathrm{nm}\times 900\mathrm{nm}$$170\mathrm{nm}\times 170\mathrm{nm}\times 900\mathrm{nm}$ around each virtual focal plane, and combined the tomogram stack into a single tomogram. Qualitative comparison validated that the RI tomograms out of a focal plane could recover image contrast using refocusing-based RI reconstruction [Fig. 2(b)]. The RI distribution root-mean-squares also increased at the out-of-focus planes, quantitatively confirming that digital refocusing improved RI contrast of reconstructed thick tissue tomograms [Fig. 2(c)].

These algorithms were implemented on a graphical processing unit (GPU) to facilitate reconstruction speed. In particular, since phase unwrapping has to be performed at every refocusing step, phase unwrapping is a particularly computationally intensive part of the process. To improve performance, we implemented a custom version of Goldstein’s phase unwrapping algorithm,⁴⁴ where the residue pairing was realized on the CPU, while residue detection, residue linking, rasterization, and phase unwrapping were performed on the GPU.

2.7.

Demonstration in Thick and Thin Tissues

To validate the capacity of the present method for tissue slides volumetric histopathology, we first measured human small intestine tissue samples. For comparison purposes, we prepared two consecutive tissue slides, one using conventional H&E staining and the other unlabeled. The XY cross-sectional images of the unlabeled slide in the optical focal plane were consistent with the conventional in-focus bright-field images of the H&E-stained slide, demonstrating the high accuracy of this method [Figs. 3(a) and 3(b)]. Overall tissue anatomy and subcellular features were clearly seen in both the conventional H&E image and RI tomogram. To demonstrate applicability, various histopathology tissue slides were imaged [Figs. 3(c)–3(e)]. Unlabeled tissue samples from human organs (pancreas, small intestines, and large intestines) were cut at a thickness of $100\mu m$$100\mu \mathrm{m}$ and analyzed using the present method. The resulting FoVs were $2\mathrm{mm}\times 1.75\mathrm{mm}\times 0.2\mathrm{mm}$$2\mathrm{mm}\times 1.75\mathrm{mm}\times 0.2\mathrm{mm}$. The capability of volumetric histopathology can be seen in the high-resolution images of subparts at various axial foci [insets, Figs. 3(c)–3(e)]. Not only can subcellular features be accessed without labeling, but the 3D architecture of tissue structures can also be investigated.

2.8.

Imaging Pathologic Tissue Slides

To validate further the present method for clinical applications, volumetric images of unlabeled thick pancreas tissues samples obtained from patients with pancreas neoplasms are shown in Fig. 4. The RI tomograms of tissues from patient #1 [Fig. 4(a)] show neoplastic epithelial cells with phenotypic neuroendocrine differentiation, characteristics of a pancreatic neuroendocrine tumor (PanNET). The back-to-back neoplastic cells with minimal stroma are easily appreciated, as is the uniformity in the shape of the round nuclei. All of these allow for the recognition of the diagnosis. For validation purposes, adjacent tissues were prepared and imaged using the conventional H&E staining method, which exhibit good agreement with the present method. The RI tomograms of the tissues from patient #2 show low-grade pancreatic intraepithelial neoplasia (PanIN) as well as normal pancreatic ducts, which show different cytoplasmic features between normal pancreatic ducts and low grade PanIN [Figs. 4(b) and 4(c)]. In the liver tissues from patient #3, the normal bile duct, intraductal papillary neoplasm of bile duct (IPNB), and cirrhotic nodule are imaged with the present method [Figs. 4(d)–4(f)]. Thus, different tumor types and a variety of precursor lesions and pathologies can be visualized.