1. Introduction

As a biomedical research tool, OCT balances the high resolution offered by optical imaging with greater depth penetration than other optical microscopy techniques [1]. OCT visualizes 2–3 mm below most tissue surfaces, spanning superficial epithelia into underlying layers. However, the resolution and contrast offered by conventional OCT systems and signal processing remain far poorer than that offered by gold-standard histopathology, preventing OCT from realizing its full potential in guiding highly-targeted biopsies for more sensitive detection of tissue pathologies and improved patient guidance. Increasing the resolving power of an OCT system to effectively visualize sub-cellular structures is non-trivial, particularly using hardware that maintains clinical compatibility [2], and finer sampling increases demands on acquisition time and susceptibility to motion artifacts. As an alternative, scattering or attenuation coefficient mapping infers sub-resolution information from the OCT intensity through signal processing [3]. OCT contrast arises predominantly from microscale structures in tissue which interact with incident light as described by light scattering theory. The way these structures scatter light may be summarized through optical properties such as scattering coefficients, which are defined principally by the diameters, refractive indices, and density of a group of scattering particles [4]. In general, the absorption coefficient in tissue at near-infrared wavelengths typically used for OCT imaging (e.g. near 1.1 µm or 1.3 µm) is small compared to the scattering coefficient [5], which is thus the primary contributor to the rate of attenuation of light in the sample. Optical attenuation coefficients are measurable directly from OCT intensity data, and in a well-behaved, homogeneous sample, are proportionate to the slope of the log-linear intensity decay over a given depthwise region of interest (ROI). Initial methods for calculating the attenuation coefficient, $\mu$, from OCT data extracted $\mu$ as a single value over an ROI spanning up to the entirety of the signal penetration depth [3], while more recent work enabled fully depth-resolved $\mu$ mapping for comprehensive 3D visualization throughout volumetric data [6].

Attenuation coefficients have been explored in several biomedical applications as potential sources of contrast to differentiate healthy and pathological tissues [7–9]. One compelling application space is in the diagnosis of dysplastic, or pre-cancerous, lesions. Gold-standard dysplasia diagnoses are made from histological evaluation of biopsied tissues, which assesses the appearance of cell nuclei. With pre-cancerous or cancerous progression, nuclei become larger and denser (e.g., higher nuclear-cytoplasmic ratio). Optical techniques such as light scattering spectroscopy have found that in epithelial tissues, nuclei are dominant scatterers, and that their sizes and densities play a substantial role in determining the overall light scattering behavior of these tissues [10]. Although most OCT systems have insufficient resolution to directly visualize nuclei, this premise may enable attenuation coefficient imaging as an attractive surrogate measurement for reporting the relative sizes and densities of nuclei to differentiate dysplastic or cancerous lesions from healthy surrounding tissue.

A number of challenges complicate the computation of attenuation coefficients that accurately report actual tissue structure. Firstly, in order to isolate the scattering behavior of a sample from the collected OCT intensity signal, it is necessary to fully remove or account for all aspects in which the OCT system itself influences that signal, particularly the confocal function imposed by the system’s focusing optics. Several prior publications have reported methods to mitigate the impact of confocal diffraction from the intensity signal [11,12], but doing so remains a critical and sometimes challenging first step in accurate attenuation coefficient quantification. Secondly, attenuation coefficients fundamentally depend on both the sizes and the concentrations of scattering particles within a region of interest (ROI); using $\mu$ alone, it is not possible to decouple these two independent particle properties. However, it is possible to derive the sample backscattering fraction as an additional optical property from OCT intensity data [13]. Unlike the scattering coefficient, the backscattering fraction, $\alpha$, behaves independently of particle concentration, and therefore, measurement of these two optical properties potentiates the discrimination of both particle diameter and concentration. In our previous work we demonstrated the extraction of a layer-resolved $\alpha$ alongside a depth-resolved $\mu (z)$ [14], potentiating the mapping of layer-resolved particle diameters and depth-resolved volume fractions based on this additional signal.

Even with two scattering properties measured properly, a third major limitation to investigating particle behavior through optical properties is the contribution of multiply scattered light to the OCT intensity signal. Both the coherence and confocal gate reject multiple scattering, but not to a complete extent [15,16]. The detection of photons with path-lengths that do not correspond to actual scatterer location in tissue degrades the accuracy of the application of light scattering theory, which only accounts for singly scattered photons [4]. Although it is acknowledged that multiple scattering events can occur in tissue, and are more likely in samples with larger and denser particles [17], conventional OCT imaging and signal processing do not discriminate singly from multiply scattered light [18]. One theoretical framework for OCT signal formation, the Extended Huygens Fresnel (EHF) model [19], estimates the contribution of multiple scattering to the collected intensity signal by accounting for not only the optical properties of the sample, but also aspects of the OCT system, including the system numerical aperture (NA) and the position of the focal plane (each of which are considered independently from their influence on the confocal function). These system-related parameters are critical in framing multiple scattering as a function of progressive beam distortion with increasing depth into a highly scattering medium, which is itself dependent on the initial beam focusing parameters. Interestingly, their role suggests that modulating system NA or focal plane position could alter the amount of multiple scattering that is detected for a given sample.

In this work, we present our results investigating the influence of modifications to OCT system focusing geometry on the amount of multiple scattering captured in OCT intensity data, as reported by optical properties. Through tissue-mimicking phantoms of known composition, we demonstrate that dilute samples of small particles exhibit consistent behavior regardless of OCT system configuration during imaging, but dense samples of larger particle size exhibit disparate optical properties, particularly when imaged using deeper focal planes and higher system NA. These trends, which were supported by complementary polarimetric measurements and theoretical model analysis [19], potentiate the usage of differential NA imaging to discriminate localized regions of multiple scattering. To assess the impact of strong multiple scattering on optical property interpretation, we additionally propose methods for relating sample $\mu$ and $\alpha$ to the physical diameters and concentrations of constituent spherical particles using Mie theory and known or estimated refractive indices. We demonstrate accurate retrieval of this information in the absence of multiple scattering, and how differential NA measurements can be used to confirm its presence. In murine brain tissue, multiple scattering was identified in blood vessels and a tumor-bearing region, but different regimes of scattering (larger versus denser particles) could nonetheless be decoupled using particle size and concentration measurements even in these areas. Overall, our methods quantitatively connect OCT intensity data to sub-resolution information about the particles giving rise to it, while contextualizing this information in a framework reporting its accuracy. This detailed microstructural assessment has substantial potential to impact various areas of biomedical imaging, including early cancer diagnosis for improved patient guidance.

2. Methods

2.1 OCT system

The custom OCT system used in these experiments was previously described in detail [14]. Briefly, the system incorporates a wavelength-swept Thorlabs VCSEL laser with center wavelength of 1.3 µm, an A-line rate of 100 kHz, and optical k-clocking. The sensitivity of the system was previously measured to be -111 dB with an axial resolution of 10 µm over an imaging range of 11 mm. All visualizations of OCT intensity data have been logarithmically compressed and scaled over a dynamic range of at least 30 dB, beginning at the system noise floor value of 75 dB. Complex fringe data were recorded with polarization-diverse balanced detection for further analysis in MATLAB.

Free space optical elements were added to the system’s sample arm to implement differential NA imaging using two different methods [Fig. 1(a)]. In one configuration, an iris diaphragm was placed to clip the expanded, collimated beam to a manually set diameter. Iris diameter settings used ranged from 12 mm, allowing the collimated beam of 11.3 mm $e^{-2}$ diameter to pass with negligible impact, down to 10 mm, 8 mm, and 6 mm aperture settings. After the iris, the beam was deflected by a galvanometer scanner with 10 mm mirrors (ScannerMAX) for raster scanning over a limited 1D scan range of 4 mm after focusing by the 10$\times$ system objective lens (Thorlabs, LSM02). This configuration was used for all phantom-based experiments.

 

Fig. 1. We used a custom sample arm design (a) to enable two different methods for differential NA imaging. The first method used an iris-based approach (b) to clip the collimated sample arm beam prior to focusing. For each of the four iris diameter settings (corresponding to NA = 0.16, 0.21, 0.26, 0.30), the Rayleigh range and focal plane position were determined experimentally throughout 2-D calibration data (c), where dashed lines represent curves fitted to smooth these results. Applying these confocal parameters resulted in linear intensity decays, as expected (d). Data plotted in panels (c) and (d) correspond to each NA setting given by the outer shading of boxes in panel (b). The second differential NA method used a pathlength-multiplexing element (PME) to simultaneously record angular-diverse signal copies (e). Confocal parameters were also determined experimentally throughout 3-D calibration data (f). Applying these confocal parameters resulted in more uniformly equivalent $\mu (z)$ measured throughout the calibration phantom, as expected (g). Scale bars = 1 mm. M, mirror; L, lens; CL, collimating lens; GS, galvanometer scanner; RM, removable mirror; OL, objective lens.

Download Full Size | PDF

The second sample arm configuration, which was accessible by removing a magnetically-mounted mirror in the free space path, generated depth-multiplexing of differential NA signal copies using a pathlength multiplexing element (PME). The PME design was based on those previously deployed in OCT imaging to achieve depth-of-focus extension [20] or angular diverse acquisition [21] and consisted of a polycarbonate slide of 5.8 mm thickness into which a central hole of 4.55 mm diameter was drilled. After passing through the PME, the initially expanded beam was demagnified to 4.5 mm for coupling with a second galvanometer scanner with 5 mm beam supporting mirrors (Cambridge Technology / Novanta), focused by another 10$\times$ LSM02 scanning lens, and raster-scanned over a larger field of view (8 mm x 8 mm). The resulting tomogram [Fig. 1(e)] featured three signal copies corresponding to different optical path lengths based on passage of the illuminating and detected light through the inner diameter of the aperture (air), the outer polycarbonate portion, or a combination of both as a result of double-pass imaging. The uppermost signal copy, having traveled twice through air, corresponds to the lowest effective NA signal, and the bottommost signal copy, having traveled twice through glass, corresponds to the highest effective NA image copy. These uppermost and bottommost image copies were used for further analysis, and the intermediate copy with combinatorial information was discarded. This configuration was used for tissue experiments, as well as their corresponding phantom calibration measurements.

2.2 Phantom imaging

For ground truth validation of calculated scattering properties and corresponding particle size and concentration values, tissue-mimicking phantoms were constructed consisting of polystyrene beads in solution. The refractive index (1.57 at 1.3 µm wavelength), sizes, and concentrations of these bead-based phantoms were all known within the given range of variation provided by the manufacturer (Polysciences). A solid phantom was made using agar (1% weight by volume) and 0.35 µm diameter polystyrene beads at a volume fraction of 0.012. This phantom was used both for calibration measurements to determine system confocal parameters, as detailed below, and as a predominantly singly scattering reference signal present in the acquisition of all other phantom data.

Additional phantoms were constructed using this solid 0.35 µm bead phantom as a base under a glass coverslip of 200 µm thickness. A 3 µl droplet of larger diameter polystyrene beads at a given volume fraction was placed on top of the coverslip to create another scattering layer [Fig. 2(a)]. For the droplet layer, beads of diameters 0.5, 1.0, 1.5, and 3.0 µm were chosen to span a diverse range of scattering anisotropies, with corresponding $g$ values equal to 0.47, 0.82, 0.89, and 0.94. Droplets were made using a serial dilution to span volume fraction values of 0.0015, 0.003, 0.006, and 0.012, or for 0.5 µm phantoms, 0.003, 0.006, 0.012, and 0.024, due to their lower overall scattering cross section. This series of volume fractions will be denoted as [$V_{f1}, V_{f2}, V_{f3}, V_{f4}$] throughout the remainder of the manuscript for simplification.

 

Fig. 2. Layered, polystyrene bead-based phantoms were constructed as droplets of variable bead size above a glass coverslip and an underlying calibration phantom (a). Intensity (b), $\mu (z)$ (c), and $\alpha$ (d) cross sectional data are shown for layered phantoms at a given focal plane position with NA = 0.21 and $V_f = 0.003$. Scatterplots with $\mu (z)$ (e) and $\alpha$ (f) results averaged throughout white box ROIs show variation in measured optical properties with system NA (along horizontal axis) and focal plane position (individual data points plotted for each NA). Theoretically calculated $\mu$ and $\alpha$ values are provided as black dashed lines, but omitted from 3 µm bead panels due to significant discrepancies with measured data. Error bars indicate the 10th to 90th percentiles of $\mu (z)$ or $\alpha$ values over the shown ROI. Scale bars = 1 mm.

Download Full Size | PDF

Each layered droplet phantom was imaged using each of the four system NA settings accessed via iris-based NA differentiation and were further translated in height with respect to the system objective lens (and thus the focal plane) four times for each NA configuration to span a 1.5 mm distance at regular intervals, with the shallowest focus located 0.5 mm below sample surfaces. In total, the four bead diameters were examined at four different volume fractions, using four different system NA and four different focal plane positions, resulted in 256 unique cross-sectional datasets acquired. Phantoms were imaged using 4 mm of lateral scan range at a lateral sampling rate of 4 µm per pixel. Out-of-plane scanning was not used, but four B-scans were acquired in the same position for each measurement and averaged to reduce the impact of speckle in the liquid droplet.

2.3 Murine brain tumor imaging

A four-week-old mouse underwent cranial window and cerebellar tumor implantation under IACUC protocol 2011N000012. The mouse was euthanized two weeks later by ketamine injection as dictated by the IACUC protocol immediately prior to OCT imaging through the cranial window using the PME-based sample arm configuration for simultaneous differential NA imaging. An 8 mm $\times$ 8 mm field of view was used to image throughout the full 7 mm-diameter cranial window, at a lateral sampling rate of 5.2 µm/pixel and an out-of-plane sampling rate of 62.5 µm/pixel. After imaging experiments concluded, the tumor-bearing brain region was excised, fixed in Formalin, and embedded in paraffin for standard hematoxylin and eosin (H&E) staining. High-resolution digital pathology images were acquired and analyzed using a built-in nuclear segmentation routine in QuPath [22] to create regional maps of nuclei sizes and densities.

2.4 Data analysis

2.4.1 OCT system calibration

A custom calibration routine was used to compensate for the impact of the confocal function on the measured OCT intensity signal, which was unique for each of the four NA and four focus positions at which phantom data were acquired. The solid calibration phantom (0.35 µm beads embedded in agar) was imaged at two different focal plane positions for each NA used in the first system configuration as manually adjusted by the diaphragm position [Fig. 1(b)]. For each scan, the focus depth ($z_f$) and Rayleigh range ($Z_R$) [Fig. 1(c)] were determined on an A-line by A-line basis using methods similar to those previously published by Dwork et al.[11] $Z_R$ is scaled by the refractive index ($n$) of the medium, equal to the group refractive index of water (1.34 at 1.3 µm system wavelength) for both the phantom studies and for the tissue experiments, as $n Z_R$. Since the reference arm position was kept constant throughout all phantom data acquisition, the absolute focus position [Fig. 1(c)] was determined by accounting for the location of the sample surface and the change in optical distance of the focus for each measurement, accounting for sample refractive index. Determined Rayleigh ranges were consistent with expectations (see Fig. S1, supplemental document), and these $Z_R$ and $z_f$ values were used for confocal function compensation in phantom data [Fig. 1(d)], using the absolute $z_f$ converted back to the optical depth of $z_f$ for each individual scan (see Fig. S2, supplemental document).

A similar routine was used to determine effective $Z_R$ and $z_f$ values for removing the confocal function throughout volumetric aperture-multiplexed data [Fig. 1(e)]. The multi-height curve fitting method [11] was used to extract confocal parameters for the uppermost, lowest-NA signal copy, and the bottom, highest-NA signal copy in another 0.35 µm calibration phantom. Confocal parameters were determined for each A-line of the volume over which tissue experiments were conducted [Fig. 1(f)]. Following $H(z)$ removal from the phantom tomograms, extracted $\mu$ values were more consistent both throughout a given scan and across differential NA acquisition [Fig. 1(g)].

For backscattering fraction calculations, the absolute incident light intensity, $L_0$, of the OCT system must be known under given acquisition conditions. Due to the methods used to vary system NA, which involved clipping or segmenting the incident beam, it was necessary to account for a different $L_0$ in each of these cases. Straightforward calibration based on the overall collected light intensity is inadequate, since $\alpha$ itself also changes with NA independently of its relationship to $L_0$. Therefore, our procedure for determining each $L_0$ value used the 0.35 µm calibration phantom data that was acquired beside the layered droplet phantom region as a well-behaved sample with presumed single scattering behavior that could reliably report $\alpha$ as a function of $L_0$ and known NA modulation. The factors used to achieve an $\alpha$ in the calibration region that varied according to expectation with NA under all acquisition conditions were applied to the rest of the tomogram to globally account for $L_0$ differences (see Fig. S3, supplemental document). For the aperture-multiplexed configuration, the 0.35 µm phantom imaged prior to tissue data acquisition was used analogously to compensate for the difference in $L_0$ based on aperture-based differential NA acquisition.

Fourier optics simulations for the resulting confocal function following aperture-based beam segmentation, as well as iris-based clipping, are provided in Figure S1 of the supplemental document. These simulations confirm that despite the deviation of input beam profiles from a Gaussian shape, their confocal functions may still reasonably approximated as Gaussian and described using previously formulated expressions [12].

2.4.2 OCT data processing

In phantom experiments, four B-scans were acquired without out-of-plane scanning and averaged to decrease speckle contrast in the liquid droplet portion. Using parameters determined for each system NA, noise floor subtraction and confocal correction were applied prior to additional median filtering for further speckle removal. The median filter used axial and lateral kernel sizes of 5 and 15 pixels, corresponding to 27 and 58 µm, respectively. In tissue experiments, depth-multiplexed differential NA data acquired with the aperture-based system configuration were processed to extract the surface of the highest and lowest image copies and isolate each of these signals. Whereas median filtering was sufficient to mitigate the impact of speckle in homogeneous phantom data, more sophisticated depeckling was required to preserve structural information in tissue samples. Accordingly, the higher and lower NA intensity tomograms were independently despeckled computationally with TNode [23] prior to noise floor subtraction and confocal correction based on calibrated confocal parameters for this setup as described in Section 2.4.1. As previously reported [14], the OCT system used for these experiments has minimal sensitivity roll-off (<1.5 dB) over the full depth range (11 mm). Further, since each PME-multiplexed signal copy is independently calibrated for its respective $L_0$, the amount of roll-off to be considered within the depth range of each signal copy (to which axial signal analysis is confined) is less than 0.25 dB. Therefore, it was not necessary to account for this additional system effect.

For both phantom and tissue data, we applied our previously published layer-based, depth-resolved algorithm to calculate depth-resolved attenuation coefficients and layer-resolved backscattering fractions [14]. Briefly, this approach automatically defines depth-wise regions characterized by different $\alpha$ or $\mu$ and analyzes the scattering properties over these effective layers independently, which prevents depth-varying $\alpha$ from impacting $\mu (z)$ calculation, a limitation in the previously developed depth-resolved model [6]. For $\alpha$ calculation, which is based on the scattering rather than attenuation coefficient, we determined $\mu _s = \mu - \mu _a$, where a $\mu _a$ of 0.11 mm$^{-1}$, the absorption coefficient of water at the 1.3 µm central wavelength, was used [24]. For summarization of results in scatter plots, data were averaged over a superficial ROI (shown as white box in figures) with error bars indicating the 10th and the 90th percentile values of the data distribution.

As a complementary metric for sample characterization, the degree of polarization (DOP) was also mapped throughout phantom data. DOP for the single input polarization state was calculated from the system’s polarization-diverse detection channels as

2.5 Scattering theory calculations

2.5.1 Forward calculations

Mie theory was used to calculate expected scattering coefficient and backscattering fraction values for particles in tissue-mimicking phantoms of known composition [4]. To estimate theoretical values for attenuation coefficients in these phantoms, $\mu _s$ values were augmented by the absorption coefficient of water at the OCT system center wavelength. Unlike $\mu _s$, $\alpha$ is directly dependent on OCT system parameters, as the system NA determines the angular extent over which backscattered light is captured with respect to all light scattered by the particle. The backscattering fraction is defined as

The Extended Huygens Fresnel model was also used to provide a theoretical context for assessment of the contribution of multiple scattering to intensity signals [19]. Intensity decays incorporating the EHF model were calculated to mimic phantom data using known particle sizes and densities to find Mie-derived $\mu _s$ and root-mean-square scattering angle, $\theta _{RMS}$, as previously described [13]. These values, as well as OCT system NA and focal plane location, were used to find local beam waist in the presence and absence of forward scattering throughout sample depth, and determine OCT intensity signals including multiply scattered light contribution (and excluding impact of the confocal function and noise floor) as formulated by Thrane et al. [19]

2.5.2 Inverse calculations

To better demonstrate how accuracy in measurement of scattering properties impacts their ability to relay physical particle information, we developed a framework for translating these scattering properties into actual particle diameters and concentrations. This framework relies on assumptions including known or accurately estimatable refractive indices of scattering particles. It is also limited to the applicability of Mie theory to spherical or ellipsoidal particles, and cannot be accurately applied to cylindrically-scattering particles, such as fibrous tissues or extracellular matrix (see Discussion). For a given refractive index, Eq. (2) was used to calculate $\alpha$ values for a range of scattering particle sizes as measured by each system NA, generating NA-specific curves of $\alpha$ versus particle diameter over the range of expected particle diameters, 0–3 µm [Fig. 5(c)]. These curves were smoothed using a median filter with kernel size based on the variability in bead size range (as stated by the manufacturer through a coefficient of variation) and the non-zero bandwidth of the OCT system (see Fig. S4, supplemental document). These curves were then used to map layer-resolved $\alpha$ values in cross-sectional phantom data to scattering particle diameter, $d$.

From layer-resolved $\alpha$-derived particle diameters ($d$) and depth-resolved $\mu (z)$ maps, it was also possible to map particle volume fraction by making further use of Mie theory. As an intermediate step, particle diameter values were used to find scattering cross section values, $\sigma _s$, which are calculated using $d$, illumination wavelength, and the refractive indices of the particle and medium (each of which were also used to map $\alpha$ to $d$). Particle volume fraction, $v_f(z)$, is then calculable as

This framework was also applied to extract scattering particle diameter and concentration from tissue data. $\alpha$ values were calculated over a larger range of potential diameters, 0–10 µm, for a range of previously estimated refractive index values for scattering particles in tissue epithelia (e.g., nuclei), 1.36–1.42 [26–28] (against a presumed medium refractive index equivalent to the group refractive index of water, 1.34, at 1.3 µm wavelength), over the appropriate NA values determined for the aperture-based angular-diverse signal multiplexing. Based on an $\alpha$-$d$ curve selected for a given particle refractive index value (see Figs. S4 and S7, supplemental document), the procedures described above to extract $d$ and $v_f(z)$ in phantoms were then also used to assess the particulate makeup of healthy and neoplastic tissue.

3. Results

3.1 Multiple scattering behavior in phantoms is modulated by system configuration

Bead-based phantoms were used to assess the reliance of scattering properties on OCT system NA and focal plane position. Layered, droplet-based phantoms [Fig. 2(a)] were constructed as described in Methods and imaged under variable sample conditions, including modulation of system NA and focal plane position. Intensity tomograms were corrected for the impact of the confocal function [Fig. 2(b)] and filtered (as described in Methods) prior to computation of depth-resolved attenuation coefficients [Fig. 2(c)] and layer-resolved backscattering fractions [Fig. 2(d)] using our previously published methods [14]. In Fig. 2, representative cross-sectional data is shown for each superficial bead size evaluated, at a given volume fraction (0.006), NA (0.21), and focal plane position. Attenuation coefficient and backscattering fraction values were averaged throughout the superficial bead ROI shown by the white boxes and compared across each set of OCT system imaging conditions for each bead at each volume fraction [Fig. 2(e, f)]. In general, it was evident that measurements of $\mu (z)$ in phantoms with smaller, less dense particle makeup were better described by the theoretical predictions of Mie theory, but also more consistent with variation in system NA and focal plane position. In phantoms containing larger and denser particles, the inconsistency of $\mu (z)$ measurements under different system NA and $z_f$ positions suggested a deficiency in the ability of the confocal function to accurately compensate beam behavior in the sample, as expected with multiple scattering. Apart from the inapplicability of standard methods for confocal function correction to multiply scattering samples, in general, the coupling of multiply-scattered light into the system increases the intensity of the tomogram compared to the expected intensity in a single-scattering regime: for this reason, measured $\mu$ values are lower than those predicted [19]. This behavior is clear in 1.5 µm beads at $V_{f4}$ and 3.0 µm beads at $V_{f3}$ and $V_{f4}$. In backscattering fraction data, $\alpha$ is expected to increase as a function of system NA in a single-scattering regime [13], but calculated $\alpha$ values were increasingly higher than expected with greater particle volume fractions, another indication of multiple scattering. In dense ($V_{f4}$) samples, we visualized the impact of increasing NA on the detected contribution of multiple scattering, as indicated throughout data as a reduction in measured particle $\mu (z)$, by comparing the average $\mu (z)$ measured across all focal plane positions for each system NA to the same metric captured by the lowest system NA (see Fig. S5, supplemental document). Particularly in 1.5 µm and 3.0 µm bead samples, we observed a progressive reduction of this ratiometric value with increasing NA.

3.2 Comparison of multiple scattering behavior in phantoms to a theoretical model

The observations above indicate that detected multiple scattering can increase with increasing system NA, a striking observation that does not align with the behavior of decorrelation tails observed in OCT angiography, which are known to arise from multiple scattering in blood vessels [29,30]. In particular, it has been well-documented that higher NA optics limit the size of decorrelation tails in the retina that appear below large blood vessels, pointing to a reduction in detected multiply-scattered light. To more closely examine these data, we compared individual tomograms of our layered phantoms with smaller [0.5 µm, Fig. 3(a)] and larger [3.0 µm, Fig. 3(b)] beads at relatively high concentrations ($V_f = 0.006, V_f = 0.012$). In these layered samples, it should be possible to visualize the border between the bottom of the droplet and the glass coverslip (see Fig. 2), except for in the case of significant multiple scattering, wherein this border is blurred or lost entirely. This effect has been previously demonstrated in multiply scattering samples, where photons having undergone forward- and side-scattering events appear below their physical location in the tomogram [15]. By visualizing the intensity tomograms without confocal function compensation, it was evident that higher NA acquisition (NA = 0.3) did not necessarily prevent the OCT system from capturing multiple scattering, and even appeared to record greater relative false intensity (e.g., tail artifacts) as compared to the lower NA (NA = 0.21) data. In high NA (0.3) tomograms, blurring at the sample surface is also visible in both axial and lateral dimensions. These strong aberrations suggest a regime in which the confocal gate itself may be disrupted. We note that due to the beam-clipping methods used to reduce system NA, there is also a reduction in incident light intensity ($L_0$) with lower NA, but in visualizing individual A-lines of normalized intensity [Fig. 3(c)], $L_0$ differences alone do not appear to account for the relatively higher amount of tail intensity in the coverslip region following the droplet border (shaded in gray) for which, in some lower concentrated cases, the signal can be seen to drop rapidly down to the noise floor as expected.

 

Fig. 3. Higher (0.012) and lower (0.006) volume fraction cross-sectional intensity data are shown at higher (NA = 0.3) and lower (NA = 0.21) NA for 0.5 µm beads (a) and 3.0 µm beads (b) where contribution of multiple scattering can be visually appreciated as intensity tails beyond the coverslip interface. A rightward ROI (blue lines shown in top row of panels a and b) is visualized as individual A-lines as intensity without confocal correction to show where the low-signal coverslip region (shaded in gray) can or cannot be visualized, as well as the magnitude of tail intensity (c). A leftward ROI (green lines shown in top row of panels a and b) is visualized as individual A-lines for the higher $V_f$ cross-sections at the full NA range, where it is possible to see, following confocal correction, a relative lack of NA-based variation in 0.5 µm beads and substantial variation in 3 µm beads (d). These trends were corroborated by the EHF model (e). Scale bars = 1 mm.

Download Full Size | PDF

In Fig. 3(c), depth profiles were plotted for an ROI (blue lines in Fig. 3(a, b), averaged over 10 A-lines, or 40 µm, laterally) where the coverslip border was still appreciable under some conditions, particularly lower NA and lower volume fractions. In Fig. 3(d), a second ROI was examined for higher concentration data where this border was consistently entirely absent, and the intensity decay could be visualized as a function of only the upper droplet layer with minimal to no contribution from the underlying calibration phantom. In this ROI, following confocal function correction, the concentrated 0.5 µm bead data still showed substantially less variation in intensity decay rate with depth (e.g., approximately equivalent $\mu$) over the full range of system NA explored (0.16 - 0.3), whereas the intensity decays for 3 µm bead data exhibited greater variation with NA. This trend was supported by OCT signal simulation using the Extended Huygens Fresnel model [Fig. 3(e)] [19], where results were predicted for each given bead size at respective NA values and focal plane position (here, located 2 mm below the sample surface) to match experimental imaging conditions, as described in Methods.

Our measurements of scattering properties in phantoms supported our hypothesis that singly-scattering samples behave independently of OCT system NA and focal plane position, whereas samples in which multiple scattering contributes more greatly to the overall intensity signal exhibit altered attenuation coefficients and backscattering fractions under these different imaging conditions. These results were also corroborated by the EHF model and direct visualization of the appearance of intensity tails in signal-free areas. Although this behavior is seemingly in contradiction with previous OCT angiography work demonstrating improved suppression of multiple scattering tails below large blood vessels with higher NA systems [30], we hypothesize that these divergent observations may result from different regimes of scattering in each respective environment. Major blood vessels behave as islands of multiple scattering behavior within the inner retina [31], which is otherwise presumed to be a primarily singly scattering tissue. This allows an unperturbed confocal gate to reject multiple scattering to an increasing extent with higher NA. In our highly-scattering phantoms, multiple scattering appears to be the consistently dominant mode of scattering behavior. We hypothesize that such conditions disrupt the confocal gate and compromise its ability to reject multiple scattering, resulting in the decreased $\mu$ and increased $\alpha$ observed in our data with increasing system NA.

3.3 Impact of multiple scattering on degree of polarization measurements

To further validate our findings, we also calculated the degree of polarization (DOP) throughout phantom data (Fig. 4). Previous work has reported that DOP is reduced in the presence of multiple scattering [32,33]. In our SNR-corrected measurements, which account for variable $L_0$ throughout differential NA data, DOP was observed to decrease with increasing particle size and concentration, as expected under these sample conditions which increasingly favor multiple scattering. For phantoms with greater bead diameters and volume fractions, DOP was also progressively reduced with increasing NA. For further validation, DOP was also compared in the calibration ROI of each phantom and found to be consistently high in the 0.35 µm bead data regardless of system NA (see Fig. S3, supplemental document). This DOP-based assessment of multiple scattering impact and its variability under different imaging conditions is an independent and complementary result to findings in intensity-derived scattering properties (Figs. 2, 3).

 

Fig. 4. SNR-corrected DOP was mapped cross-sectionally for 0.5 µm beads (a), 1.0 µm beads (c), 1.5 µm beads (e), and 3.0 µm beads (g) as overlays on signal intensity. Results are summarized in corresponding scatterplots (b, d, f, h) as averaged throughout the representative ROI shown with error bars indicating the 10th and 90th percentile values of the data distribution. Scale bars = 1 mm.

Download Full Size | PDF

3.4 Impact of multiple scattering on optically-derived physical properties

Optical properties accurately measured in singly scattering samples are directly linked to physical sizes and concentrations of the constituent scattering particles in those samples. Since the backscattering fraction is independent of particle density, with known particle and medium refractive indices, it is possible to relate $\alpha$ to the particle diameter, and to further retrieve particle volume fraction based on this $\alpha$-determined diameter and $\mu (z)$. To assess the feasibility of this analysis framework, we used the phantom data summarized in Fig. 2 to calculate particle physical properties based on sample optical properties.

To map $\alpha$ values to particle diameters ($d$), we first used differential scattering cross sections, $\sigma _s(\theta )$, to determine the expected $\alpha$ under our various system NA [Fig. 5(a)] for particles of a continuous range of diameters [Fig. 5(b)] and known refractive indices of polystyrene (particle, 1.57) and water (medium, 1.34) at the OCT system central wavelength. Diameter-$\alpha$ curves were smoothed to generate a mostly-monotonic relationship between $\alpha$ and $d$. Based on these curves, layer-based $\alpha$ values calculated throughout phantom data [Fig. 5(c)] were mapped to particle diameter, $d$ [Fig. 5(d)]. Using $d$ and $\mu (z)$, we also retrieved depth-resolved particle volume fractions ($v_f(z)$) [Fig. 5(c, d)]. Diameter and volume fraction values mapped throughout phantom data were averaged over the superficial ROI shown by the white boxes for summary in scatter plots for calculated particle diameters [Fig. 5(e)] and volume fractions [Fig. 5(f)], with error bars indicating the 10th and the 90th percentiles of parameter values over the analyzed ROI. Analogously to scattering property measurements, $d$ and $V_f$ values for smaller and less concentrated beads samples were generally closer to their physical ground-truth values, whereas larger and more concentrated beads (particularly of 3 µm in diameter) exhibited more divergence from these expected values, and also greater variation with system NA. In particular, the higher $\alpha$ reported for more multiply-scattering samples, particularly at higher NA, resulted in erroneously small diameters, which in turn reduced the overall variation in particle volume fraction. Nonetheless, values retrieved for samples exhibiting even moderate multiple scattering (as assessed by ratiometric $\mu (z)$ methods: see Fig. S5, supplemental document) were often still within reasonable range of their ground truth values, an encouraging result for the application of this type of assessment to OCT data.

 

Fig. 5. Angle-dependent particle scattering cross sections were analyzed at each system NA (angles subtended for each NA are represented by dashed vertical lines in panel a) to report trends in $\alpha$ with particle diameter, $d$, for bead-based phantoms of known refractive index (b). Optical properties [$\alpha (z)$, $\mu (z)$] measured in phantoms (c) were translated into corresponding physical properties (d) as described in Methods. Results are summarized for the $\alpha$-derived particle diameter $d$ in each sample (e) and resulting volume fraction $v_f$ (f) with ground truth values provided as dashed lines. Error bars indicate 10th and 90th percentiles of the data distribution. Scale bars = 1 mm.

Download Full Size | PDF

3.5 Exploiting differential NA imaging to map tissue microstructure

Based on these encouraging results, we applied our analysis framework to tissue data to assess its realistic performance in a biomedical application. A murine brain was imaged freshly post-mortem through a cranial window using the aperture-based OCT system configuration described in Methods for depth-multiplexed, simultaneous differential NA imaging [Fig. 6(a)]. Within the window ROI, a cerebellar tumor had been implanted and subsequent histology confirmed its location post-imaging [Fig. 6(b)]. The cancer cell line used was specifically generated to recapitulate typical nuclear-cytoplasmic ratios (NCR) observed in human glioblastoma [34]. Cell nuclei were segmented throughout the histologically processed brain slice to generate a map of relative densities of cell nuclei, which were determined to be much higher in the tumor ROI [Fig. 6(c)]. More quantitative metrics characterizing nuclei diameter and NCR are provided in Fig. S6 of the supplemental document. Intensity data [Fig. 6(d)] were analyzed as described in Methods, and attenuation coefficient data were extracted [Fig. 6(e-h)] following confocal correction. In general, $\mu (z)$ values were found to be similar for each system NA in healthy cerebellar ROIs but some variation was visible in areas of suspected tumor burden (blue dashed ROI), and also near blood vessels (red dashed arrows). In en face intensity data, the tumor was visible as a region with relatively higher intensity, and a major blood vessel was visible as a region with relatively lower intensity [Fig. 6(i)]. In the corresponding $\mu (z)$ map [Fig. 6(j)], both the tumor and blood vessel had elevated attenuation coefficient values, but in the backscattering fraction map, higher $\alpha$ in the tumor area contrasted with lower $\alpha$ in the vasculature [Fig. 6(k)].

 

Fig. 6. Brain imaging in a mouse was performed through a cranial window (a) to access a cerebellar tumor (b), which was subsequently processed histologically and analyzed in terms of density of detected cell nuclei areas (c). From the aperture-based differential NA imaging volume (d), $\mu (z)$ was analyzed independently for higher and lower NA signal copies (e-h) revealing the tumor-bearing region (blue dashed ROI) and blood vessels (red dashed arrows) which were also visible in the en face intensity map (i) at a location 100 µm below the tissue surface. Corresponding en face $\mu (z)$ (j) and $\alpha$ (k) maps, calculated from the high NA signal channel, were used to calculate the multiple scattering index (l), $V_f$ (m), and $d$ (n). Scale bars=1 mm.

Download Full Size | PDF

In assessing relative differences in $\mu (z)$ captured at higher and lower NA as a function of multiple scattering contribution, taking a direct ratio of $\frac {\mu _high(z)}{\mu _low(z)}$ as determined in phantom data produced significant artifacts in the noisier tissue data. As a more robust metric, to highlight areas of potential multiple scattering as features with a high overall scattering coefficient and also a strong difference in $\mu (z)$ for higher ($\mu _{hi}(z)$) and lower ($\mu _{lo}(z)$) system NA, we defined a Multiple Scattering Index, $\mu _{hi}(z)|\mu _{hi}(z) - \mu _{lo}|$, which was found to be elevated in both the tumor and blood vessel [Fig. 6(l)]. In phantoms, this metric also produced trends that recapitulated the purely ratiometric differential NA behavior (see Fig. S5, supplemental document). In tissue, this finding was consistent with both the dense particle makeup of the tumor and the known multiple scattering behavior of red blood cells as imaged with OCT [31]. Finally, $\alpha$ and $\mu (z)$ values were used to determine the volume fractions [Fig. 6(m)] and diameters [Fig. 6(n)] of scattering particles in tissue based on the high-NA signal copies, using an estimated particle refractive index of 1.38.

These results demonstrate the ability of OCT-based scattering properties to retrieve physical properties of particles in tissue with relative accuracy even in the presence of multiple scattering, including the decoupling of concentration and diameter for diversely scattering particles. However, based on the observed impact of multiple scattering on optically-derived physical properties in phantoms of known composition (Fig. 5), $d$ and $v_f$ values determined in areas of high multiple scattering index are expected to diverge from ground truth values. These discrepancies may be appreciated in both the brain tumor and vasculature. The blood vessel should report a scattering particle concentration comparable to typical murine hematocrit levels, which can vary between 0.4-0.6 [35], but optically-determined $v_f$ were much lower, and comparable to that of surrounding cerebellar tissue. The optically-determined $d$ within the tumor region were similarly lower than those of the surrounding tissue area. However, more intensive assessment of tumor histopathology reported nuclei in this area with diameters comparable to and larger than those of the healthy cerebellum (see Fig. S6, supplemental document). Further, while the $v_f$ determined in the tumor region was substantially higher than that of the surrounding tissue, it was still more than an order of magnitude below corresponding NCR values calculated from the histology. These findings indicate that while our proposed analysis may be effective in decoupling scattering behavior driven by particle sizes and concentrations, particularly in areas determined to be significantly impacted by multiple scattering, the optically-determined $d$ and $v_f$ in these areas should be interpreted as only as relative values with a large margin of error.

For our assessment of $d$ and $v_f$ throughout tissue data, we used the higher NA signal copy due to the finer image quality provided by this higher resolution data. However, we recognize that based on our evaluation of the impact of multiple scattering as captured by different system NA on calculations of $d$ and $v_f$ throughout phantom data, it could also be advantageous to use the lower NA signal copy, based on the potentially lesser contribution of multiple scattering in this channel. Although this presents a possible trade-off between image quality and measurement accuracy, we envision that the information from the two channels could be further combined such that in tissue areas with a high multiple scattering index, the NA channel with lower contribution of multiple scattering (i.e., higher $\mu (z)$ and lower $\alpha$) would be used for $d$ and $v_f$ calculation, whereas in areas with low multiple scattering index, the high NA channel would be used consistently. In the future, we hope to determine a threshold of multiple scattering index based on evaluation of further tissue samples at which to apply this multi-NA post-processing strategy, which further motivates our differential NA imaging platform for more informative tissue microstructure characterization.

4. Discussion

Attenuation coefficients and other optical properties are attractive for characterizing biological tissues due to the quantitative information that they provide beyond the OCT intensity signal alone. However, several challenges complicate the usage of optical properties to accurately and intuitively provide sub-resolution microstructural information profiling tissues of interest. Previous work has shown that it is possible to measure optical properties in homogeneous phantoms of known composition to a high degree of accuracy compared to ground truth values [13]. Further methods have been developed to extend such analysis to more complicated (e.g., layered) samples under more flexible imaging conditions [6,11,14]. In this work, we leveraged these and additional methods to measure attenuation coefficients and backscattering fractions in layered phantoms under a variety of imaging conditions and assessed their accuracy as a function of system geometry and expected multiple scattering contribution. We have further demonstrated the extension of this approach to a tissue sample containing healthy and pathological regions.Specifically, we chose to use our methods to characterize a murine brain tumor, due to the compelling rationale of identifying structural elements with the potential to discriminate cancerous from healthy tissue (e.g., elevated nuclear-cytoplasmic ratio). However, we note that these methods are not limited to brain tumor, but are generally applicable to any cancerous or dysplastic tissues, as well as a number of previously identified applications, such as characterization of atherosclerotic plaques [7], burn wounds [36], and optic nerve pathologies [9].

Based on the extraction of both layer-resolved $\alpha$ and a depth-resolved $\mu (z)$, we have proposed methods to retrieve scattering particle diameter and volume fraction directly from these optical properties. This is a useful technique for two reasons. Firstly, compared to the usage of attenuation coefficient to profile tissue microstructure alone, our approach effectively decouples the independent contributions of both particle diameter and concentration to $\mu (z)$. Secondly, outputting particle physical properties rather than optical properties yields more fundamentally intuitive sample information with greater interpretability within existing clinical frameworks (e.g., percent of enlarged nuclei above a certain size threshold as a hallmark of dysplasia or neoplasia).

However, we note a number of limitations motivating future work to continue improving our approach. The initial step in our physical particle property processing pipeline relies on accurate $\alpha$ measurement to relate particle diameter, and this information is in turn used in calculating particle volume fraction. Unlike $\mu (z)$, $\alpha (z)$ requires concrete knowledge of the OCT system incident light intensity, $L_0$ (though as we have previously noted, even without a precise $L_0$, relative $\alpha$ measurements can still yield valuable sample contrast [14]). Less work has been done previously in fine-tuning methods for accurate $\alpha$ calculation compared to those for $\mu (z)$ calculation. Various assumptions that we used, particularly in smoothing $\alpha$–$d$ curves (justified by angular distribution within measurements, manufacturer-reported variability in particle size, and the wavelength range of the OCT laser source) merit further investigation. Additionally, our automated methods for layer-based $\mu (z)$ and $\alpha (z)$ processing [14] can be subject to error in the case of imaging artifacts such as strong specular reflections and are generally more difficult to implement in thinly layered samples. This posed a challenge in some of our measurements of 1.5 µm bead phantoms, which suffered from greater lateral droplet spread and thus a thinner sample for analysis, and contributed to the greater uncertainty observed for these measurements [Fig. 2(e, f), Fig. 5(e, f)]. The ability of strong multiple scattering to distort beam focus and create artifacts [see surface aberrations in Fig. 3(a,b)] could also impact our optical property extraction methods. Efforts are ongoing to continue to improve these image processing routines and increase their robustness in samples with challenging geometries.

In tissue, cell nuclei can exceed 10 µm in diameter, placing them in the upper range of sizes analyzable by Mie theory and OCT. However, their reported refractive indices are much lower than that of polystyrene, decreasing their relative scattering cross section values and placing them more comfortably in a regime where our analysis is applicable. In phantom measurements, particularly due to the high refractive index of polystyrene, we were limited to exploring particle sizes and densities up to $d = 3$ µm and $v_f = 0.024$, respectively, before observing significant degradation in optical and physical property quantification accuracy. The expansion of methods for phantom fabrication to better mimic realistic in vivo environments would be useful in further experiments comparing the scattering properties measured in a sample to known ground truth values. We also relied on the known refractive index of polystyrene to inform inverse calculations in phantoms (e.g., optically-deriving $d$ and $v_f$). Reported refractive indices of nuclei and nucleoli vary in literature, with typical estimates falling between 1.35 and 1.45 for near-infrared illumination wavelengths near 1.3 µm [26–28]. We demonstrate in Figure S7 of the supplemental document that varying the estimated refractive index of cell nuclei over this range impacts resulting $d$ and $v_f(z)$ measurements using a given refractive index as an input value, but maintains overall trends in particle size and concentration. Additionally, previous work demonstrating that tissue can be considered to have a fractal distributions of scattering particles with refractive indices varying over a given range further justifies the smoothing of $\alpha$–$d$ curves, [37,38] and is another aspect that we hope to explore further in the future. Without precise knowledge of refractive index or tissue fractal dimension, $d$ and $v_f(z)$ measured in tissue could be considered as relative particle size and density indices rather than their direct physical values, acknowledging the limitations of these assumptions yet still has the capacity to provide valuable insight into tissue behavior.

The methods presented here for translating optical properties into physical properties are applicable only in the case of elliptically-shaped scattering particles which meet the criteria for the application of Mie theory [4]. In tissue, these criteria apply to many cellular organelles such as nuclei, but not to fibrous and extracellular matrix components such as collagen and elastin. Using polarization-sensitive (PS-) OCT, [39,40] it is possible to identify such fibrous areas based on their birefringence, which could enable the exclusion of such fibrillar structures with birefringence above a given threshold from Mie-based processing. The anisotropic scattering behavior of fibrillar structures in tissue has previously been assessed using cylindrical scattering theory [41], which is more complicated than the isotropic scattering of spheres due in part to necessary consideration of angular orientation of individual scatterers with respect to the imaging system [42]. Nonetheless, similarly quantitative assessment of $\mu$ and $\alpha$ throughout fibrillar tissue regions based on cylindrical scattering theory remains an interesting future direction.

Another key element of this work was in exploring the ability of differential NA and focal plane measurements to capture differing levels of multiple scattering in imaged samples. In phantoms, knowledge of theoretical Mie-based $\mu$ and $\alpha$ helped to determine the presence of multiple scattering through its reduction of measured optical properties compared to these ground truth values. These predictions helped to confirm our hypothesis that in predominantly singly scattering samples, $\mu$ and $\alpha$ measurements would be unaltered by NA or focal plane modulation, whereas in multiply scattering samples, changing these parameters also altered measured $\mu$ and $\alpha$. Although this result is in line with intensity behavior predicted by the EHF model [19], our finding that higher NA imaging conditions actually captured more multiply scattered photons was contradictory to prior work in OCT angiography [30] and general expectations of confocal gating. This finding is potentially explained by the underlying theory of the EHF model, which determines the contribution of multiple scattering to the intensity signal based on distortion of the illumination beam with depth in a sample with strong forward scattering [19]. Additionally, previous work in dynamic light scattering has reported insensitivity to the NA of the transition to a multiple-scattering regime in media with large scattering anisotropy [29]. However, more experimental work or advanced theoretical modeling [43] corroborating this finding in phantom and tissue samples is necessary to determine whether the mechanism at play is indeed the disruption of the confocal gate to the extent that it fails to reject multiply scattered photons.

We also note that all of our imaging conditions featured a focal plane that was located within samples at sub-surface locations ranging from 1 mm to 2 mm. We chose these locations based on the general advantage of better image contrast that accompanies focusing within a sample, thus evaluating our methods within a realistic framework applicable to many imaging applications. However, it is appreciable in the data shown that NA and focal plane location each separately impact multiple scattering behavior (Fig. 2, Fig. 5), with deeper focal planes generally resulting in greater multiple scattering contribution. It is therefore possible that multiple scattering behavior with NA could differ particularly if focal planes were located at more superficial locations within the sample with respect to the system $Z_R$ or even above sample surfaces, offering another topic for further investigation.

It is useful to know where multiple scattering significantly impacts collected OCT intensity within an imaged volume [18], both to appreciate where more quantitative analysis may be rendered inaccurate and as its own potential form of contrast. Our differential NA imaging platform offers this capability, as well as the potential to explore the relative strength of multiple scattering through our Multiple Scattering Index. Each of the methods we have proposed to accomplish differential NA imaging can be achieved through relatively simple alterations to standard OCT imaging hardware [Fig. 1(a)]. The aperture-based approach is advantageous for tissue imaging, especially in vivo, due to simultaneous signal acquisition for precise image copy registration. One disadvantage of this approach is the reduction in imaging range available for each individual signal copy, which cannot exceed 1/3 of the full imaging range. PME thickness should also be optimized to suit this criterion, and separate each signal copy by this axial extent. Since most tissues and biological samples of interest fully attenuate incident light over a depth of less than 3 mm, the baseband imaging range of our OCT system (11 mm) is sufficient to effectively multiplex this information without overlap. The thickness of our PME was chosen to enable a signal copy separation of 1.6 mm, which eased requirements on element construction and alignment, and was sufficient for the experiments presented here. However, it is possible that this extent of separation could be insufficient in some biomedical applications analyzing samples with lower scattering. In some areas of our imaged brain tissue, there appeared to be a small amount of overlap between the middle signal copy and the bottom signal copy (high NA). This effect was more pronounced in areas of lower attenuation coefficient, corresponding to the healthy brain regions, whereas the highly-scattering tumor-burdened region fully attenuated the intensity signal prior to the beginning of the following signal copy. In cases of overlap, the tail of the intermediate NA signal was found to be approximately 1–2 orders of magnitude lower than the intensity of the high NA signal in the overlap region. Considering that the impact of the overlapping signals on $\mu$ is only measurable in terms of a change in slope, the quantitative impact of this minor overlap is negligible despite its pronounced appearance in logarithmically-scaled B-scans [(e.g., Fig. 6(d)]. Further, given that the measured scattering properties in the healthy brain regions appeared to be well-behaved (i.e., low multiple scattering index) compared to the tumor-burdened region, we do not believe that any cross-term and high NA signal overlap in our system design was sufficient enough to impact the scattering property results. Alternatively, focal plane modulation could also be deployed to obtain multiple scattering contrast, but is less easily achieved through simultaneous acquisition. Although either NA-based or $z_f$-based differential imaging requires further hardware or data collection beyond standard OCT imaging conditions, the additional information gained stands to greatly enhance the interpretation of optical property measurements and derived particle physical properties. Throughout both phantom and tissue data, the impact of multiple scattering can be observed on both measured optical properties and corresponding optically-derived physical properties; these measurements in regions with strong multiple scattering require less quantitative and more nuanced interpretation. In future work, we hope to demonstrate the feasibility of differential NA imaging not only in benchtop OCT systems, but also in catheter-based devices to potentiate the extension of our work to in vivo settings.

5. Conclusion

In this work, we have aimed to characterize the accuracy of optical scattering properties calculable from OCT intensity data–namely, attenuation coefficients and backscattering fractions–and improve the interpretability of the information that they convey. To do so, we have measured well-characterized tissue-mimicking phantoms under varying OCT system conditions, which has given us greater insight into how modulating system parameters, such as NA, can impact the contribution of multiple scattering to the intensity signal and also the reliability of intensity-derived optical properties. Further, we have formulated methods to translate optical properties into physical properties based on relationships defined by Mie theory. Leveraging our differential NA imaging methods, we have also characterized the impact of multiple scattering on the determined sizes and concentrations of scattering particles in both phantoms and in tissue ex vivo. Although our methods for particle size and concentration determination show promise in enhancing the information implicit in optical properties, the accuracy of this information remains contextualized within the degree of multiple scattering impact as measurable by our differential NA imaging methods. Many interesting future directions exist for this work, including further validation of observations, more quantitative interpretation of scattering information within multiply scattering sample regions, and extension of our approach to additional areas of clinical imaging interest. Overall, we hope that our methods provide a useful approach for maximizing the utility of OCT-based optical property measurements, particularly for microstructural tissue characterization in biomedical settings.