2.1.

Sample Preparation

The samples were acquired as 5-mm biopsies during routine pathological preparation of the specimen directly after the surgery. The biopsies were snap frozen in liquid nitrogen and subsequently stored at $-80°\mathrm{C}$. A $5\text{-}\mu \mathrm{m}$ frozen slice was used for H&E preparation according to the standard protocol and the rest of the tissue was used for imaging. Tissue samples were thawed at room temperature and placed in a FluoroDish FD35 Petri dish filled with formalin in such a manner that the imaging was conducted on the consecutive slide to the cut for H&E surface enabling matching of the H&E image with the CARS and/or other measurements. The sample was fixated at the bottom of the Petri dish by putting a transparent block of PDMS on top of the sample.

2.2.

Experimental Setup

The setup consists of a high-power aeroPULSE fiber laser (NKT Photonics)—a passively mode locked fiber laser with a center wavelength of 1032 nm. It has a pulse length of 6 ps and a repetition rate of 80 MHz, with a nominal power output of 10 W. Part of the 1032 power is used as the Stokes beam in the CARS process, whereas the rest is frequency doubled to 516 nm. This light is then used to pump the Levante optical parametric oscillator (OPO) (APE, Germany), which converts part of it into a signal and an idler beam (${\omega }_{\text{pump}}={\omega }_{\text{signal}}+{\omega }_{\text{idler}}$). The idler is not used in this setup and is discarded in a beam dump. The signal beam is used as the pump in the CARS process, and hereafter is referred to as “pump beam.” The wavelength of the pump beam can be tuned from 690 to 990 nm and used in the range of 775 to 807 nm to generate CARS range of 2700 to $3200{\mathrm{cm}}^{-1}$. A small part of this pump beam is split off and directed into a power meter, which is used to monitor the average output power. These data are used later to correct for power fluctuations during the spectral measurements. Another part of the pump beam is split off and coupled into a spectrometer (Ocean Optics HR2000, the Netherlands), which is used to verify that the desired wavelength is generated. If this is not the case, various parameters on the OPO can be adjusted, such as cavity length, Lyot filter position, and crystal temperature. The Stokes beam can be modulated with an acousto-optic modulator (Crystal Technology 3080-197, USA), and is spatially and temporally overlapped with the pump beam through the use of a delay stage and an appropriate dichroic mirror. The overlapped beams proceed into an Olympus FV300 laser scanning microscope which consists of a scanning box with galvo mirrors, and an Olympus IX71 inverted microscope. The microscope is equipped with a filter wheel, with different combinations of filters and dichroic mirrors. Furthermore, a set of objectives (Olympus 60X/1.2NA water immersion and Olympus 20X/0.75NA air objective) are installed allowing to create an image of $235.7\times 235.7\mu {\mathrm{m}}^2/\text{frame}$ or $707\times 707\mu {\mathrm{m}}^2/\text{frame}$, respectively.

For signal detection, a photo multiplier tube (PMT)—Hamamatsu R3896—is placed in the epi-direction. An optional extra (narrowband) filter can be placed directly in front of the PMT (in addition to the filter in the filter wheel) to reject background light from the SHG signal.

2.3.

Homodyne Detection

The TPFE signal partially overlaps spectrally with the CARS signal. The majority of the TPFE can be collected using a spectral filter but the minor overlapping part presents a substantial background to the smaller CARS signal. To separate the contribution of the TPFE from the CARS signal, homodyne detection is used, in which the Stokes beam in the CARS process is modulated (on/off). As a result, the CARS signal, which depends on the presence of the Stokes, is modulated as well whereas the TPFE that is generated predominantly by the pump beam remains constant. By detecting only the modulated part of the signal around the CARS wavelength, we separate the CARS portion. We separately verified that the Stokes beam alone does not generate any significant TPFE. The Stokes beam is modulated by an acoustic optical modulator, with a square wave at 1 MHz generated by an SFG-2110 function generator (GW Instek, Taiwan). The resulting signal from the PMT is demodulated on a HF2LI Lock-in amplifier (Zurich Instruments, Switzerland). A low-pass filter at roughly the same frequency as the pixel sampling frequency is applied, to average the demodulated signal over a single pixel.

2.4.

Data Acquisition

An overview of SHG, TPF, or CARS (2850 or 2940) image acquired with the $60\times$ objective was recovered by stitching the individual frames. Measured with $256\times 256\text{pixels}$ resolution these images feature a spatial resolution of $0.91\mu \mathrm{m}$. For the hyperspectral scan, a single frame was measure repeatedly, while changing the CARS wavelength over the range of 2700 to $3200{\mathrm{cm}}^{-1}$. A combination of a stitching scan and hyperspectral scan (hyperstitch) was made on selected regions of interest.

To accurately reconstruct the Raman spectrum, the measured intensity was calibrated with respect to the output power of the OPO and Stokes beam power, thus obtaining a spectrum that is independent of input power. The data were also corrected for the constant offset on the ADC, of about 50 mV and variation of the PMT voltages were used for signal detection. The stitched images were corrected for the uneven illumination (bright in the center and less bright at the edges) by fitting a fourth order 2-D polynomial to the overall averaged intensity of all frames and dividing the individual frames by this shape.

3.1.

Tissue Fixation

In different studies, it was shown that formalin, as well as liquid nitrogen fixation, can generate some additional peaks in the fingerprint region while altering the intensity of existing peaks which would make the data analysis and results translation to the nonfixed tissue more complex. To examine whether a similar effect happens in the high-wavenumber region, we did measurements using a model tissue—pork. We used pork as fresh, formalin-fixed, and fresh frozen samples imitating regular pathology tissue preparation. For all samples, a CARS spectrum was acquired in a high-wavenumber region 2700 to $3200{\mathrm{cm}}^{-1}$ (see Fig. 1). Our results show that there is no significant difference in spectra acquired from samples treated in different ways. Moreover, our measurements indicated a clear differentiation of different types of pork tissue, namely fat and meat. From Fig. 1(c), it is clear that fat tissue has a dominant peak at $2850{\mathrm{cm}}^{-1}$, corresponding to ${\mathrm{CH}}_2$ stretch, while in the meat tissue Fig. 1(b), the signal from OH stretch (water) is more pronounced. Furthermore, we tested a combined effect of the fixation methods: liquid nitrogen fixed tissue was thawed at room temperature and then fixed in formalin. When CARS spectra were acquired from these samples, no difference was observed when compared to fresh or formalin or liquid nitrogen fixed samples. These results enabled us to develop a protocol for patient tissue handling and use that is independent from the surgery planning. This protocol was implemented as described in M&M. In short, a 5-mm biopsy was acquired directly after surgery and snap frozen in liquid nitrogen. A frozen section H&E was prepared and the consecutive slice was used for the CARS measurements at room temperature in formalin. Such imaging of the consecutive planes enabled matching of the H&E and CARS images (see Fig. 2).

Fig. 1

Spectra for different fixation methods: (a) sample with meat (above) and fat (below), (b) enlarged section of meat, scalebar is $500\mu \mathrm{m}$, (c) enlarged section of fat, scalebar is $500\mu \mathrm{m}$, and (d) CARS spectra for different fixation methods on meat (blue) and fat (green). The error bars indicate standard deviation over the measured region.

Fig. 2

Three examples of matching consecutive tissue slices. The images in the left of each pair are measured with CARS and in the right are stained with HE. The false color in the left uses blue for CARS at $2850{\mathrm{cm}}^{-1}$, red for SHG, and green for TPF. The circled areas are matched. In the first pair of images, the circled area contains healthy cells. In the middle pair of images, the areas circled in red indicated suspicious cells (unpon further inspection labelled as healthy), yellow are healthy cells, and the green areas indicate inflamed cells. The image in the right is all designated as tumor cells.

3.2.

Tissue Spectra

Matching the H&E slides to the microscopic/spectroscopic measurements is not a trivial challenge. Securing consecutiveness of the observation planes and orientation was a great aid; however, tissue deformation could not be eliminated. We created stitched images of the samples using CARS at $2850{\mathrm{cm}}^{-1}$ to highlight the lipids, SHG to highlight the connective tissue, and TPF to highlight cells and connective tissue (i.e., elastin). The composite of these images (a representative example is shown in Fig. 3) was used to find a morphological match between the sample and its corresponding H&E. Due to major tissue deformation, it appeared impossible to make a full sample match; however, the regions of interest for further examination could be selected. For each selected region of interest (ROI), CARS measurements were performed from 2700 to $3200{\mathrm{cm}}^{-1}$. In these data sets, three main tissue groups were defined: fat, connective tissue, and cells; and spectra for each group were reconstructed (see Fig. 4).

Fig. 3

Different signals from the same region: (a) CARS signal at 2980 and $2930{\mathrm{cm}}^{-1}$, respectively, (b) these two overlaid in blue and cyan, respectively, (c) the SHG and TPF signal, and (d) three signals overlaid where SHG is shown in red and TPF in green and the CARS at $2850{\mathrm{cm}}^{-1}$ is shown in blue. The squares are caused by stitching different images together.

Fig. 4

CARS spectra for three tissue types: (a) the combined image, (b) the CARS image at $2850{\mathrm{cm}}^{-1}$, and (c) the average spectra for the blue (fat), green (connective tissue), and red (cells) regions.

Notably, stitched images exhibit uneven illumination, which complicates the data analysis. Several elaborate methods have been proposed over the years to compensate for this nonuniform illumination problem postacquisition.^46,47 However, these are relatively complicated and difficult to implement. In this study, we first subtracted the constant background on the PMT. Then the average intensity distribution over all the individual frames of the images was calculated and smoothed by fitting a fourth-order two-dimensional polynomial to it. The frames were corrected by dividing them by this average intensity distribution, which yields the results shown in Fig. 4. This correction was applied to all collected measurements. Although it might seem more appropriate to measure a homogenous sample and use that as a reference for the intensity distribution, we found that the distribution showed day to day variation. We tentatively assign these to variation in alignment and differences in imaging depth for different samples. The averaging procedure that we use takes no extra time and reduces the stitching errors to an acceptable level.

In total, 16 samples from 14 patients were examined of which 4 were biopsies acquired from the healthy tissue and 12 samples were acquired from cancerous tissue. In total, 396,646 spectra were collected and split into five distinct classes (fat, cancer cells, cancer connective tissue, healthy cells, and health connective tissue). To reduce noise, a block of $4\times 4\text{pixels}$ (an area of $3.7\times 3.7\mu \mathrm{m}$) was averaged to produce one spectrum. Figure 5 presents the averaged CARS spectra collected from all the tissue samples, averaged over areas that were assigned by pathology on the adjoining slice as indicated in Figs. 2 and 4. From these spectra, it is clear that the fat tissue has the most distinct profile. As expected, the $2850{\mathrm{cm}}^{-1}$ peak is most pronounced making this spectrum easily distinguishable from the spectra of other tissue types. When evaluating spectra from cells and connective tissue, the differences are much less pronounced, however, there is a difference in cells’ spectra collected from healthy-tissue biopsy versus cancerous tissue biopsy (see Fig. 5 blue versus purple line). The same is true when comparing the spectra collected from healthy connective tissue versus connective tissue in cancerous tissue. In accordance with the results previously reported in the literature,^39,41 the water/fat ratio (signal at $3150{\mathrm{cm}}^{-1}/\text{signal}$ at $2845{\mathrm{cm}}^{-1}$) in the cancerous tissue is lower than in the healthy tissue.

Fig. 5

Averaged CARS spectra for different regions of interest where the annotation comes from the matched pathology.

Additionally, the intensity distribution of the SHG and TPF signals was examined (see Fig. 6). As expected, the SHG signal for cells has lower intensity than for connective tissue and there is no distinct difference in SHG distribution between cancerous and healthy cells. In the case of connective tissue, there is more signal of lower intensity for the healthy tissue, while cancerous tissue exhibits a more even distribution. For the TPF signal, the connective tissue does not show much difference between cancerous and healthy tissue. In the case of cells, the healthy cells have a peak in the low-intensity region while the cancer cells have a more pronounced peak at the high-intensity region. Despite these distribution differences, the SHG and TPF signals alone do not contain sufficient discriminative power for tissue classification.

Fig. 6

Histograms showing the relative frequency of the SHG and TPF intensities (digitized PMT signals, where 1 $\mathrm{V}=4096$) for the different tissue classes. The horizontal axis is consistent so that the SHG and TPFE signal levels are on the same scale and can be compared. Fat is not included as it does not produce SHG or TPF.

3.3.

Machine Learning Analysis

For a more quantitative analysis, we addressed a machine learning analysis. In this study, we used two methods: (1) recursive feature elimination (RFE)⁵⁰ and (2) feature importance analysis.^51,52 In both cases, we aimed to identify the wavenumbers with the highest discriminative power as well as the minimum number of wavenumbers required for the best result.

3.3.1.

SVM-RFE

Due to the limited number of healthy tissue samples that were measured during this research, there are comparatively few spectra of healthy cells and fat. Thus for fat and healthy cells, SVMs were trained with all available spectra. For the other tissue classes, a random subset of 30,000 spectra for each tissue class was taken for the training of the SVMs.

The SVM-RFE feature selection algorithm repeatedly builds an SVM model to classify the input data and rejects the features with the lowest weights. Following this algorithm, 10 most important wavenumbers were identified (2943, 2892, 3192, 3093, 2994, 2700, 2829, 3042, 2748, and $3144{\mathrm{cm}}^{-1}$). To avoid adjacent wavenumbers, a minimum spectral distance of 40 inverse centimeters was required, which corresponds to the width of the Raman peaks in this region.

As an SVM is inherently a binary classifier, a single SVM can only be trained to separate two classes. For more than two classes, the problem has to be broken down into cascading binary problems. Each binary problem can then be solved by a single SVM, and a majority vote of all the SVMs decides the final classification. For training purposes, SHG and TPF data were also included for some of the tissue classes resulting in a total number of 10 SVM classifiers (Table 1).

Table 1

The 10 SVMs used in the classification scheme, along with their respective input data. Not all available data are given to all SVMs, as this proved detrimental to the scheme’s performance.

SVM	Input class 1	Input class 2	Input data
1	Cancer cells	Cancer connective tissue	SHG + CARS
2	Cancer cells	Fat	TPF + CARS
3	Cancer cells	Healthy cells	TPF + CARS
4	Cancer cells	Healthy connective tissue	SHG + CARS
5	Cancer connective tissue	Fat	SHG + CARS
6	Cancer connective tissue	Healthy cells	SHG + CARS
7	Cancer connective tissue	Healthy connective tissue	SHG + CARS
8	Fat	Healthy cells	TPF + CARS
9	Fat	Healthy connective tissue	SHG + CARS
10	Healthy cells	Healthy connective tissue	SHG + CARS

A sweep was performed over the number of CARS wavenumbers included in the SVM training set. As shown in Fig. 7, the highest Matthews correlation coefficient (MCC) of 0.64 was achieved with a feature set that contains the SHG signal and four CARS wavenumbers. Increasing the number of features did not result in a better performance instead the performance decreases slightly before stabilizing. The optimal feature set exhibits an accuracy of 0.73; sensitivity of 0.75 and specificity of 0.92.

Fig. 7

(a) Influence of the number of parameters (including SHG and TPF) on the classification and (b) row-normalized confusion matrix for a data set with six features (SHG, TPF, and CARS 2892, 2943, 3093, and 3192).

The types of errors that the algorithm is likely to make can be seen from the confusion matrix in Fig. 7. As expected, the performance on the healthy tissue classes is worse than that of cancer classes, likely due to the small number of healthy samples in the dataset. Even though fat has the lowest number of spectra of any class in the dataset, it still had the highest performance. This is owed to the highly distinguishable CARS spectrum of fat, making it easy to identify even when relatively few examples are available to learn from. Moreover, we can see that the algorithm is most likely to make mistakes between two classes of the same tissue type, i.e., it mistakes cancer cells for healthy cells and cancer connective tissue for healthy connective tissue. This means that good separation is achieved between cells and connective tissue as a whole.

To calculate the predicted majority class per sample, the best performing algorithm (the one with SHG, TPF, and CARS 2892, 2943, 3093, and 3192) was used to make a prediction on unseen data. Afterward, the number of pixels that were identified as cancer (either cancer cells or cancer connective tissue) was summed, as are the pixels that were identified as healthy (either healthy cells or healthy connective tissue). The final sample classification is given by whichever class (i.e., cancer or healthy tissue) has the most pixels. Out of the 16 samples in the dataset, 15 were correctly classified. Only one sample in the dataset was misclassified as cancer tissue, while containing tissue that was classified as healthy by pathology. The feature importance analysis exhibited result similar to the SVM-RFE. The best performance was achieved already when four CARS wavenumbers were combined with the SHG data. The accuracy and sensitivity were 0.87, whereas specificity was as high as 0.95; and MCC, 0.87.

3.3.2.

Image reconstruction

To visualize the achieved tissue discrimination results, we applied the trained classifiers to a raw dataset. This visual interpretation of the algorithm’s output can also be compared with H&E annotated slices to check fidelity. As pointed out, the best performance for both machine learning algorithms was achieved with a feature set of SHG, TPF, and four CARS wavenumbers, thus this set was used for classification prior to image reconstruction. The sample data for the reconstructed image were excluded from the training set. Figure 8(a) shows an example of healthy tissue; the connective tissue is generally classified correctly, but the performance on the cells is not as accurate. A lot of pixels were classified as cancer cells. As discussed before, this misclassification can be expected as the algorithm suffers from the shortage of data on normal cells. In the cancer tissue sample [Fig. 8(b)], the cancer cells and cancer connective tissue are classified with good accuracy. Some patches of healthy connective tissue (orange) remains, however, in general the image is conclusive. In both cases, the reconstructed image had clear separation and easily recognizable patterns matching the H&E image.

Fig. 8

Image classification for (a1)–(a3) a healthy and (b1)–(b3) cancer samples. Column (1) shows the H&E stained counter slice, column (2) shows the color coded predicitions. Orange, healthy connective tissue; green, healthy cells; violet, cancer cells; light blue, cancer connective tissue; and yellow, fat. Column (3) shows an overlay of the SHG, TPF, and CARS image.