Microscopic assessment of biopsied and resected tissues is central to understanding the underlying pathophysiology and clinical states of many patients. Because most microscopic specimens are translucent to visible light, conventional histological methods achieve image contrast using exogenous labels, such as haematoxylin and eosin (H&E) or immunohistochemical labeling.^1–3 Although these procedures allow pathologists to easily identify cellular processes due to high specificity, they have several limitations, including heavy use of chemicals during sample preparation, which can induce artifacts and distort samples;⁴ long sample preparation time, which can delay pathological diagnosis; and labor-intensive histochemical processes.⁵ Importantly, the staining and visualization of thicker (>5μm) tissues using traditional histology is limited by strong light absorption and scattering. Although recent advances in structured illumination microscopy⁶ and light-sheet microscopy⁷ have alleviated many of those disadvantages, label-free microscopy is still attractive because it might further reduce sample preparation.

A variety of diverse label-free microscopy methods have been employed to circumvent these challenges. Nonlinear microscopy methods are widely used⁸ because chemical specificities are retained without exogenous labels via characteristic vibrational modes^9,10 or nonlinear susceptibilities^11–14 of biomolecules. Nevertheless, the long dwell time to collect weak nonlinear signals impedes high-speed image scanning and rapid identification of pathological regions. Optical coherence tomography also shows promising results in histology,¹⁵ but imaging of subcellular features is limited due to speckle noise and spatial resolution. Quantitative phase imaging (QPI) is an alternative candidate for next-generation label-free histological methods, because it accelerates imaging speed by exploiting holographic imaging and using elastically scattering light information as an intrinsic contrast agent for biological specimens.¹⁶ Nevertheless, conventional QPI techniques provide only two-dimensional (2D) information of thin tissue samples.^17,18

Optical diffraction tomography (ODT) is a three-dimensional (3D) QPI technique that achieves label-free volumetric histopathology. From transmitted-light holographic measurements at various angles, similar to x-ray computed tomography, ODT reconstructs the 3D refractive index (RI) distribution of biological cells.^19,20 Because RI values depend on the number of intracellular biomolecules, including proteins and lipids, ODT allows label-free quantitative 3D morphological mapping of biological specimens and has been widely utilized to advance our understanding of the physiology of various live cells.^21–25 However, several technical issues, such as a small holographic field retrieval field-of-view (FoV) and image degradation due to multiple light scattering, have hindered 3D RI-based histopathological analyses. When using this algorithm, tissue slices with thickness of up 100μm could be reconstructed with high contrast.

Here, we demonstrate multiscale label-free volumetric histopathology of thick tissue slides. We addressed the issues of limited FoV and image degradation by constructing a long-working-distance ODT optical setup and developing a robust multiscale ODT reconstruction and stitching algorithm, which consider optical aberration due to thick tissues. The experimental setup enabled mesoscopic imaging of various pathologic tissues over a millimeter-scale FoV with submicrometer resolution.

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 wavelength=457nm, 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 (NA)=0.6, Olympus] was of 0.2mW/mm2 when using a camera exposure of 2 ms. The scattered field from the sample was collected by the objective lens (UPlanSAPO20X, 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^{27 [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μ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.

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μ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 ∼100μ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.

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μm thickness, deparaffinized using three consecutive 10 min xylene bath, and finally mounted between two number zero coverslips using Permount mounting medium (Permount, RI=1.52, Fisher Chemical). Thin 3-μm-thick sections were mounted on 1-mm-thick slide glass with Acrymount mounting medium.

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^{39 [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.

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μ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.⁴³

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μm, reconstructed volumetric RI tomograms 9μm thick with a voxel size of 170nm×170nm×900nm 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.

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μm and analyzed using the present method. The resulting FoVs were 2mm×1.75mm×0.2mm. 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.

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.

We present multiscale volumetric imaging of unlabeled tissue samples. We imaged millimeter-scale, unstained, 100-μm-thick tissues at a subcellular 3D resolution, which enabled the visualization of individual cells and multicellular tissue architectures in a variety of tissues. By systematically stitching 3D RI tile data obtained with a long-working-distance ODT setup, we demonstrated an FoV of 2mm×1.75mm×0.2mm with a spatial resolution of 170nm×170nm×1400nm, and, throughout this extended imaging volume, image degradation due to sample-induced aberration was significantly reduced using the holographic refocusing strategy. The present method can be understood as extremely high content imaging. The number of optical modes, calculated as the reconstructed volume divided by the size of the point spread function, was ∼20giga voxels.

Considering the label-free volumetric imaging capability exhibited across a variety of tissue samples, this method could potentially be used for rapid cancer diagnosis during intraoperative pathology consultations or small biopsy samples. In addition, the present digital pathology approach has considerable advantages over conventional imaging methods, particularly in terms of measuring seamless multiscale volumetric images; conventional methods only generate one or a few focused plane images. In the current system, the time needed to measure 11×11 sections (1.75mm×1.75mm) and produce volumetric reconstructions is 40 and 135 min, respectively. However, this can be further expedited by implementing a high-speed translational stage, a fast image sensor, and more powerful computing power.

To make the method more applicable, more research is needed on the imaging of fresh tissues and the histological interpretation of RI information. Because RI tomograms do not generate the same data as H&E staining, they would be utilized as complementary information. We envision that this generic approach could have far-reaching applications in histopathology and cytometry, possibly in conjunction with newly emerging machine learning methods for segmentation and cell-type classification.