2.1.

Generic Model of the Optical Coherence Tomography Signal

The most common approach to date to quantify optical properties from OCT images is by the means of nonlinear least squares (NLLS) curve-fitting of an OCT signal model to the measured data. If the signal model is parameterized properly, such that no mutual dependence exists between the fit parameters, the covariance matrix obtained from the NLLS fit procedure yields uncertainty estimates on the fitted parameters such as standard deviation of multiple fits for which the region of interest (ROI) is changed by 10%. We adopt the model by Xu et al.,¹³ which assumes only single backscattered light contributes to the OCT signal. Light to and from the backscatter location at depth $z$ is attenuated due to scattering and absorption. This is accounted for by the attenuation coefficient ${\mu }_t$, such that ${\mu }_t\mathrm{d}z$ is the probability of interaction in path length $\mathrm{d}z$ and $1-{\mu }_t\mathrm{d}z$ is the transmission probability. The intensity incident at the backscattering site is thus reduced by a factor $(1-{\mu }_t\mathrm{d}{z}_1)\times (1-{\mu }_t\mathrm{d}{z}_2)\times \dots \times (1-{\mu }_t\mathrm{d}{z}_N)$ compared with the intensity impinged on the tissue. If the sample is homogeneous and the product ${\mu }_t\mathrm{d}z$ is small, we write $\underset{n\to \infty }{\lim }{[1-{\mu }_t(z/n)]}^n={\mathrm{e}}^{-{\mu }_{t}z}$. Several preclinical and phantom studies demonstrated that, for most tissues, the OCT signal amplitude versus depth $z$ is best described by this single exponential decay curve, combined with descriptions for the confocal PSF^14,15 and sensitivity roll-off in depth for frequency domain OCT systems. The resulting expression for the mean-squared OCT depth signal after noise subtraction for $z\ge z_0$ is^16,17

The expression for the sensitivity roll-off in depth $H(z)$ is given in terms of sampling and optical resolution (in wavenumbers) in the following equation:²¹

The last term of Eq. (1), ${\mu }_{B,\mathrm{NA}}\exp [-2{\mu }_{\mathrm{OCT}}(z-z_0)]$, contains the optical properties of the sample, but is influenced by system properties as well. The backscattering coefficient within the NA (${\mu }_{B,\mathrm{NA}}$) clearly depends on the confocal properties of the OCT system. Throughout this paper, we assume ${\mu }_{B,\mathrm{NA}}$ only differs per sample, but is constant in depth for homogeneous phantoms for each individual sample. The decay constant ${\mu }_{\mathrm{OCT}}$ depends on the scattering coefficient ${\mu }_S$ and absorption coefficient ${\mu }_A$. For weakly scattering samples, with negligible contributions from multiple scattered light, ${\mu }_{\mathrm{OCT}}={\mu }_S+{\mu }_A$. The experimentally obtained value may be influenced by the contribution of multiple scattering to the OCT signal, determined by the confocal properties of the system and on the angular scattering characteristics of the sample. Both a lower NA and samples with highly forward directed scattering combined with a high scattering coefficient, corresponding to values of the scattering anisotropy $g$ close to unity, lead to larger contributions of multiple scattered light in the detected signals.¹² Thus, for these situations, detection of multiple scattered light leads to a decreased decay of the OCT signal in depth, i.e., ${\mu }_{\mathrm{OCT}}<{\mu }_S+{\mu }_A$. In that case, the following equation is applicable:

2.2.

Optical Properties of the Phantom

Consider an isotropic DRM consisting of identical hard spherical particles with known number density $\rho$ [# $\text{particles}/{\mathrm{m}}^3$]. The particles cannot overlap in volume. For the individual particles, Mie theory²³ yields the differential scattering cross section ${\sigma }_S(\theta )$ (${\mathrm{m}}^{-2}{\mathrm{sr}}^{-1}$), where $\theta$ is the scattering angle. The total scattering cross section is found by integrating over the angular coordinates while taking into account the phase differences in the scattered fields from the different particles (which are determined by their distinct positions $r$)

Fig. 1

A graph of the structure factor at volume fractions $f_v=0.001$, 0.2, and 0.5 as a function of $Dq=(kD)*2\sin (\theta /2)$. The curve was calculated using the algorithm described by Tsang et al.,²⁴ as the numerical Fourier transform of the Percus–Yevick pair correlation function (PCF).

The scattering coefficient (the product of particle density and scattering cross section, ${\mu }_S=\rho {\sigma }_S$) is then calculated as

2.3.

Multiple Scattering Model

For a high value of the sample’s scattering coefficient (${\mu }_S$) and a scattering anisotropy ($g$) close to 1, an increased contribution of multiple scattered light to the OCT signal will lead to a reduced measured ${\mu }_{\mathrm{OCT}}$. The most comprehensive model to date describing this effect is based on the EHF principle introduced by Schmitt and Knüttel¹¹ and Thrane et al.¹² Here, in the paraxial approximation, the mean squared OCT signal—excluding the effects on attenuation of the confocal PSF and sensitivity roll-off—is given by a contribution of three terms: the single-backscattered field [as in Eq. (1)], the multiple (forward) scattered field, and a coherent cross term between these two fields as

Here, $w_H(\mathrm{z})$ is the local beam waist in the absence of forward scattering (e.g., of the single-backscattered beam):

3.1.

System Parameter Calibration

All recorded data were corrected for system-induced attenuation, i.e., caused by the confocal PSF and sensitivity roll-off as described in Sec. 2. Due to the low NA and high spectral resolution of our system, both are slowly decaying functions of depth and can in practice be approximated by a single exponential decay with attenuation constant ${\mu }_{\mathrm{CAL}}$. Thus, Eq. (1) can be modified to

3.2.

Sample Preparation

Monodisperse silica beads (Kisker Biotech, Steinfurt, Germany) with measured mean diameters of 0.47, 0.70, 0.91, and $1.60\mu \mathrm{m}$ were suspended in distilled water for a monodisperse concentration series. For an accurate measurement of the size distribution, the silica beads were imaged using transmission electron microscopy (TEM) (Philips CM-10) (Fig. 2). The TEM images were acquired and analyzed following the protocol described in previous work.¹⁰ The obtained values of the mean diameter and standard deviation of the bead sizes are shown in Table 1. To prevent aggregation of the beads in suspension, 0.03 mM of sodium dodecyl sulfate was added to all samples.²⁶ The 0.47- and $0.70\text{-}\mu \mathrm{m}$ sized beads were diluted to 10 samples with volume fractions ranging from 0.02 to 0.2. The 0.91- and $1.60\text{-}\mu \mathrm{m}$ sized beads were diluted to 5 samples with volume fractions from 0.02 to 0.1. Suspensions with a higher volume fraction of silica were not studied due to aggregation of the beads. All samples were 4 to 6 times alternatingly vortexed and sonicated for 15 min.

Table 1

Mean and standard deviation of the size distribution of monodisperse silica beads measured with transmission electron microscopy.

Fig. 2

Histogram of the size distribution of monodisperse silica beads (PSI 1.0) measured with transmission electron microscopy.

3.3.

Optical Coherence Tomography Data Analysis

We used a custom-written code (LABVIEW 2013, National Instruments) for the analysis of the acquired OCT data. The attenuation (${\mu }_{\mathrm{OCT}}$) was obtained from the data by NLLS fit of a single exponential decay model $y=\mathrm{AMP}*\exp (-\mu z)+y_0$ to the OCT signal in depth. The data were first averaged over at least 1000 A-lines. The selection of the fit was done manually where the start region was chosen to approximately be $50\mu \mathrm{m}$ behind the cuvette boundary. The end point of the fit was chosen such that within the ROI, upon visual inspection, the signal appeared to follow a single exponential decay curve. An offset added to the fit model was fixed at the average noise level directly at the backside of the cuvette. The amplitude and decay constant $\mu$ were independently varied parameters. The fit was repeated 100 times, during which the boundaries of the ROI were varied within 10% of their operator-selected values. The mean value of the 100 measurements was taken as the final $\mu$ value, with uncertainty estimated as the standard deviation of the fits. The obtained attenuation $\mu$ was corrected for system parameters, as described in Sec. 3.1, yielding ${\mu }_{\mathrm{OCT}}=\mu -{\mu }_{\mathrm{CAL}}$.

The backscattering coefficient, ${\mu }_{B,\mathrm{NA}}$, was determined from the amplitude of the backscattered OCT signal^18,19 in a superficial ROI starting at $z_{\mathrm{ROI}}=0.03\mathrm{mm}$ below the glass–sample interface. The ROI is chosen as such to exclude contribution of multiple scattered light. The region extends only 0.026 mm (4 pixels) in depth to minimize the effects of attenuation, which was verified by ensuring that the OCT amplitude speckle contrast in the ROI was 0.52,²⁷ corresponding to fully developed speckle, which can only be achieved if no appreciable attenuation in depth takes place in the ROI.²⁸ By Eq. (11), in this ROI, we have $⟨A^2(z_{\mathrm{ROI}})⟩=\alpha \exp (-{\mu }_{\mathrm{CAL}}{z}_{\mathrm{ROI}}){\mu }_{B,\mathrm{NA}}$ which features only system-dependent parameters. This allows us to define a scaling factor between measured amplitude and backscatter coefficient ${\mu }_{B,\mathrm{NA}}=s{f}_{Ø}⟨A^2(z_0)⟩$, where the scaling factor is determined separately for each particle diameter using the highest volume fraction $f_{v\mathrm{MAX}}$: $s{f}_{Ø}={[{\mu }_{B,\mathrm{NA}(\mathrm{MIE}-\mathrm{PY})}/⟨A^2(z_0)⟩]}_{Ø,f_{v\mathrm{MAX}}}$.

3.4.

Multiple Scattering

A practical approach for assessing the influence of multiple scattering for confocal reflectance microscopy and OCT was proposed by Jacques et al.²² Based on Monte Carlo simulation, a reference grid is constructed, allowing for the mapping of scattering the coefficient and scattering anisotropy to the amplitude and decay rate and vice versa.^22,29 We adopt this approach, but rather than relying on Monte Carlo simulations, we use the analytical EHF model proposed in Sec. 2.3 to simulate the OCT signal at a range of values of the scattering coefficient and anisotropy and, subsequently, fit the single exponential decay to the simulated signal to retrieve the corresponding predicted value of ${\mu }_{\mathrm{OCT}}$. This procedure is illustrated in Fig. 3.

Fig. 3

Simulated extended Huygens–Fresnel (EHF) model for optical coherence tomography (OCT) (green dashed line), experimental OCT curve (black solid line) fitted with a single exponential decay within the ROI for samples of 0.1 volume fraction: (a) $0.47\text{-}\mu \mathrm{m}$ beads and (b) $1.60\text{-}\mu \mathrm{m}$ beads. The simulated single scattering ($\mathrm{Mie}+\mathrm{PY}$) curve (blue dotted line) is plotted for reference. The exponential decay of the fits on the experimental curves is (a) ${\mu }_{\exp }=5.0{\mathrm{mm}}^{-1}$ and (b) ${\mu }_{\exp }=10.8{\mathrm{mm}}^{-1}$. After correction of system-induced attenuation (${\mu }_{\mathrm{CAL}}=1.1{\mathrm{mm}}^{-1}$), the decay rates match closely with the decay of the simulated curves: (a) ${\mu }_{\mathrm{sim}}=3.6{\mathrm{mm}}^{-1}$ and (b) ${\mu }_{\mathrm{sim}}=9.7{\mathrm{mm}}^{-1}$.

3.5.

Calculations

To calculate predictions for ${\mu }_{\mathrm{OCT}}$ and ${\mu }_{B,\mathrm{NA}}$, we have written a Labview code (Labview 2013, National Instruments). We implemented Mie calculations from Ref. 30 and calculated the differential scattering cross section for each particle by integrating over the size distribution measured with TEM (Table 1 and Fig. 2). The refractive index of the silica beads was fixed at 1.425 (within the range provided by the manufacturer) and converted to 1309 nm using the dispersion curve from Malitson³¹ ($n=1.419$ to 1.449). The refractive index of the medium was assumed to be equal to that of water (1.324).³² The Percus–Yevick structure factor $S(\theta )$ was calculated in Ref. 24, using volume fraction and mean particle diameter (Table 1) as input parameters. By combining the Mie-calculated ${\sigma }_S(\theta )$ with this structure factor, we calculate ${\mu }_S$ [Eq. (6)], ${\mu }_{B,\mathrm{NA}}$ [Eq. (7)], and rms scattering angle ${\theta }_{\mathrm{rms}}$ for each sample ($\mathrm{Mie}+\mathrm{PY}$). The EHF calculations were based from Ref. 12 and took the system properties center wavelength, focal length, and Rayleigh length as input parameters, as well as the medium refractive index and the ${\mu }_S$ and ${\theta }_{\mathrm{rms}}$ ($\mathrm{MIE}+\mathrm{PY}$) for each sample. To quantify the accuracy of the prediction, for each model, the standard error of the estimate ($s$) was calculated as $\surd \sum {({\mu }_{\mathrm{OCT},\text{EXPERIMENTAL}}-{\mu }_{\mathrm{OCT},\text{PREDICTION}})}^2/N$, where $N$ is the number of data points (volume fractions). The lower $s$ corresponds to a better predictive value. An uncertainty estimate on $s$ was calculated using a bootstrapping method,³³ we calculate $N$ values of $s$ from an ($N-1$) sized dataset by leaving out each value of the original dataset once. The standard deviation of the resulting sequence of $s$-values is taken as an uncertainty estimate (Table 2). A similar analysis was performed for the backscattering coefficient (Table 3).

Table 2

The standard error of estimates is calculated to evaluate the match between the experimental data and the theoretical predictions for μOCT by the three models (Mie, Mie–PY, and Mie–PY–EHF). Uncertainties are estimated using a bootstrapping method.

Table 3

The standard error of estimates is calculated to evaluate the match between the experimental data and the theoretical predictions for μB,NA by the two models (Mie and Mie–PY). Uncertainties are estimated using a bootstrapping method.

5.1.

Concentration-Dependent Scattering Versus Multiple Scattering

For higher scattering coefficients and low values of ${\theta }_{\mathrm{rms}}$ ($\sim g$ close to 1), multiple forward scattered light has a large contribution to the OCT signal. Corresponding DRM for calibration purposes requires high volume fraction of the scatterers. However, for high volume fractions, the bulk optical properties of the DRM cannot be extrapolated from the optical properties of a single scatterer simply through scaling by volume fraction. Rather than adding scattered intensities, scattered fields (amplitude and phase) from the individual particles in the DRM should be added (taking into account their interference) to calculate the optical properties. This situation is commonly referred to as “concentration-dependent scattering” as the bulk optical properties now depend nonlinearly on volume fraction. The effect of concentration-dependent scattering is that the scattering coefficient ${\mu }_S$ scales as a lower-than-linear rate with volume fraction [${\mu }_S<(f_V/V_P){\sigma }_{\mathrm{SCAT}}$]; Sec. 2.2), leading to a lower ${\mu }_{\mathrm{OCT}}$ compared with the “concentration-independent scattering” case. Increased detection of multiple-forward scattering leads to further reduction of ${\mu }_{\mathrm{OCT}}$. Therefore, care must be taken to separate both effects in data analysis. Here, we predict ${\mu }_{\mathrm{OCT}}$ by first calculating ${\mu }_S$ and ${\theta }_{\mathrm{rms}}$ for our DRM samples by using Mie theory (scattering properties of individual spheres) and the Percus–Yevick structure factor (concentration-dependent scattering). These values are then input to the EHF equation that describes the OCT signal in the multiple scattering regime. ${\mu }_{\mathrm{OCT}}$ is derived by fitting a single exponent to the simulated data. The measured values for ${\mu }_{\mathrm{OCT}}$ correspond closely to these predictions (Fig. 4). For samples with a low amount of forward scattering (low $g$-value, Fig. 5), no contribution of multiple scattered light is observed. Although the EHF model is limited by the paraxial approximation, its prediction does not deviate from the concentration-dependent scattering curve and is in good agreement with the experimental data, which suggests a broad applicability of the EHF model. For the case of weakly scattering media with moderate scattering anisotropies, the simpler Mie–Percus–Yevick model may be used for accurate prediction of ${\mu }_{\mathrm{OCT}}$. The backscattering coefficient, ${\mu }_{B,\mathrm{NA}}$, is best predicted by Mie–Percus–Yevick theory for all particle sizes. Even though no statistical difference is found between the prediction based on Mie theory or Mie–Percus–Yevick theory, only the latter model yields accurate predictions for both parameters. The physical explanation is that, as the volume fraction increases, the increment in scattering coefficient decreases [e.g., the integral part of Eq. (6) is a decreasing function of volume fraction] while at the same time increasingly more light is scattered to the backward direction [e.g., the integral part of Eq. (7) remains almost constant]. Effectively, the backscattering coefficient scales linearly with volume fraction in both models.

We conclude that our approach correctly separates the contributions of concentration dependent and multiple scattering, substantiated by the correct prediction of the backscattering coefficient ${\mu }_{B,\mathrm{NA}}$ from the same calculations (Fig. 6). We measure ${\mu }_{B,\mathrm{NA}}$ at depths in the order of 1 mean free path to minimize attenuation of the signal within the ROI (minimum mean free path is ${\mu }_S^{-1}=30\mu \mathrm{m}$ for $Ø=1.60\text{-}\mu \mathrm{m}$ particles, $\text{volume fraction}=0.1$). Moreover, we have recently demonstrated correct prediction of scattering coefficients ${\mu }_S$ by Mie–Percus–Yevick theory, using transmission-based low-coherence interferometry.¹⁰ In that experiment, the contribution of small-angle multiple forward scattering is greatly reduced by probing only the ballistic light transmitted through the sample.

The weighting of the differential scattering cross section with the structure factor (concentration-dependent scattering) also leads to more light being scattered into the backward direction. Ultimately, these two contrary effects result in an almost linear increase of the backscattering coefficient ${\mu }_{B,\mathrm{NA}}$ with volume fraction. Clearly, this result is not general but depends on the exact scattering phase function of the used sample/tissue and the confocal configuration of the OCT system as predicted by the presented theory.

To eliminate contributions of multiple scattering, Xu et al.¹³ propose to confine the ROI of the fit between $I_{\max }$ and $I_{\max }/\mathrm{e}$. As can be seen from the EHF equations [Eqs. (8)–(10)], the multiple scattering starts to contribute directly from $z=z_0$, or $\mathrm{\Delta }z=(z-z_0)=0$. The nonscattered beam attenuates with $\exp (-{\mu }_s\mathrm{\Delta }z)$ and the (small angle, forward) scattered beam grows in intensity with [$1-\exp (-{\mu }_s\mathrm{\Delta }z)$]. This latter beam, however, broadens upon propagation [Eq. (10)], so that coupling back into the single-mode fiber of the OCT system will become more difficult with increasing $\mathrm{\Delta }z$. The contribution of multiple scattering, therefore, is a balance between these two factors. At small values of $\mathrm{\Delta }z$, the power in the scattered beam is small but, since the spread (increased divergence of the scattered beam) is small, it is more effectively coupled back. In contrast, for larger values of $\mathrm{\Delta }z$, the power in the scattered beam is increased, but due to the increased spread, the coupling efficiency is worse. At very large depths, the highly increased magnitude of the scattered beam starts to dominate the single-backscattered contribution.

Two solutions proposed in literature are commonly adopted to take into account the multiple scattered light contribution to the OCT signal: the EHF model described here and a Monte Carlo derived model introduced by Jacques et al.²² We have adopted the EHF model in this work, but adjusted for the confocal PSF for a diffuse reflection, where we assume that the lateral phase of the probing beam is lost upon reflection. In practice, this is accounted for by setting the Rayleigh length in the medium to 2 times the value measured for specular reflection (Sec. 2.3).^15,20 In the other approach, an empirical equation is derived by Monte Carlo simulations to approximate the contribution of multiple scattering.^22,29 Levitz et al.²⁹ use this approach and translate ${\mu }_{\mathrm{OCT}}$ and amplitude to ${\mu }_S$ and $g$ based on calibration of the peak signal of a known interface (scattering phantom with known optical properties). The strength of our approach is that the variation of any parameter in the model may be changed easily to study its influence on the retrieved optical properties (of course, within the range of validity of the model itself, e.g., dominant forward scattering). The empirical model by Jacques et al.²² would require a new set of simulations if, for example, the detection NA changes. In concept, however, both methods are equivalent in that the mapping function $f({\mu }_S)$ in Eq. (1) is studied either using an analytical model or Monte Carlo simulations.

5.2.

Limitations and Clinical Implications

Our model calculations assume applicability of Mie–Percus–Yevick calculations to account for concentration-dependent scattering and the EHF model to account for multiple small-angle forward scattering. Input for the Mie–Percus–Yevick calculations is size and refractive index of the silica spheres making up the DRM. In our present calculations, we use the mean diameter of the size distribution as measured by TEM. Alternatively, an effective “scattering diameter” can be computed by weighting the size distribution with the corresponding scattering cross sections which yield only small deviations from our present calculations ($<1\%$ in ${\mu }_{\mathrm{OCT}}$). By averaging over all particle sizes, we can draft a rule of thumb expression for the concentration-dependent scattering coefficient ${\mu }_S$ as a function of volume fraction $f_V$ and the concentration-independent scattering cross section ${\sigma }_{S,\mathrm{MIE}}$ of particles with volume $V_P$:

Since the specified uncertainty of the refractive index of silica beads is 2%, we have selected the refractive index ($n=1.425$) that describes all measured data most accurately. To study how the refractive index affects the calculations, we have varied the optimal refractive index with 1% (not shown here). This resulted in a maximum uncertainty in ${\mu }_{\mathrm{OCT}}$ of 20% for $0.70\mu \mathrm{m}$ spheres at 0.2 volume fraction (${\mu }_{\mathrm{OCT}}8.9\pm 1.8{\mathrm{mm}}^{-1}$).

An alternative approach to determining ${\mu }_{\mathrm{OCT}}$ would be direct fitting of the EHF model [Eq. (8)] to the OCT signal in depth to obtain ${\mu }_S$ and ${\theta }_{\mathrm{rms}}$. However, this is not always possible because ${\mu }_S$ and ${\theta }_{\mathrm{rms}}$ (via $w_S$) are competing parameters: change in one of the parameters can often be compensated by a change in the other leading to different sets of (${\mu }_S$, ${\theta }_{\mathrm{rms}}$) with equivalent statistical goodness of fit. The resulting curve fit may, therefore, not always converge to realistic values.¹⁵ A promising venue to be explored in future studies is the combination of OCT with other (fiber-based) technologies such as single-fiber reflectance spectroscopy³⁸ that can provide independent parameterization of the scattering phase function—and thus in principle of ${\theta }_{\mathrm{rms}}$.

We calculate the amplitude directly from the OCT data and do not use the amplitude from the fitted A-scans. A practical concern is that in many cases of interest, such as epithelial tissues, the number of data points is limited for an NLLS curve fit. With a common OCT resolution of $10\mu \mathrm{m}$, and for example a layer thickness of $50\mu \mathrm{m}$, only five unique data points are available, which may lead to inconclusive fitting results on ${\mu }_{\mathrm{OCT}}$ and ${\mu }_{B,\mathrm{NA}}$. We scaled the measured OCT signal amplitude to the backscattering coefficient, ${\mu }_{B,\mathrm{NA}}$, within the system’s NA using a scaling factor that was calculated at the highest volume fraction. Ideally, a single calibration factor could be used for scaling.¹⁰ In (clinical) practice, however, the scaling is influenced by parameters that may change between measurements, for example, when switching disposable probes with slightly different coupling efficiencies. In practice, it is, therefore, desirable to calibrate each probe/measurement individually using phantoms with known scattering properties. In this work, it required calibration for each volume-fraction series, per particle diameter.

To obtain quantitative measurement of OCT-derived parameters ${\mu }_{\mathrm{OCT}}$ and ${\mu }_{B,\mathrm{NA}}$, it is crucial to calibrate for system parameters (confocal PSF and sensitivity roll-off). However, the contribution of multiple scattering is challenging to quantify. From Fig. 4, we conclude that multiple scattering starts to have a large ($>10\%$) effect on the measured ${\mu }_{\mathrm{OCT}}$ for samples with a high anisotropy, $g>0.8$, and high scattering coefficient, ${\mu }_S>10{\mathrm{mm}}^{-1}$. In this case, it is advisable to consider multiple scattering to map ${\mu }_{\mathrm{OCT}}$ to ${\mu }_S$. For samples with lower anisotropy factor and scattering cross-section, values of ${\mu }_S$ are within 10% variation of ${\mu }_{\mathrm{OCT}}$.

The role of concentration-dependent scattering is absent for tissue, since this only relates to upscaling optical properties of individual scatterers to those of a solution of scatterers at a given volume fraction. In general, both “an individual scatterer” and “volume fraction” cannot be defined for biological tissue. Instead of a DRM, tissue is arguably better described as a continuous random medium,³⁹ where statistical properties, such as variance, correlation length, and fractal dimension of the spatial refractive index fluctuation, take the role of refractive index contrast, correlation length of Percus–Yevick PCF, and particle size that describes the DRM. Linking these different sets of parameters to OCT measured optical properties will be subject of further research, deploying increasingly hybrid samples (duo and polydisperse bead mixtures of known size distribution).⁴⁰

The exceptional biological specimen for which scattering properties may indeed be described by the Mie–Percus–Yevick formalism, is whole blood. The main scatterers (red blood cells) can be regarded as spherical when sufficient averaging over orientation is possible (validating the use of Mie theory).⁴¹ The physiological volume fraction ($f_V\sim 0.45$) is sufficiently high to warrant the use of the Percus–Yevick structure factor. In a recent theoretical study annex literature review, we found that ${\mu }_s\sim {(1-f_V)}^2(f_V/V_{\mathrm{RBC}}){\sigma }_{\mathrm{SCAT}}$ as opposed to ${\mu }_s=(f_V/V_{\mathrm{RBC}}){\sigma }_{\mathrm{SCAT}}$ for concentration-independent scattering.³⁷ Specifically, the predicted scattering coefficient and scattering anisotropy $g$ at 800 nm are ${\mu }_s=71{\mathrm{mm}}^{-1}$, $g=0.9812$; ${\mu }_A=0.38{\mathrm{mm}}^{-1}$ and ${\mu }_s=63{\mathrm{mm}}^{-1}$, $g=0.9820$; ${\mu }_A=0.47{\mathrm{mm}}^{-1}$ for oxygenated and deoxygenated blood, respectively. The major contribution of multiple forward scattered light in this case resulted in the highly decreased contribution of scattering properties to the total OCT attenuation coefficient, e.g., measured ${\mu }_{\mathrm{OCT}}$ of 5.5 and $5.8{\mathrm{mm}}^{-1}$ for oxygenated and deoxygenated blood, respectively.⁴² Accurate separation of scattering and absorption contributions to ${\mu }_{\mathrm{OCT}}$ may enable quantitative OCT-based measurement of ${\mu }_A$ (spectra) leading to quantification of localized oxygen saturation. Previously, we raised the question of whether quantitative measurements of attenuation by blood coefficients by OCT are feasible⁴³ and concluded that better modeling of the concentration-dependent OCT signal was needed. With the theoretical results presented in Ref. 37 and the experimental validation of the influence of light scattering described in this contribution, accurate determination of attenuation by blood, and subsequent extraction of oxygen saturation, seems feasible and warrants further experimental validation studies.