Assessment of photoacoustic tomography contrast for breast tissue imaging using 3D correlative virtual histology

Current breast tumor margin detection methods are destructive, time-consuming, and result in significant reoperative rates. Dual-modality photoacoustic tomography (PAT) and ultrasound has the potential to enhance breast margin characterization by providing clinically relevant compositional information with high sensitivity and tissue penetration. However, quantitative methods that rigorously compare volumetric PAT and ultrasound images with gold-standard histology are lacking, thus limiting clinical validation and translation. Here, we present a quantitative multimodality workflow that uses inverted Selective Plane Illumination Microscopy (iSPIM) to facilitate image co-registration between volumetric PAT-ultrasound datasets with histology in human invasive ductal carcinoma breast tissue samples. Our ultrasound-PAT system consisted of a tunable Nd:YAG laser coupled with a 40 MHz central frequency ultrasound transducer. A linear stepper motor was used to acquire volumetric PAT and ultrasound breast biopsy datasets using 1100 nm light to identify hemoglobin-rich regions and 1210 nm light to identify lipid-rich regions. Our iSPIM system used 488 nm and 647 nm laser excitation combined with Eosin and DRAQ5, a cell-permeant nucleic acid binding dye, to produce high-resolution volumetric datasets comparable to histology. Image thresholding was applied to PAT and iSPIM images to extract, quantify, and topologically visualize breast biopsy lipid, stroma, hemoglobin, and nuclei distribution. Our lipid-weighted PAT and iSPIM images suggest that low lipid regions strongly correlate with malignant breast tissue. Hemoglobin-weighted PAT images, however, correlated poorly with cancerous regions determined by histology and interpreted by a board-certified pathologist. Nuclei-weighted iSPIM images revealed similar cellular content in cancerous and non-cancerous tissues, suggesting malignant cell migration from the breast ducts to the surrounding tissues. We demonstrate the utility of our nondestructive, volumetric, region-based quantitative method for comprehensive validation of 3D tomographic imaging methods suitable for bedside tumor margin detection.

www.nature.com/scientificreports/ Cytological examination and frozen sectioning are the current gold standards, but both methods are timeconsuming and have suboptimal accuracy 6,7 . Therefore, the clinical problem is two-fold as positive margins are detected after the procedure leading to relatively high reoperative rates. Additionally, intraoperative histological analysis is a destructive process that limits post-operative histological evaluation. There is a clear clinical need to develop nondestructive intraoperative methods for rapid tumor margin detection to reduce reoperative rates and improve patient outcomes.
In the past decade, numerous reports have described prospective methods to improve intraoperative breast tissue characterization 6,8 . The problem is that a single technique cannot rapidly distinguish tumor margins with high sensitivity and deep tissue penetration. Radio frequency spectroscopy and ultrasound imaging can provide real-time imaging to reduce intraoperative time, but neither provide compositional or molecular information that can enhance cancer tissue detection [9][10][11][12][13] . Near-infrared fluorescence imaging has shown potential in providing in vivo information during breast cancer removal [14][15][16] , but requires exogenous contrast agents that still must be optimized for targeting efficiency, toxicity, and specificity, leading to regulatory approval challenges. Likewise, diffuse optical tomography has also been developed for in vivo breast cancer characterization but lacks the necessary resolution required to assess tumor margins 17,18 . Other optical techniques, such as diffuse reflectance imaging, optical coherence tomography, and Raman spectroscopy can provide high-resolution images to better characterize breast tumor margins, but these methods generally have long image acquisition times and suboptimal penetration depths [19][20][21][22][23][24][25] .
Photoacoustic tomography (PAT) has emerged as a promising label-free method to improve breast margin characterization by providing compositional information with superior penetration depth than conventional optical techniques [26][27][28][29] . The photoacoustic effect relies on interactions between pulsed laser light and tissue chromophores that generate acoustic waves acquired to reconstruct an image. The wavelength of pulsed laser light is tuned to visualize tissue composition, such as hemoglobin and lipid distribution at 1100 nm and 1210 nm light 26,30 . When combined with conventional ultrasound imaging, the user can visualize both breast tissue structure and composition. It is hypothesized that PAT can differentiate non-cancerous lipid-rich and abnormally vascularized hemoglobin-rich cancerous tissue 8,[31][32][33] . Multiple reports support that dual-modality photoacoustic and ultrasound imaging may improve intraoperative breast cancer detection [26][27][28] . More recently, a conflicting study suggested that ex vivo PAT imaging of hemoglobin distribution correlates poorly with ultrasound-determined tumor location 34 . Comparing gold-standard histology, clinical imaging, and photoacoustic data is inherently challenging due to image co-registration issues. Thus, these conflicting PAT reports may be partially due to a lack of robust co-registration methods to accurately compare histology or gold-standard clinical imaging results with photoacoustic imaging datasets at the same location.
Inverted selective plane illumination microscopy (iSPIM) is an advanced optical imaging technique and variant of light-sheet microscopy that can help address image co-registration challenges for the validation of volumetric label-free imaging modalities. iSPIM uses an excitation light sheet to illuminate a plane of tissue, inducing a fluorescence emission collected using a detection objective. This approach allows users to acquire high-resolution volumetric datasets with imaging rates 30 × faster compared to conventional microscopy techniques 35 . Although iSPIM has limited depth penetration for intraoperative use, combined with tissue clearing techniques, it can provide large volume datasets with histological contrast and resolution that can be co-registered to lower resolution label-free techniques [36][37][38][39][40][41][42] . These high-resolution volumetric datasets can provide structural information regarding ducts and blood vessels, lipid distribution, and cellular density. Taken together, iSPIM can help coordinate PAT image co-registration with gold-standard histology, and additionally provide independent volumetric highresolution histological confirmation of sources of contrast present in 3D label-free PAT images.
Here, we present a novel multimodality imaging workflow to quantify the utility of PAT for breast biopsy characterization ( Fig. 1). To do this, we developed a method using high-resolution iSPIM imaging to facilitate image co-registration between PAT-ultrasound datasets and gold-standard histology. We chose rapid, high-resolution PAT and ultrasound as a proposed method for depth-resolved tumor margin characterization due to their ease of integration that provides both tissue structure and label-free composition. Additionally, we use co-registered iSPIM volumes for independent 3D histological validation due to its rapid, large field of view imaging capabilities while being compatible with conventional sample mounting. We take advantage of our nondestructive volumetric imaging datasets and also quantify tissue composition and distribution in multiple regions throughout the breast tissue sample. Finally, we describe the current challenges regarding clinical multimodality breast tissue imaging and potential solutions using emerging engineering approaches.

Methods
Tissue acquisitions and preparation. Six breast cancer tissues from distinct patients were acquired from the Indiana University Simon Cancer Center (IUSCC) Tissue Bank. All six samples were imaged and three heterogenous biopsies were chosen for quantitative analysis. Samples were excised from female patients who underwent lumpectomy procedures and were subsequently diagnosed with infiltrating/invasive ductal carcinoma. Tissues were initially frozen followed by fixation in 4% paraformaldehyde for 24 h. These samples were then submerged in 0.1% paraformaldehyde and stored in a 4 °C refrigerator until imaged.
Photoacoustic tomography. A custom designed PAT system was used to acquire 3D compositional information from the breast cancer biopsies ( Supplementary Fig. S1). This system consisted of a high-resolution small animal ultrasound system (Vevo2100, FUJIFILM VisualSonics) with a 128 element, 40 MHz central-frequency probe coupled with a Nd:YAG optical parametric oscillator laser (Surelite EX, Continuum). Nd:YAG laser delivered a 10 Hz, 5 ns pulse ranging from 670 to 2300 nm to the tissue through a fiber optic bundle. The fiber optic bundle had an opening diameter of 1.0 cm and rectangular terminals of 12 × 2 mm 2 . These rectan- www.nature.com/scientificreports/ gular terminals were attached to the sides of the ultrasound transducer to couple light delivery and ultrasound-PAT signal acquisition. Synchronization of pulse delivery and signal acquisition was performed using a function generator (33220A, Agilent) that delivered a 10 Hz, 10 µs transistor-transistor logic signal to the laser and ultrasound system. PAT and ultrasound co-registration was achieved using a delay generator (DG535, Stanford Research System). For this study, we used 1100 nm light to produce hemoglobin-specific contrast, 1210 nm light to produce lipid-specific contrast, and 1400 nm light as an off-resonance control.
To minimize tissue movement during 3D scans, the breast cancer biopsies were fixed in place using 1% lowmelt agarose in a silicon petri dish. Deionized water was used as an acoustic coupling media to eliminate bubble artifacts that are commonly observed with ultrasound gel. A MATLAB algorithm was then used to control a 3D step-wise motor to translate the PAT probe across the breast cancer biopsy in the y-direction and acquire both PAT and ultrasound images. The 3D motor was moved 0.193 mm across the sample and then paused to acquire 10 PAT images. This process was repeated for the anterior and posterior sides of the breast biopsies. Once the entire sample was imaged, the data was exported and processed using median averaging of 10 images for every location.
Inverted selective plane illumination microscopy. The ASI iSPIM system has been described in our previous publications and is also illustrated in Supplementary Fig. S2 41 . This system was equipped with 488 nm and 647 nm laser sources (Omicron), fiber-coupled laser scanner (ASI) to generate a scanned light sheet, multiimmersion cleared tissue objectives (CTO, ASI/Special Optics) and CMOS cameras (ORCA-Flash 4.0, Hamamatsu). Two immersion objectives were aligned perpendicular to each other, with each having a 45° angle to the horizontal plane 41 . Both objectives switched roles between illumination and detection, making isotropic resolution possible by deconvolving dual view stacks 40 . However, in this experiment a single view was utilized. The field of view was 0.8 mm by 0.8 mm when the refractive index of immersion solution was 1.46. When large samples were imaged, 3D imaging could be realized by using a motorized stage with sub-micron repeatability to reconstruct samples in three dimensions.
After ultrasound-PAT imaging, the breast tissue samples underwent X-CLARITY polymerization and electrophoresis accelerated clearing (XCLARITY, Logos Biosystems) for one week. All cleared samples were stained in 50 µM DRAQ5 solution diluted in phosphate buffered saline overnight (Biostatus, Ltd), and then in 2 mg/ mL Eosin Y solution diluted with 80% ethanol for 30 min (Sigma). After the staining procedure, samples were rinsed three times with deionized water for one minute, X-CLARITY mounting solution for 10 min, and fresh X-CLARITY mounting solution for one additional hour. This process helps match the sample's refractive index to the mounting solution. Samples were fastened on the bottom of the imaging chamber with silicone gel (Dowsil) and then imaged with the iSPIM system. During imaging, the step size between two light sheets was set as 2 μm, overlap between two strips was 15%, and the layer step size was 450 μm. A 488 nm laser source was used for illumination of the Eosin channel and a 647 nm laser source was used for the DRAQ5 channel. Since the iSPIM system is able to image the breast samples cleared by X-CLARITY within a depth of 2 mm, e imaged both the anterior and posterior sides of the breast tissue sample. After imaging, the image stacks of the light sheet strips were reconstructed in a customized MATLAB program 43 , then stitched with the FIJI stitching plugin 44 , converted to pseudo-color images similar to hematoxylin and eosin (H&E) stained slide 45 , and visualized in Amira 40,46 .

Figure 1.
Overview of multimodality workflow used to characterize human invasive ductal carcinoma breast tissues. Qualitative breast tissue analysis was first performed by co-registering PAT, ultrasound, and iSPIM images with histology using overall morphology, internal structures, and composition. Quantitative breast tissue characterization was performed by verifying PAT-ultrasound and iSPIM co-registration using lipid-weighted threshold images, as well as co-registering high-resolution iSPIM images with histology. We then used iSPIM to facilitate co-registration between PAT-ultrasound and histology to assess how tissue composition changes in cancerous and non-cancerous regions. www.nature.com/scientificreports/ Histology. After completion of ultrasound, PAT and iSPIM imaging, the breast cancer biopsies were sent for standard histological sectioning. Tissues were paraffin embedded and cut with a slice thickness of 4 μm for histological sections. We used a step size of 200 μm for qualitative comparison and 30 μm for quantitative analysis. The slides were stained using H&E and then underwent blinded review by a board-certified pathologist to identify cancerous and non-cancerous regions.
Image processing for quantitative analysis. All volumetric datasets underwent image thresholding to assess the percent and spatial composition of the breast biopsies. We first matched H&E histology sections with ultrasound, PAT, and iSPIM datasets based on the overall tissue morphology, internal structures (e.g., ducts), and composition. We then created iSPIM and PAT analysis regions consisting of one imaging slice matched with histology, as well as two imaging slices adjacent to the histology-matched imaging slice. A total of five analysis regions were created for each biopsy, consisting of three images that occupied a biopsy length of 0.386 mm. The PAT images were cropped to the imaging depth of the iSPIM images prior to post-processing. Both datasets initially underwent edge-aware local contrast enhancement to improve morphological contrast and identification of structural edges. We then applied image binarization using Otsu's method to identify PAT-specific lipid and hemoglobin pixels 47 . The iSPIM datasets underwent conventional thresholding using the image histogram to identify a threshold point in each region. The threshold points were averaged between regions and applied to the dataset to discriminate between lipid and stromal tissue, as well as DRAQ5-specific pixels that represent nuclei. Finally, we created tissue composition maps by averaging the binarized biopsy images to topologically visualize tissue composition and calculate percent composition in each region.
We then assessed how breast tissue composition differs in cancerous and non-cancerous regions. To do this, we assigned non-cancerous (presumed healthy) and suspicious cancerous ROIs to our tissue composition maps. ROI size was defined as 1/3 of the 2 mm imaging depth, or 0.66 × 0.66 mm. We first placed 20 ROIs in suspicious cancerous breast tissue determined by a board certified pathologist. We then placed 30 ROIs in non-cancerous breast tissue with lipid, hemoglobin, or nuclei content. ROIs were scored 1-3 using the compositional positive pixels, calculated by summing the number of binarized images that expressed lipid, hemoglobin, or nuclei contrast for a given pixel. ROIs with no hemoglobin or lipid contrast were scored as 1, composition positive pixels between one or two were scored as 2, and composition positive pixels greater than or equal to two composing at least 50% of ROI were scored 3. Nuclei ROI with no cellular content were scored as 1, composition positive pixels between one and two were scored 2, and composition positive pixels greater than or equal to two were scored 3.
Statistical analysis. Linear regression statistics was used to evaluate the correlation between iSPIM and PAT-derived lipid composition, as well as the correlation between iSPIM nuclei and PAT hemoglobin composition. Shapiro-Wilk normality test was performed on all lipid, hemoglobin, and nuclei composition datasets. Datasets that failed normality test underwent log(y) or 1/y transformations. One way analysis of variance (ANOVA) with a Tukey post-hoc test was then performed to determine statistical significance in lipid, nuclei, and hemoglobin composition amongst the biopsies and analysis regions. Kruskal-Wallis test with Dunn's multiple comparisons test was performed on ROI analysis data. The diagnostic ability of lipid, nuclei, and hemoglobin composition was quantified by calculating the area under the curve from receiver operating characteristics (ROC) curves. Statistical significance was considered at p < 0.05. Data are shown as mean ± standard deviation.

Results
Qualitative multimodality image comparison. Ultrasound and PAT images provided morphological and compositional information that was comparable to histology ( Fig. 2a-c). Ultrasound images clearly resolved overall tissue morphology, including specific internal structures such as ducts ( Fig. 2d-f). PAT images allowed for visualization of hemoglobin-rich regions via 1100 nm light ( Fig. 2g-i) and lipid-rich adipose tissue via 1210 nm light ( Fig. 2j-l). We see in sample 1 and 3 that lipid-rich tissues match well between histology and 1210 nm PAT images. Further, in sample 2 and 3 we see the ability of PAT to resolve large blood vessels (Fig. 2b,h), as well as smaller vessels that may feed the ductal carcinoma (Fig. 2c,i). Perfect co-registration in breast sample 2 was difficult between the lipid signal and traditional histology sections due to slight differences in the sectioning/imaging planes and varying step sizes. Adjacent histological sections to breast sample 2, however, reveal greater lipid-rich adipose tissue. Interestingly, breast sample 1 did not produce blood contrast surrounding the ductal carcinoma.
Comparison between iSPIM pseudo-color image and H&E-stained histological slides revealed similar compositional information (Fig. 3a,b). Texture-based volumetric rendering (TVR; Fig. 3a inset) was used to identify iSPIM images that match H&E histology. We were able to clearly identify nuclei structure and density, as well as spatially-resolve breast biopsy tissue composition. Indeed, iSPIM and H&E histology both revealed regions with adipose tissue (Fig. 3c,f), pseudoangiomatous stromal hyperplasia (Fig. 3d,g), and ductal carcinoma in situ (Fig. 3e,h).
Threshold validation and regional compositional analysis. Image thresholding for breast biopsy composition adequately distinguished stroma (Fig. 4a,b), adipose tissue (Fig. 4e,f), nuclei (Fig. 4c,d), and hemoglobin-rich regions (Fig. 4g,h). We qualitatively identified comparable ultrasound, PAT, and iSPIM images using overall biopsy morphology and internal structures (e.g., ducts). A total of 15 regions were quantitatively analyzed, with each region comprising three image slices. This resulted in 45 total post-processed images. Coregistered images first underwent image thresholding to calculate the percent of the biopsy composed of stromal and lipid-rich adipose tissues. Linear regression revealed strong correlation in lipid composition between iSPIM and PAT data (Fig. 4i, R  www.nature.com/scientificreports/ artifacts that were manually segmented and not included in the analysis. Interestingly, there was no correlation between hemoglobin and nuclei content in these biopsies (Fig. 4j, R 2 = 0.03, p = 0.56), suggesting that nuclei dense regions do not necessarily coincide with hemoglobin-rich regions. We quantified breast tissue composition and averaged imaging slices within each region (Fig. 5a) to topologically visualize lipid-rich adipose and fibrotic stromal tissue distribution (Fig. 5b,c), as well as nuclei (Fig. 5d) and hemoglobin distribution (Fig. 5e). Pooled iSPIM and PAT lipid analysis results revealed average lipid composition of 30.6 ± 9.0% in breast sample 4, 10.3 ± 4.2% in breast sample 5, and 13.8 ± 9.1% in breast sample 6 ( Fig. 5f). Hemoglobin contrast also varied greatly amongst different regions and samples, revealing 70.7 ± 5.8% in breast sample 4, 12.9 ± 6.6% in breast sample 5, and 30.0 ± 20.6% in breast sample 6 (Fig. 5g). Both lipid and hemoglobin composition were statistically greater in breast sample 4 than breast sample 5 and 6 (p < 0.05). Average nuclei content comprised 2.1 ± 1.3% of all three breast samples (Fig. 5h), but localization of cellular contrast varied greatly between regions and samples.
Quantitative multimodality imaging and histology comparison. ROI analysis allowed us to compare breast sample lipid, hemoglobin, and nuclei composition in non-cancerous and cancerous regions (Fig. 7a-f). A board certified pathologist identified 20 tissue regions with suspicious cancerous features. Regions    (Fig. 7h). Non-cancerous hemoglobin ROIs scored significantly higher at 2.54 ± 0.79 compared to suspicious cancerous regions at 1.65 ± 0.88 ( Fig. 7h; p < 0.05). Finally, ROI analysis scores were then used to generate ROC curves (Fig. 7i), revealing area under the curve of

Discussion
This study presents a multimodality imaging workflow to characterize breast tissue composition using a volumetric region-based quantitative method. The breast biopsies obtained for this study were classified as invasive ductal carcinoma, where cancer growing in the milk duct begins to infiltrate into surrounding stromal and adipose tissue. We first qualitatively assessed how PAT lipid and hemoglobin-contrast and iSPIM lipid and cellular distribution matched with gold-standard histology in three breast tissue biopsies. Dual-modality PAT and ultrasound allowed us to visualize the breast tissue sample morphology and differentiate stromal tissue from lipid-rich adipose tissue using 1210 nm light (Fig. 2j-l). PAT also provided hemoglobin contrast using 1100 nm light that correlated with regions with arterioles or large arteries (Fig. 2h,i). We chose a single near-infrared wavelength for lipid and hemoglobin imaging to minimize image acquisition time and maximize penetration depth. High-resolution iSPIM provided detailed, histology-like volumetric images that allowed visualization of nuclei-rich areas such as ducts and blood vessels and tissue composition, including adipose tissue (Fig. 3c,f), stroma (Fig. 3d,g), and ductal carcinoma in situ (Fig. 3e,h). We experienced challenges, however, when comparing volumetric datasets with traditional histology sections, due to variable image acquisition planes, histological artifacts, and substantial changes in breast tissue composition between 200 μm histology sections (Fig. 2b,k), supporting the need for iSPIM 3D histological imaging as a replacement for traditional histology sections in this context. In three additional breast tissue biopsies, we created a systematic approach to compare multimodality imaging datasets and quantify lipid, stroma, nuclei, and hemoglobin distribution. PAT and iSPIM imaging slices obtained from the same volumes were first co-registered using overall tissue sample morphology, internal structures, and lipid distribution. High-resolution iSPIM images were also compared to gold-standard histology, allowing us to verify and use iSPIM as a 3D replacement for traditional histological validation. The three breast biopsies used for quantitative analysis were divided and analyzed into five subregions to topologically visualize and quantify Figure 7. Comparison of lipid, hemoglobin, and nuclei composition in cancerous and non-cancerous breast tissues. A board certified pathologist first diagnosed suspicious tissue regions using H&E histology images (a). Ultrasound (b), PAT lipid (c), iSPIM lipid (d), PAT hemoglobin (e), and iSPIM nuclei (f) were used to assign cancerous (black box) and non-cancerous tissue (white box) ROIs. Each ROI was then assigned a score using predetermined rules. ROI analysis suggests that lipid rich regions correlate with non-cancerous breast tissue and stromal-rich regions correlate with cancerous tissues (g). Hemoglobin and nuclei content, however, did not positively correlate with regions with increased breast cancer prevalence (h). ROC area under the curve analysis suggests that lipid and hemoglobin contrast can help distinguish cancerous and non-cancerous region (i). *Denotes p < 0.05 between groups. Error bars represent standard deviation. www.nature.com/scientificreports/ tissue composition throughout the tissue. This approach allows the user to quickly visualize how lipid, stroma, hemoglobin, and nuclei are distributed throughout breast biopsy sub-regions compared to a representative histology image. Image thresholding was performed to extract breast tissue composition in multiple tissue regions. Each thresholding method was tailored to a subset of PAT and iSPIM imaging slices and then applied to the entire breast tissue sample dataset. We found a positive correlation between iSPIM and PAT lipid composition (R 2 = 0.78, p < 0.05), suggesting our thresholding method applied to label-free PAT imaging can be used to segment lipidrich regions. Interestingly, we found no correlation between PAT hemoglobin composition and iSPIM nuclei composition (R 2 = 0.03, p = 0.56), suggesting that iSPIM and PAT can provide unique compositional information to improve tissue characterization. In fact, the regional analysis showed that nuclei content comprises approximately 2.1 ± 1.3% of all three biopsies, while hemoglobin contrast varied from 12.9 ± 6.6% to 70.7 ± 5.8%. These results suggest iSPIM and PAT images can be co-registered using tissue morphology and composition to obtain lipid, stroma, hemoglobin, and nuclei distribution.
After validating our multimodality comparison method, we assessed tissue composition in cancerous and non-cancerous regions. We first assigned ultrasound, PAT, and iSPIM image ROIs in suspicious cancerous regions confirmed by a board-certified pathologist. We then placed ROIs in non-cancerous regions with lipid, hemoglobin, and nuclei contrast. ROIs were placed in the same location on the PAT, ultrasound, and iSPIM image and then scored using a predefined set of rules. Overall, our lipid-rich ROI scoring was not statistically different between iSPIM and PAT images, corroborating that PAT contrast at 1210 nm is lipid-specific and that we are comparing similar regions between volumetric datasets. Our ROI scoring and ROC area under the curve analysis suggest that lipid-rich regions are more likely to express non-cancerous phenotypes and regions with diminished lipid content are more likely to express cancerous phenotypes. This is consistent with histopathology reports that show that malignant cancerous cells migrate from the breast duct to the surrounding stroma, thus transforming ductal carcinoma in situ to invasive ductal carcinoma 48 . Cancer cell migration into adipose tissue may also decrease lipid content in the infiltrated region. Recent reports using PAT and broadband near-infrared mammography study found that reduced lipid concentration was an indicator of malignant breast tissue 33,34 . Therefore, lipid-weighted PAT images may enhance positive breast margin detection and allow clinicians to quantify cancer aggressiveness.
Pathologists routinely assess breast tumor grade using a modified Bloom-Richardson system to semi-quantitatively assess nuclear atypia 49 . Here, we visualized nuclei distribution throughout the breast tissue sample using image segmentation and then compared how nuclei distribution correlates with non-cancerous and cancerous breast tissue regions. Our ROI analysis suggested that nuclei distribution was similar in both non-cancerous and cancerous breast tissue regions. This data may indicate cancerous cells infiltration from the ducts into the surrounding non-cancerous tissues, but healthy control breast biopsies are needed to confirm this hypothesis. ROC area under the curve analysis also suggests limited diagnostic ability of iSPIM nuclei contrast alone. The diagnostic ability of nuclei weighted images may be improved by combining our method with automated segmentation techniques to quantify heterogeneous cancer nuclei characteristics such as enlarged and irregularly shaped nuclei throughout the breast biopsy 50 .
Interestingly, non-cancerous regions also showed greater PAT hemoglobin-weighted contrast compared to cancerous regions. This result contradicts the previously proposed hypothesis that regions with increased PAT hemoglobin contrast may indicate high-vascularized cancerous regions [31][32][33] . Excised breast biopsies lack dynamic blood flow causing blood vessels to collapse and altering blood distribution within the tissue sample. Additionally, hemoglobin distribution is likely altered through the tissue handling process 8,34 , which in this study includes washes after tissue resection, fixation, agarose embedding, and imaging with multiple modalities. Hemoglobin redox state also changes in excised tissues, causing oxygenated and deoxygenated hemoglobin to slowly convert into methemoglobin 51 . Changes in hemoglobin oxygenation and redox state will alter hemoglobin absorption spectra and subsequent PAT contrast 42 . Our results are also supported by Kosik et al., who reported that 690-nm deoxyhemoglobin-weighted PAT images correlated poorly with ultrasound-determined cancer regions in freshly excised biopsies 34 . These factors highlight a need to establish breast biopsy preservation guidelines to preserve hemoglobin distribution and redox state after tissue resection.
The dual-modality PAT-ultrasound imaging used in this study can be further improved for intraoperative breast tumor analysis. We used high-resolution ultrasound to facilitate iSPIM and PAT co-registration; however, ultrasound elastography can be adapted to provide tissue stiffness measurements that may further improve diagnostic accuracy 52-54 . Bae et al. recently used shear-wave elastography to measure mean and maximum stiffness and elasticity ratio in 228 breast lymph nodes obtained from 55 patients 52 . Their results suggest that metastatic breast lymph nodes are significantly stiffer than non-metastatic breast tissue lymph nodes. The elasticity ratio was also higher in metastatic lymph nodes compared to non-metastatic lymph nodes. This pilot study highlights the potential of ex vivo shear wave elastography to enhance intraoperative breast cancer characterization.
PAT image quality is dependent on the light penetration, as photons are more likely to be scattered and absorbed when traveling deep inside tissue. Therefore, heterogeneous light distribution in tissue limits PAT image contrast quantification as superficial tissue absorbs more light and generates greater PAT signal than deep tissue. To our knowledge, PAT transducer frequency and illumination geometry have not been optimized for PAT signal generation in the first 2 mm of the breast biopsy. Previous PAT studies characterize cancerous breast tissue using 2, 12, and 21 MHz center frequency transducers 26,32,55 , while our study used a 40 MHz transducer. Similar to ultrasound, PAT imaging that relies on higher frequency acoustic waves acquisition provides superior resolution, but diminished penetration depth compared to lower frequency acoustic waves acquisition. Additionally, the photoacoustic signal is the product of the physical properties of the tissue known as the Grüneisen coefficient, tissue-specific optical absorption coefficient, and light fluence. The optical absorption coefficients and Grüneisen coefficients have has been characterized for many different issues, but there is still room to optimize light www.nature.com/scientificreports/ illumination geometry to maximize light fluence into the region of interest 56,57 . Therefore, transducer frequency and illumination geometry studies should be performed to obtain the greatest penetration depth and signalto-noise ratio in the first 2 mm of the breast biopsy tissue. Optimized photoacoustic illumination geometries coupled with deep learning algorithms can enhance image contrast in low-fluence regions and sparse data 58,59 . Light penetration also affects iSPIM penetration depth. More effective tissue clearing methods, including Ethyl Cinnamate clearing (ECi), clearing-enhanced 3D (Ce3D), and clear, unobstructed brain/body imaging cocktails (CUBIC), can help reduce light scattering and increase the signal-and-noise ratio in deep tissue regions. These methods have been previously used for human kidney, lung, and prostate tissue 21,60,61 . Enhanced image contrast from optically cleared tissues may also serve as controls to show that the wavelength used during PAT provides contrast from the correct tissue type 62 . Furthermore, the limited iSPIM field of view requires that multiple sub-volumes are stitched to reconstruct 3D datasets. iSPIM stitching introduces image artifacts that limit quantitative analysis. The stitching artifact could be diminished by using stitching algorithms or acquisition techniques such as electronic confocal slit detection or structured illumination, to enhance image contrast in deeper tissue regions 41 . Together, these methods enable more comprehensive supplements to help validate PAT for breast tumor margin detection.
Multimodality approaches inherently require additional manpower to acquire, post-process, and interpret volumetric datasets, which can quickly increase intraoperative margin assessment. Artificial intelligence has the potential to significantly improve multimodality intraoperative workflow by deciding which modalities are needed to assess breast margins, improving image quality, and automating complex volumetric image analysis. Li et al. calculated a theoretical PAT image acquisition rate of approximately 36 cm 2 /min using two wavelengths through a 10 Hz laser, with a 15 mm × 15 mm illumination area and a 200 μm translational step size 26 . On the other hand, iSPIM has a theoretical image acquisition rate of 0.008 cm 2 /min after sample staining and clearing. In the future, it may be possible to use PAT for initial breast margin assessment and then use iSPIM imaging to acquire high-resolution datasets in smaller rapidly-cleared regions where artificial intelligence-determined margin assessment is unclear. Moreover, manual thresholding and image segmentation of large volumetric datasets was labor-intensive in this study and exacerbated by imaging artifacts. Machine learning algorithms can improve image quality by improving image contrast, decreasing background noise, and identifying image artifacts 58,63,64 . Artificial intelligence is actively being developed to analyze mammography and pathological images 65,66 . Pacilè et al. recently performed a multireader, multicase retrospectively study using 14 radiologists and artificial intelligence to assess 240 digital two-dimensional breast cancer mammography images 65 . The authors of this study found that artificial intelligence improved breast cancer diagnostic performance without prolonging the radiologists' workflow. Data obtained from ultrasound, PAT, and iSPIM can be input into deep learning algorithms to differentiate positive and negative breast tissue margins. However, this approach requires quantitative methods to determine where the cancer is located within the biopsy, as presented in this article.

Conclusion
We developed a systematic method that uses iSPIM to rigorously compare volumetric ultrasound and PAT with gold standard histological sections and demonstrate the use of virtual histological imaging as an independent validation tool for label-free volumetric imaging techniques. Diagnostic features, such as lipid and stromal distribution and hemoglobin and nuclei localization, were extracted from multiple breast tissue regions and topologically visualized. Our results suggest that PAT and iSPIM lipid contrast can differentiate cancer-rich stromal tissue and non-cancerous lipid-rich adipose tissue. Nuclei data showed similar cellular content in cancerous and non-cancerous tissue, potentially indicating cancer infiltration from the milk duct to the surrounding tissue. Finally, hemoglobin-weighted PAT images suggest that non-cancerous regions contain more hemoglobin contrast than cancerous regions, conflicting with previously published reports. Together, a multimodality approach has the potential to enhance breast tumor margin characterization compared to a single imaging technique. Our systematic method comparing multiple modalities can also advance artificial intelligence methods aimed at differentiating positive and negative breast tumor margins.