3.1.

Second-Harmonic Generation Microscopy and Spatial Light Interference Microscopy Imaging Systems

The optical setups for the two imaging modalities are illustrated in Fig. 1. While both have been discussed in detail in previous publications,^17,51,56 here, we briefly describe their operating principles.

Fig. 1

(a) Optical setup of the SLIM system, built as a module attached to the output port of a commercial phase contrast microscope. The inset shows an H&E image of the TMA slide, as well as H&E and SLIM images of one of its cores. Scale bar: $200\mu \mathrm{m}$. (b) Optical setup of the SHGM system. The inset shows an H&E image of the TMA slide as well as H&E and SHGM images of one of its cores. C, condenser; Ob, objective; TL, tube lens; IP, image plane; M, mirror; BS, beam splitter; SP, short pass; BP, band pass; SLM, spatial light modulator.

SLIM is built as a module at the output port of a commercial phase contrast microscope (Zeiss Observer Z1). The conjugate image plane at the output port of the microscope is imaged onto a CCD camera using a $4f$ lens system (comprising lenses ${\mathrm{L}}_1$ and ${\mathrm{L}}_2$). At the Fourier plane of the first lens ${\mathrm{L}}_1$, a spatial light modulator controls the phase of the unscattered portion of light with respect to that of the scattered portion in increments of $\pi /2$. The camera acquires four intensity images at four different phase modulations $(0,\pi /2,\pi ,-\pi /2)$ and, from these images, the phase of the imaging field $\varphi (x,y)$ can be computed, as detailed in Ref. 51. Throughout our experiments, we used a $40\times /0.75$ NA phase contrast objective. Using a slide-scanning platform, developed in-house using C++, we imaged the entire TMA slide consisting of 24 cores at $0.5\text{-}\mu \mathrm{m}$ transverse resolution. Subsequent stitching of individual tiles and segmentation of individual core images were performed using algorithms developed in-house.⁴² The acquisition rate for the SLIM system was $0.08\mu \mathrm{s}/\text{pixel}$.

Figure 1(b) illustrates the SHGM system. A Ti:Sapphire laser was used to produce 70-fs pulses at an excitation wavelength of 780 nm and a repetition rate of 80 MHz. The beam was scanned onto the sample and focused by a $20\times /0.8$ NA air-illumination Zeiss condenser. The transmitted beam was collected in the forward direction by a $40\times /0.9$ NA Zeiss objective and then passed through two filters. The first filter was a 680-nm short-pass filter ($680\mathrm{nm}/\mathrm{SP}\text{-}25$) for blocking the laser light, and the second was a 390-nm band-pass filter ($390\mathrm{nm}\pm 18$ to 25 nm) for selecting the second harmonic signal. In contrast to the wide-field CCD detector used in SLIM, SHGM uses a single-point photodetector. More details regarding the system have already been published in Refs. 17 and 56. The acquisition rate for SHGM system was $7\mu \mathrm{s}/\text{pixel}$.

3.2.

Tissue Microarray

The breast TMA used in this study was purchased from US Biomax Inc. (Serial # T088b). The TMA was received from the manufacturer with all human subject information deidentified. Neither the authors of this work nor their institutions were involved in tissue collection. The TMA is comprised of 24 cores from six different cases with eight benign/normal and 16 malignant cases. The malignant cases included cores corresponding to three different stages of breast cancer: IIa, IIb, and IIIa. The insets of both Figs. 1(a) and 1(b) show a hematoxylin and eosin (H&E) stained bright-field image (henceforth referred to as “H&E image”) of the TMA. H&E, SLIM, and SHGM images of one of the cores in the TMA are also shown.

Figure 2 shows H&E, SLIM, and SHG images of three benign and three malignant cores, each at a different disease stage, as indicated. The SLIM images generate uniform contrast across the cores, including in cellular structures. SHGM images, on the other hand, generate contrast only in areas where collagen fibers are present. Note that the SHG signal is more sparse in the malignant cores versus the benign cores due to the higher fraction of epithelial cells in the former, associated with tumor invasion into surrounding stroma. Row 2 of Fig. 2 shows the square of the quantitative phase map $\varphi (x,y)$ (${\mathrm{SLIM}}^2$) for each core.

Fig. 2

Morphological comparison between (a) benign and (b) malignant cores from patients at three different stages of disease as captured by SLIM, SHGM, and H&E images. The square of the quantitative phase map in SLIM is referred to as ${\mathrm{SLIM}}^2$ and is proportional to $[{\chi }^{(1)}{]}^2$.

Figure 3 provides a comparison between morphologies of stromal and epithelial tissue within one core, as revealed in H&E, SHGM, and SLIM images. As shown in Figs. 3(a) and 3(g), both SLIM and SHGM reveal qualitatively similar stromal structures. The situation is different in areas with epithelial tissue, where SLIM images [Fig. 3(c)] contain centrosymmetric structures, such as epithelial cells, which are absent in SHGM images [Fig. 3(i)]. The presence of these cellular structures is confirmed by the H&E stain, which shows cell nuclei in purple color [Fig. 3(l)]. In areas that are a mixture of cells and collagen fibers, quantifying fiber alignment and orientation from SLIM images requires segmentation of epithelial cells, as discussed in detail in Sec. 4.

Fig. 3

(a–c) SLIM, (g–i) SHGM, and (j–l) H&E images of both stromal and epithelial/stromal mixed tissue in a TMA core. (a,d,g,j) Left-most and (c,f,i,l) right-most columns show zoomed-in versions of the regions indicated on the core in the central column (b,e,h,k). As is evident from a comparison with the H&E images, the SLIM images show contrast in both stromal and epithelial regions whereas the SHGM images only generate contrast in stromal regions. (d–f) Segmented SLIM images were obtained from the SLIM images by numerically removing the cellular structures and preserving the collagen fibers.

4.1.

Decoupling Isotropic from Anisotropic Signals in Spatial Light Interference Microscopy Images

Prior to measuring collagen fiber orientation, the SLIM and SHGM images were registered and an image segmentation algorithm was used to remove isotropic structures from SLIM images, leaving behind the highly anisotropic collagen fibers. No information from the SHGM data was used for the segmentation of SLIM images.

To perform image registration, the SLIM images were first downsampled, from $\mathrm{12,000}\times \mathrm{12,000}$ to $9216\times 9216\text{pixels}/\text{core}$, to match the sampling of the SHGM images. These images were then coregistered using ImageJ by selecting control point pairs for each image and matching them using transformation techniques (primarily translation and rotation).

The image segmentation algorithm used for removing isotropic cellular structures from SLIM images is schematically illustrated in Fig. 4. Specifically, for each core image, we first computed the response to the Leung–Malik (LM) filter bank comprising gradient filters at 50 different orientations.^57,58 Each filter in the bank computed the directional image gradient using the first derivative of a Gaussian, oriented at $\alpha =\frac{m\pi }{50}$, where $m=0,1,2,\dots 49$, from the horizontal axis in the image. The resulting stack of responses represented both the magnitude and direction of the local gradient or anisotropy in a pixel’s vicinity for each TMA core. The magnitude of these gradients was then normalized by the highest value in the stack and was summed along the stack (along the $z$-axis in Fig. 4). The resulting image represented the isotropy map of the tissue with higher values representing locally isotropic structures (background and cells) and lower values representing locally anisotropic structures (collagen fibers). This isotropy map was then low-pass filtered (Gaussian kernel of $15\times 15\text{pixels}$, which is slightly larger than one epithelial cell) and the Otsu’s thresholding method was used to find the grayscale level separating anisotropic from isotropic pixels.⁵⁹ The Otsu thresholding method dynamically sets the threshold by minimizing the intraclass variance and maximizing the interclass variance under the assumption of a biomodal class histogram distribution. Since this method establishes a threshold based on the natural separation of isotropic and anisotropic pixels in the data, it can be applied to arbitrary SLIM images without prior training. After setting the isotropic pixels to zero, we employed another round of low pass filtering (Gaussian kernel, $5\times 5\text{pixels}$) followed by thresholding to remove any remaining background pixels and obtain the label map. During this process, the maximum gradient along $z$ for each pixel was also used for detecting the remaining background pixels and marking them for removal. The final segmentation map was then computed using the label map as a mask. The segmentation algorithm was coded in MATLAB, and a prepublished LM filter generator was used for the purpose.⁶⁰

Fig. 4

The segmentation algorithm for removing isotropic structures from SLIM phase images. (a) SLIM phase map before segmentation (color bar in radians). (b) Response to LM filter bank. (c) Normalized isotropy map. (d) Binary label map. (e) Segmented SLIM image (color bar in radians). (f) and (g) Zoomed-in portions of the SLIM and segmented SLIM images, respectively, showing an epithelial stromal boundary. The segmentation algorithm removes isotropic structures and preserves anisotropic signals associated with collagen fibers.

As apparent from comparing the SLIM and segmented SLIM images in Figs. 4(f) and 4(g), the algorithm leads to oversegmentation, and sometimes collagen fibers that are isotropic at length scales equal to or smaller than the size of one epithelial cell are segmented out. In some cases, the segmentation algorithm may break up thin or twisting fibers. As shown in Sec. 4(c), these imperfections in segmentation are subdominant, and the local collagen orientation remains similar between SLIM and SHGM images. Furthermore, as demonstrated in Refs. 4, 20, and 55, biomarkers for prognosis are based on average orientation of collagen over spatial scales that are longer than those at which segmentation errors occur in our images. Thus, for the relevant clinical applications, sufficient information remains in the segmented SLIM images to measure fiber-orientation based biomarkers. Another important consideration for clinical applications is the relative orientation of tumor adjacent collagen fibers with the tumor boundary. While SLIM images can be segmented to digitally remove the epithelial cells, they clearly resolve epithelial structures, allowing determination of the tumor boundary orientation with respect to the collagen fibers.

4.2.

Fourier Analysis

Fourier analysis was carried out on SHGM and segmented SLIM images to extract collagen fiber orientation. Each image was sectioned into subimage regions using grids ($16\times 16$, $32\times 32$, or $64\times 64$), and the localized orientation per subimage (which we call $\theta$) was determined using the Fourier analysis technique outlined in Ref. 55. Regions having a mean direction above a chosen threshold were referred to as anisotropic regions, while those under this threshold were labeled isotropic. To highlight preferred orientation in each subimage, quiver plots were superimposed on anisotropic regions to give an orientation map. Histograms of $\theta$ for anisotropic regions ($\text{number of bins}=32$) and bar plots showing the isotropic and anisotropic region counts were also generated from these data. The $\theta$ histograms were further normalized to obtain orientation probability densities $P_{\mathrm{SHG}}(\theta )$ and $P_{\mathrm{SLIM}}(\theta )$. Figure 5 shows the comparable results obtained for a selected pair of segmented SLIM and SHGM images.

Fig. 5

The Fourier analysis procedure for computing collagen fiber orientation probability densities. (a) SLIM image. (b) Segmented SLIM image. (c) Orientation map of SLIM image. (d) SHGM image. (e) Orientation map of SHGM image. (f) Probability density of collagen fiber orientation $\theta$ computed from the SLIM orientation map in (c). (g) Bar chart showing the number of isotropic and anisotropic regions in segmented SLIM image. (h) Probability density of collagen fiber orientation $\theta$ computed from the SHGM orientation map in (e). (i) Bar chart showing the number of isotropic and anisotropic regions in SHGM image. The Pearson’s correlation $\rho$ between the probability densities for the two modalities is also shown.

4.3.

Comparison Between Second-Harmonic Generation Microscopy and Spatial Light Interference Microscopy Signals

Figure 6 shows the fiber orientation probability densities $P_{\mathrm{SHG}}$ and $P_{\mathrm{SLIM}}$ for three different cores, extracted from their respective SHGM and segmented SLIM images. As shown, the shapes of the density functions obtained from the two modalities are qualitatively similar for each of the cores. To obtain a quantitative measure of this similarity, the following procedure was used. The cross-correlation between $P_{\mathrm{SHG}}$ and $P_{\mathrm{SLIM}}$ was first obtained, and the circular lag corresponding to maximum cross-correlation was computed. The two densities were then shifted relative to one another by this lag, so any errors due to overall rotation between the two images are minimized. This alignment procedure is required for a fair comparison between the two densities because errors in registration of the images from the two modalities can cause one density to be slightly shifted with respect to the other. After alignment, the Pearson’s correlation coefficient $\rho$ between the two densities was computed as a quantitative measure of their similarity. The Pearson’s correlation measures the similarity of any two random variables $X$ and $Y$ and is defined as follows:

Fig. 6

Comparison between the fiber orientation probability densities and bar charts counting number of isotropic and anisotropic regions for three different cores: SHGM (left column) and SLIM (right column). This comparison is shown (a and b) for a benign core, (c and d) malignant stage IIa core and (e and f) malignant stage IIIa core. The similarity of $P_{\mathrm{SHG}}$ and $P_{\mathrm{SLIM}}$ is measured for each case using the Pearson’s correlation coefficient $\rho$, as indicated. Scale bar: $200\mu \mathrm{m}$.

In Eq. (7), the operator $E[]$ refers to the expected value, and $\mu$ and $\sigma$ are the mean and standard deviation, respectively, of the random variable in question.⁶¹ The correlation coefficient $\rho$ has values over the interval $[-\mathrm{1,1}]$ with $-1$ referring to perfect negative correlation and 1 corresponding to perfect positive correlation. In our analysis, $X$ and $Y$ refer to the probability densities $P_{\mathrm{SHG}}$ and $P_{\mathrm{SLIM}}$.

The three cores in Fig. 6 belong to three different disease stages and, therefore, correspond to three different morphologies. As shown, the benign core shows the highest positive correlation between $P_{\mathrm{SHG}}$ and $P_{\mathrm{SLIM}}$ with $\rho =0.93$, which decreases to 0.82 for the stage IIa malignant core. The lowest correlation is seen for the stage IIIa core, which was computed as 0.71. This trend can be accounted for by the fact that a core from a patient at an advanced stage of disease is more likely to contain large amounts of epithelial tissue. A large proportion of epithelial tissue suppresses the SHG signal and results in greater image segmentation errors in SLIM images due to the smaller amount of collagen involved, resulting in lower agreement between the two. In addition to the probability densities, the bar charts showing the number of isotropic and anisotropic cells counted by the Fourier analysis procedure are also very similar for both imaging modalities, as demonstrated in Fig. 6.