1.

Introduction

Light that has traveled through tissue carries information about tissue properties, such as structure and biochemical composition. The distance that light travels determines the scale on which the tissue is investigated, as optical properties are averaged over the sampling volume. To detect localized, small-scale changes, techniques with small sampling volumes are, therefore, appropriate. Consequently, short source–detector separations are required to ensure most detected photons will have scattered only once or a few times. In this letter, we will investigate the reflectance measured in this so-called subdiffusive regime, specifically for overlapping source–detector geometries, relevant to, e.g., single fiber spectroscopy measurements.

In the diffusive regime, photon direction is randomized and reflectance is independent of the exact shape of the phase function [$p(\theta )$, the probability distribution of scattering angles]. In contrast, in the subdiffusive regime, the photon direction is not fully randomized; therefore, measurements are sensitive to the phase function.^1,2 Advances have been made to derive analytic models for subdiffusive reflectance.^3–5 However, for overlapping source–detector geometries, challenges remain regarding incorporation of the tissue phase function.^3,4

To model light transport, solutions to the radiative transport equation (RTE) involve expanding the radiance into a series of $i$ spherical harmonics and the phase function into $i$ Legendre polynomials. The latter are weighted with their moments $g_i$, where $g_1$ is commonly referred to as the scattering anisotropy. For $i=1$, the diffusion approximation to the RTE is obtained, with the similarity relation ${\mu }_s(1-g_1)={\mu }_s^{*}(1-g_1^{*})$, which expresses that tissues with a different scattering coefficient ${\mu }_s$ and a different scattering anisotropy $g_1$, but equal ${\mu }_s^{\prime }={\mu }_s(1-g_1)$ yield the same reflectance. However, reflectance at short source–detector separations is inadequately described for $i=1$. Improvement is possible by increasing $i$, giving additional similarity relations. For $i=2$, Bevilacqua and Depeursinge⁶ suggested an alternative form of these similarity relations using the parameter $\gamma$

We consider subdiffusive reflectance as the sum of a diffusive component and a semiballistic component: the diffusive component depends on ${\mu }_s^{\prime }{d}_{\det }$, and, for the semiballistic component, only photons that have experienced a single backscattering event in combination with an arbitrary number of small-angle forward scattering events are considered. To model this contribution, not all details of the phase function are needed. Rather, we consider the phase function within the NA of the detector. The NA characterizes the range of angles over which a detector can accept incoming photons. We propose a theoretically derived parameter ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, which is related to the NA of the detector and the integral of the phase function over the corresponding angular interval in the backward ($p_{\mathit{NA},\text{backward}}$) and forward ($p_{\mathit{NA},\text{forward}}$) directions, as

2.

Theoretical Derivation of ${\mathit{R}}_{\mathit{p}\mathit{NA}}$

We consider semiballistic photons that are backscattered once in combination with an arbitrary number of forward directed scattering events and assume that these photons will only be detected if subsequent scattering events occur at angles smaller or equal to the acceptance angle ${\theta }_{NA}$ (Fig. 1), where ${\theta }_{\mathrm{NA}}=\arcsin (\mathrm{NA}/n_{\text{sample}})$ and $n_{\text{sample}}$ is the refractive index of the sample. We introduce the probabilities $p_{\mathit{NA},\text{backward}}$ and $p_{\mathit{NA},\text{forward}}$ as

Fig. 1

Semiballistic scattering: a combination of any number of forward scattering events $N$ in combination with a single backscatter event. We assume photons will only be detected if all scattering angles $\theta$ are smaller than or equal to the acceptance angle ${\theta }_{\mathrm{NA}}$ of the detector. The scattering angle $\theta$ is determined with respect to the previous scattering direction as indicated by the dashed lines.

Summing Eq. (6) over all possible numbers of forward scattering events, $N$ from 0 to $\infty$, we obtain

Thus, the combined probability of one backscatter event and an arbitrary number of forward scattering events, each occurring within ${\theta }_{\mathrm{NA}}$, is given by integrating over all possible pathlengths, 0 to $\infty$

3.

Methods

For our Monte Carlo simulations, we modified the software of Prahl et al.¹⁵ (which was the core programming later used in Monte Carlo model of steady-state light transport in multilayered tissues (MCML)¹⁶ software), to allow the use of arbitrary phase functions using the method of Zijp and ten Bosch.¹⁷ The modified code was benchmarked against standard MCML using the Henyey–Greenstein phase function. We modeled geometries where the source and detector geometry were equal and overlapping, where photons were launched from a location based on a uniform distribution across the source with an angle from a uniform angular distribution within the NA. We performed simulations using 15 modified Henyey–Greenstein (mHG),⁶ 144 double Henyey–Greenstein (dbHG),¹⁸ 8 modified power of cosines (MPC),⁶ and 46 Reynolds McCormick¹⁹ phase functions (RMC, which is equivalent to the Gegenbauer kernel phase function) employing the parameters specified in Table 1 and applying the restrictions $g_1\ge 0.5$ and $g_2<0.9$ to exclude biologically unreasonable phase functions.

Table 1

Parameters employed in the selection of phase functions.

For each set of phase functions, simulations were performed for three values of ${\mu }_s^{\prime }{d}_{\det }:0.1$ (${\mu }_s^{\prime }=10{\mathrm{cm}}^{-1}$, $d_{\det }=100\mu \mathrm{m}$), 1 (${\mu }_s^{\prime }=100{\mathrm{cm}}^{-1}$, $d_{\det }=100\mu \mathrm{m}$), and 9 (${\mu }_s^{\prime }=100{\mathrm{cm}}^{-1}$, $d_{\det }=900\mu \mathrm{m}$) in combination with an NA of 0.22 or 0.5. An absorption coefficient of $0.1{\mathrm{cm}}^{-1}$ and refractive indices inside and outside the sample of 1.35 and 1, respectively, were used.

4.

Results

Figure 2 shows the simulated reflectance versus ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, $\sigma$, and $\gamma$ for an NA of 0.22 and 0.5, and three different values of ${\mu }_s^{\prime }{d}_{\det }$ (0.1, 1, and 9). Reflectance correlates with all three parameters for ${\mu }_s^{\prime }{d}_{\det }$ of 0.1 and 1. This correlation is far less pronounced for ${\mu }_s^{\prime }{d}_{\det }=9$. To compare ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, $\sigma$, and $\gamma$, we determined the spread in each of these three parameters for a chosen reflectance ($\pm 10\%$) relative to the total range of the parameter (Table 2). For example, for ${\mu }_s^{\prime }{d}_{\det }=0.1$, NA = 0.22 and a reflectance of 0.001 ($\pm 10\%$), $\gamma$ ranges from 1.43 to 1.66. The total range of $\gamma$ is 0.68 to 2.43, so the variability is calculated as 0.23/1.75 = 0.13. For all combinations of NA (0.22 and 0.5) and ${\mu }_s^{\prime }{d}_{\det }$ values (0.1 and 1), we determined these variability values for three reflectance values and found the variability to be lowest for ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ compared to $\sigma$ and $\gamma$.

Fig. 2

Simulated reflectance versus ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, $\sigma$, and $\gamma$ for (a)–(c) $\mathrm{NA}=0.22$ and (d)–(f) $\mathrm{NA}=0.5$. Symbols indicate ${\mu }_s^{\prime }{d}_{\det }$ values, and colors indicate phase function types. Note the log scales for both the reflectance and ${\mathit{R}}_{\mathit{p}\mathit{NA}}$.

Table 2

Variability of RpNA, σ, and γ for μs′ddet=0.1 and μs′ddet=1, defined as the spread in RpNA, σ, and γ values for a chosen reflectance (±10%) relative to the total range of each parameter.

5.

Discussion

We have theoretically derived the parameter ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ to model subdiffusive light scattering for overlapping source–detector geometries, and we have shown with Monte Carlo simulations that the reflectance depends on ${\mathit{R}}_{\mathit{p}\mathit{NA}}$. In comparison to $\sigma$ and $\gamma$, ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ improves prediction of the reflectance. Our findings indicate that the reflectance does not depend on the details of the phase function—knowledge of all Legendre moments is not required—but only on the magnitude of the phase function within the acceptance angle of the detector (in the backward and forward directions).

Currently, the parameter $\gamma$ is widely used to model subdiffusive reflectance. Since ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ improves prediction of the measured reflectance, incorporating ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ into subdiffusive models is expected to improve their reliability and thereby the estimation of other optical properties, such as ${\mu }_s^{\prime }$. This improvement is significant because subdiffusive measurements have the potential to detect small-scale tissue changes.

For higher values of ${\mu }_s^{\prime }{d}_{\det }$, the reflectance becomes more diffusive and, therefore, depends less on ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, $\sigma$, and $\gamma$. The contribution of diffuse photons is a remaining challenge for describing subdiffusive scattering. In our simulations, we kept the absorption coefficient constant. For subdiffusive reflectance, path lengths are short; therefore, the effect of ${\mu }_a$ will be minor. While the absorption affects the diffuse contribution to the reflectance, it has a minor effect on semiballistic photons. To implement ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ for subdiffusive measurements, models have to be developed incorporating ${\mathit{R}}_{\mathit{p}\mathit{NA}}$, the diffusive reflectance, and absorption for specific measurement geometries, such as single fiber reflectance spectroscopy.

6.

Conclusion

In conclusion, the theoretically derived parameter ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ predicts subdiffusive light scattering for overlapping source–detector geometries. Consequently, the reflectance does not depend on the details of the entire phase function, but on the phase function within the acceptance angles of the detector. Since ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ improves prediction of the measured reflectance compared to $\sigma$ and $\gamma$, the use of ${\mathit{R}}_{\mathit{p}\mathit{NA}}$ is expected to improve derivation of optical properties from subdiffusive measurements.