1. Introduction

A variety of optical techniques are available for characterizing light propagation in diffusive media. Among them, spatially-resolved reflectance measurements have been under significant developments for now over more than two decades and were successfully used in the medical field for characterizing tissues and detecting pathologies [1]. In principle, such techniques consist in locally injecting light into a turbid media and then measuring the diffuse light reflectance profile R(ρ) at the media surface where ρ stands for the injection/detection point separation. A post-treatment of R(ρ) with analytical or numeric models then allows deducing the optical coefficients characterizing light absorption and scattering within the turbid media.

In the past, our laboratory has developed an optical instrumentation measuring spatially-spectrally resolved reflectance profiles R(ρ) at short source-detector separation ρ [2]. Light injection/detection is carried out with an optical probe made of 10 optical fibers packaged in a 2mm-diameter rigid or flexible tube. Such small-dimension optical probes are insertable through the instrumental channel of an endoscope, which has allowed to conduct a clinical in-vivo study of gastritis in the human stomach [3]. Additionally, positioning a small dimension probe at the surface of a tissue, away from major large blood vessels, allows localized optical measurements of optical coefficients within seemingly homogeneous tissue volumes. The ability of targeting limited tissue volumes with an optical probe carrying short source-detector separation can be somehow quantified. Thus, in the case of a steady state diffusion regime in an homogenous diffusive media, it has been established [4, 5] that the mean depth reached by the propagating photons is approximately 0.5(ρ / μ_eff)^0.5 with μ_eff = (3μ_a(μ_a + μ′_s))^0.5. Considering typical μ_a and μ′_s values for biological tissues of 0.1mm⁻¹ and 1mm⁻¹ (λ = 670nm), respectively, and taking the largest source/detector separation ρ for our probes as 1.5mm, we evaluate the mean depth reached by the photons to 0.8mm (photons detected from closer fibers < ρ reach shallower depths). We then conclude that the optical setup previously mentioned carrying R(ρ) at short source-detector separation essentially probes a layer < 1mm immediately below the media surface, thus providing the measurements of optical coefficients within volumes of tissues of ∼ 1mm³. Such spatially selective optical measurements are especially useful for studying local biological processes or histological features. It has allowed measuring optical coefficients changes resulting from epithelial alterations on the mouse back [6]. A last argument supporting the development of probes with such small source-detector geometry is that it would contribute to enhance the sensitivity to high-angle scattering, what may aid to the assessment of tissue morphology [7].

For source/detector separations around one transport mean free path ( ${\mu }_t^{-1}$), measurements [8] and Monte-Carlo simulations [9–11] of R(ρ) profiles are in inadequacy with the values computed using the diffusion approximation of the transport equation. Improved analytical solutions have been proposed to predict the reflectance close to the source in cases of strong absorption [12, 13]. Among them, the P3 approximation approach clearly improved R(ρ) predictions in the short source-detector distance range [14]. Following Ishimaru’s [15] and Groenhuis work [16], Bevilacqua proposed an alternative approach [9] and demonstrated that in geometries where the source-detector distance ρ is in the range 0.5 < ρ μ′_s < 5, a third optical coefficient called γ has to be introduced to accurately describe R(ρ) near the source, according to Wyman’s similarity relations [17]. This phase function coefficient γ can be expressed as:

The γ coefficient, along with a complete theoretical framework on analytical phase functions proposed by Bevilacqua, was found suitable for computing optical coefficients based on R(ρ) measurements at small source-detector distances [18], nevertheless the understanding of the physical/biological meaning of the γ parameter still remains to be clarified. Physically, it is worth remembering that Rayleigh scattering, which is induced by particles much smaller than the light wavelength, is described by the phase function p_Rayl(cosθ) = 3(1 +cos²θ)/4 for non-polarized light. It is then purely of a second order (i.e. quadratic) in cosθ and thus only affects g₂ with g₁ = 0. On the other hand, the scattering due to objects with dimensions close to the wavelength is essentially forward-oriented with typical g₁ > 0.8. Now considering a simple turbid media containing these two types of scatterers and assuming a constant concentration of the larger scatterers, a variation of the concentration of the smaller ones would only affect the 1 – g₂ numerator in Eq. (1). This demonstration case with only two populations of scatterers suggests that γ can be interpreted as the ratio of the relative concentration of small size scatterers compared to the larger ones. However, this simple and intuitive explanation, as seductive as it appears, may not hold for explaining the physical significance of γ when measured on turbid media encountered in real world that generally contain a more complex distribution of scatterers’ sizes. This paper intends to address this issue by studying the dependance of γ(λ) on the fractal distribution of scatterers’sizes.

2. Theory

2.1. Mie’s formalism

For a given wavelength, the Mie formalism provides an analytical expression of the scattering phase function of a sphere noted p_sph(θ,x,n_sph,n_med), where x = πd/λ is the size parameter with d the sphere diameter, n_sph is the real or complex refractive index of the microsphere and n_med the real refractive index of the interstitial fluid [19]. Any phase function p(cosθ) can be profitably expanded in a series of Legendre polynomials P_m(cosθ):

In the present context, this applies to p_sph and allows computing its first and second moments (g_1,sph and g_2,sph). The Mie formalism also provides an analytical expression for the scattering cross-section σ_s,sph(x, n_sph, n_med) [20]. Knowing the concentration of microspheres ρ_sph(d), the scattering coefficient μ_s,sph is then obtained:

The linear correlation between μ_s,sph and σ_s,sph holds at low particles concentrations. At high concentrations near-field optical interactions between neighboring particles become significant and can dramatically impact macroscopic optical properties. The concentrations used in the present study never exceeded 0.1% volume.

2.2. Effective optical coefficients

Let us consider a diffusive medium containing M different types of spheres, each with a different diameter d_i and concentration ρ(d_i). The formalism introduced by Gelebart et al. [21] allows computing the effective phase function moments g_m,eff and scattering coefficient of μ_s,eff(λ) that characterize the scattering encountered by photons propagating through such a media:

Introducing the expressions of g_1,eff and g_2,eff into the definition of the γ coefficient expressed in Eq. (1), one obtains an effective γ_eff for a distribution a sphere diameters:

From μ_s,eff and g_1,eff computed for any discrete distribution N(d_i) of spheres diameters, the expression of the effective reduced scattering coefficient μ′_s,eff is thus given by:

2.3. Fractal size distribution

In tissues, both structural ordering and heterogeneities including cell membranes (∼ 0.01μm), vesicles and mitochondria (∼ 0.1μm) and nuclei and cells (∼ 10μm) contribute to the scattering signal. Consequently, biological tissues contain “scattering centers” with typical sizes covering a large range. Within this framework, Gelebart et al [21] attempted to model phase functions measured from biological tissue with a fractal size distribution of scatterers [22, 23] such that the diameters of scattering centers d are distributed according to:

It then appears that A is proportional to μ_s,eff but has not effect on the g_m,eff and γ_eff values, which means that the anisotropy factor g_1,eff and γ_eff parameters are functions of the fractal dimension α only.

3. Numeric study

3.1. Mie scattering simulations

The scattering phase function of spherical polystyrene microspheres in a water environment was computed according to the Mie theory [19] with a custom C++ routine fully validated by comparison with a web-based calculator [27]. Computations were issued for 500 discrete values of microsphere diameter d_i equally distributed ∈ [0.001, 10]μm. The scattering cross section σ_s,sph (Fig. 1) and phase function p_sph(θ) (not shown) of the polystyrene microspheres were computed for each sphere diameter d_i using complex refractive index values in the spectral range [500, 800]nm measured for latex [28] and water [29]. Using Eq. (3), g_1,sph(d_i), g_2,sph(d_i) were also derived as illustrated in Fig. 2.

 

Fig. 1 Scattering cross section according to Mie’s formalism [19] plotted here in a normalized form σ_s,sph/πd² (right) as a function of the normalized sphere diameter d/λ

Download Full Size | PDF

 

Fig. 2 Expression of g_1,sph (left) and g_2,sph (right) as a function of the normalized sphere diameter d/λ derived from p_sph as expressed in Eq. (3)

Download Full Size | PDF

3.2. Effective optical coefficients

The effective coefficients g_1,eff(λ), g_2,eff(λ) and γ_eff(λ) corresponding to such 500 discrete values of microsphere diameter were computed using Eq. (11), Eq. (12) with the σ_s,sph(d_i), g_1,sph(d_i) and g_2,sph(d_i) calculated using the Mie formalism as inputs (see Fig. 1 and 2). This calculation was carried for 61 values of λ equally distributed ∈ [500, 800]nm and for 101 values of α also equally distributed ∈ [0, 10]. Figure 3 shows the spectral dependence of coefficients g_1,eff and γ_eff.

 

Fig. 3 Numerical expression of g_1,eff (left) and γ_eff (right) as a function of the fractal dimension α as defined in Eq. (11) and (12). The dashed lines indicate the α and γ_eff values corresponding to experimental g₁ values measured on brain and liver [24, 25], thus defining an experiment-based physiological range for these coefficients.

Download Full Size | PDF

Figure 3 indicates that both g_1,eff and γ_eff are insensitive to wavelength changes in the presently simulated spectral range. This theoretical prediction is confirmed through experimental measurement of the γ achieved on biological tissues with spatially-spectrally resolved reflectance measurements [2, 6], with a variation < 20% between 500nm and 900nm.

An extensive review of the optical properties measured ex vivo on various biological tissue/organs and among different animal species was provided by Cheong et al [30], who reported a range of values for the anisotropy factor g₁ between 0.68 (liver) and 0.98 (epicardium).

Under the assumption that a fractal distribution of microspheres’ sizes adequately mimics the optical properties of a biological tissue sample, we can, based on the results of the numerical simulation of g_1,eff(α) presented in Fig. 3-left, deduce that the lower value of the physiological g₁ range (= 0.68) corresponds to α ≅ 4 for λ ∈ [400, 800]nm. This defines a higher bound for the possible α values of biological tissue. Now considering the γ_eff(α) curves in Fig. 3-right, one finds a corresponding γ_eff ≅ 1.5. The same reasoning does not apply to the highest g₁ values (= 0.98) measured on biological samples as g_1,eff reaches a maximal value ≅ 0.95 for α = 0 (Fig. 3-left). Considering Eq. (11), it appears that g_1,eff is a linear combination of g_1,sph as computed with the Mie theory, with linear coefficients < 1. The g_1,sph values for latex microspheres of diameter < 10μm being inferior to 0.95 (Fig. 2-left), it results that the simulated g_1,eff values never reach the upper physiological values of g₁ (= 0.98). In the present simulation context, the upper limit of g_1,sph arise from the large difference between n_med ≅ 1.33 and n_sph ≅ 1.58. Thus, the latex refractive index used for computing the microspheres scattering appears especially high if compared to the values for cell components reported in Table 1.

Table 1. Refractive index of cell components reported in literature

View Table

Therefore, a new simulation of g_1,eff(α) and γ_eff(α) was carried for λ = 650nm (g_1,eff and γ_eff are not wavelength dependant) and five different values of n_sph = 1.58, 1.50, 1.45, 1.40, 1.35 and n_med = 1.33 (see Fig. 4).

 

Fig. 4 Numerical expression of g_1,eff (left) and γ_eff (right) as a function of the fractal dimension α as defined in Eq. (11) and (12). The dashed lines indicate the α and γ_eff values corresponding to experimental g₁ values measured on brain and liver [24, 25], thus defining an experiment-based physiological range for these coefficients.

Download Full Size | PDF

In the left graphic of this figure, it appears that g_1,eff reaches the 0.98 value for n_sph = 1.35 and 1.4. It should be noted here that the 1.38 index of refraction measured on cells [35] falls into that specific n_sph range. This suggests that a fractal distribution of microspheres with a “physiological” index of refraction (instead of the refractive index of latex) is more suited for studying/predicting the optical properties of biological tissue samples. Pursuing our reasoning, two values of α corresponding to g₁ = 0.98 are deduced from the left graph in Fig. 4, equal to 2.3 and 3 (n_sph = 1.4 and 1.35, respectively). This then allows, considering the right graph of Fig. 4, deducing two values of γ_eff = 2.4 and 2.1 that can be interpreted as the higher bound of the γ range that can be expected from biological tissue samples. Following the same reasoning, it appears that the lower bound of the physiological g₁ range (= 0.68) intersect the curves simulated for n = 1.35 and 1.4 at α = 4.5 and 4.6, respectively. These two α values correspond to γ ∼ 1.3. This development can be summarized by saying that numerical simulations of g_1,eff and γ_eff for a fractal distribution of microspheres presented in this paragraph allow predicting, for a range values of “physiological” g₁ ∈ [0.68, 0.98], corresponding ranges for α ∈ [2.3, 4.6] and for γ ∈ [1.3, 2.4]. Additionally, simulations on Fig. 4-right show that α and γ coefficients are nearly linearly correlated within that physiological window, which suggests that γ provides a measure of the fractal power describing the sizes’ distribution of the scattering centers within the diffusive tissue.

Published experimental results support the range predictions resulting from these numerical simulation. Firstly, spatial frequency analysis of phase-contrast images of the mouse liver tissue provided an estimation of the fractal dimension around 3.6 [36], which falls within the predicted range of α values. Secondly, various studies conducted in vivo on biological tissues using spatially-resolved reflectance spectroscopy reported γ values within the wavelength range investigated here (500 to 800nm). Thus, the γ values measured on brain tissue were ∈ [1.7, 2.2] [18], ∈ [1.9, 2.2] for the stomach mucosa [2] and ∈ [1.45, 2.2] for the mouse skin [6]. In details, the study on different brain tissues reported γ values between 1.7 and 1.9 for healthy cortex, optic nerve, white matter and medulloblastoma and 2.2 for cerebellar white matter with scar tissues [18]. In the two other studies, the authors reported significant variations between γ(λ) curves measured on control and altered tissues. [2, 6]. These results also indicate that γ values for healthy tissue fall in the physiological range of γ_eff values computed from the fractal theory ∈ [1.3, 2.4], thus validating the use of a fractal model for describing the concentration of the scattering particles in a healthy tissue. On the contrary, altered tissue (scars, chemical reagent) are often characterized by γ values out of the physiological range, which suggests that fractal distribution varies significantly for altered tissues.

4. Experimental validation

An experimental investigation of the relationship between α and γ was conducted in an attempt to validate the numerical prediction presented above. Four turbid solutions were produced, each corresponding to a distinct α value within the α range [2.5,4]. These solutions were prepared by mixing 6 types of polystyrene microspheres, with distinct diameters (Ø = 0.08, 0.19,0.36, 0.77, 1.03, 9.98μm, Polybead®, Polysciences Inc.). The microspheres were immersed in distilled water and the total volume of the solution set to 7ml. For each solution, the concentrations of the 6 microsphere types were adjusted to comply with the fractal concentration power law for the corresponding α and provide a μ′_s,eff value of 1mm⁻¹ at 650nm. The relative error made on the total volume of these 5 different solutions amounts to 0.5%.

The method used for experimentally measuring the optical coefficients of these diffusive microsphere solutions, namely μ_a,exp(λ), μ′_s,exp(λ) and γ_exp(λ) consisted in immersing a multi-fiber optical probe into the diffusive media and recording the spatially-resolved reflectance with a custom spectrometric instrumentation [2]. Twenty recordings of the spatially-resolved reflectance were carried on each turbid solution. The μ_a,exp(λ), μ′_s,exp(λ) and γ_exp(λ) were then retrieved from the experimental reflectance data following a procedure developed by Thueler [2] (Fig. 5 and Fig. 6).

 

Fig. 5 Optical coefficients μ′_s,eff and μ′_s,exp as numerically simulated (right) and measured with spatially-resolved reflectance spectrometry (left). Data are plotted as a function of λ for α ∈ [2.5, 4]. Relative errors (standard deviation/mean) on μ′_s,exp amount to approximately 4%.

Download Full Size | PDF

 

Fig. 6 Optical coefficients γ_eff and γ_exp (lower row) as numerically simulated (right) and measured with spatially-resolved reflectance spectrometry (left). Data are plotted as a function of λ for α ∈ [2.5, 4]. Relative errors (standard deviation/mean) on γ_exp amount to approximately 6%.

Download Full Size | PDF

The analysis of the reflectance data confirmed that μ_a,exp(λ) is negligible over the visible spectral range in the custom diffusive solutions made of polystyrene microspheres with α ∈ [2.5, 4]. Figure 5 and Fig. 6 show a comparison of the experimental values μ′_s,exp and γ_exp with the μ′_s,eff and γ_eff values obtained from a numerical simulation carried according the fractal approach presented in the previous section and using as inputs the six different sizes’ of microspheres and the selected α values.

In both μ′_s,eff and μ′_s,exp graphs in Fig. 5, the reduced scattering coefficient decays exponentially with the increasing wavelength λ. Such a behavior resembles the wavelength dependence of the scattering cross section predicted by the Rayleigh scattering theory, with σ_s decaying as λ⁻⁴. We compared the μ_eff resulting from our numerical simulations carried with a λ⁻⁴ trend line (Fig. 8) and observed that the slope of μ_eff(λ) curves increases with α but remains inferior to the trend line slope for any α value. This suggests that the higher the fractal dimension α, the closer to the Rayleigh theory the turbid solution scattering behaves. Physically, this hypothesis makes sense as a turbid solutions with a high α contains a large concentration of small diameter scatterers relatively to larger ones, which means that most of the scattering events would result from the interaction of light with small (relatively to λ) particles. Numerical simulations also indicate that the slope of μ′_s,eff decreases with α and this prediction is experimentally confirmed for α ∈ [2.5, 3.5], while μ′_s,exp does not vary significantly when increasing α from 3.5 to 4.

 

Fig. 8 Optical coefficients μ_s,eff computed as a function of λ for a range of α values indicated in the embedded legend. A trend line indicates the λ⁻⁴ behavior expected from Rayleigh’s theory.

Download Full Size | PDF

The numerical simulations carried with 500 discrete values of microsphere diameter d_i (Fig. 3) showed that both g_1,eff (α) and γ_eff(α) do not vary for λ ∈ [500, 800]nm. The γ_eff simulations carried with only 6 values of d_i as well as the experimental measurements of γ_exp, showed moderate (< 5%) oscillations of this parameter around its mean value (Fig. 6), which essentially confirmed this trend. This result is important because it indicates that the behavior of γ relatively to the fractal dimension α is wavelength independent. It also appears looking at Fig. 6 that these moderate oscillations are present in both simulated and experimental results with a local peak around 700 nm. This comforts us with the idea that the effective coefficient γ_eff provides a good estimate of the experimental coefficient γ. For reasons that may be due to the calibration of our instrumentation, the absolute values for μ′_s and γ obtained through computing and experiments are not perfectly matched. Thus, γ_exp for α = 2.5 reaches a local maxima of 1.86 at 690 nm while the corresponding γ_eff is 2.05 (10 % difference). For α = 4, the local maxima at 690 nm amounts to γ_exp = 1.69 and γ_eff = 1.6 (6 % difference). When it comes to the reduced scattering coefficient, we observe that μ′_s,eff values at 500 nm are 20% smaller than μ′_s,exp values for all α. As for γ, one can observe a common feature in both experimental and numerical curves, namely a fast slope increase around 650–700 nm.

Figure 5 shows that the experimental and numerical values μ′_s,exp(λ,α) and μ′_s,eff(λ,α) are highly correlated (r = 0.989). Regarding the γ_exp(λ,α) and γ_eff(λ,α) values, we obtained a correlation r = 0.867. This result demonstrate the validity of the formalism introduced by Gelebart et al. [21] in predicting the effective optical coefficients μ′_s,eff and γ′_eff for a turbid solution with a fractal distribution of microspheres’ sizes.

Based on the results of the numerical simulation presented in the previous section of a diffuse solution with 500 distinct microspheres’ diameter, it was shown that the γ_eff coefficient is linearly correlated to α for a range of α ∈ [2.5, 4] corresponding to a range a measured physiological g₁ values (Fig. 4). The optical reflectance measurements obtained in a custom solution with 6 microsphere diameters between 0.08 and 9.98μm also shows a mean linear correlation r = −0.954 (p = 0.016) between the γ_exp and α for all wavelengths in the [500,800] range (Fig. 7). Furthermore, the γ_eff predicted from a simulation with these 6 microspheres’ diameter is also correlated to α with r = −0.998 (p = 0.0001) (Fig. 7). This experimentally confirms that the γ coefficient deduced from the analysis of diffuse reflectance at short source-detector distance is closely related to the fractal power α governing the corresponding distribution law of scatterers’ sizes.

 

Fig. 7 Measured optical coefficients γ_eff and γ_exp as a function of α for a discrete set of λ ∈ [500, 800]. Relative errors (standard deviation/mean) on γ_exp amount to approximately 6%.

Download Full Size | PDF

Due to the lack of basic understanding regarding the physical interpretation of γ, we believe that carrying fundamental studies involving measurements on known, well controlled, diffusive samples (here the microsphere solutions) was a mandatory step in the process of acquiring the ground knowledge prior to conducting further studies on more sophisticated/realistic media such as biological samples. The present study was focused on simulating μ′_s,eff and measuring μ′_s,exp for diffusive media containing a single fractal distribution of microspheres. However, evidence have been found that some objects encountered in nature appear to be multifractal [37]. In a next step, it would be valuable to investigate the γ properties using a statistical fractal approach, which would require extending the theoretical framework [21] and refine our experimental protocol for the preparation of multifractal diffusive solutions of microspheres.

While a large body of literature supports the idea that tissues are self-organizing systems with fractal dynamics [22, 23] and that the fractal model used in the present study was quite successful at matching tissue phase functions [21, 26], one should remain careful when speculating that the properties of γ demonstrated for microspheres solutions generally apply to more complicated diffusive structure such as biological samples. Thus, the present model does not address optical effects due to specific histologic features, for instance the anisotropic organization of cellular structures resulting in optical birefringence as encountered in the cornea or retina.

5. Conclusion

This paper investigated both numerically and experimentally the physical meaning of the γ parameter introduced by Bevilacqua to improve the calculation accuracy of reflectance profiles close to the light emission source [9]. Considering a fractal distribution of spherical scatterer diameters, the initial step consisted in computing the scattering related coefficients using the Mie formalism for each diameter and then, based on theory, deducing the effective coefficient μ′_s,eff andγ_eff for a set of fractal dimensions over the visible/nearIR spectral range. γ_eff(λ) for such fractal distributions proved to be almost flat, though it is a combination of oscillating functions. Furthermore, a correlation was found between γ_eff and the fractal power α, in the form of a hyperbolic tangent. This correlation was shown to be wavelength and refractive index independent for the range of fractal powers usually found in tissues.

Turbid solutions containing 6 different types of microspheres obeying a fractal size distribution were produced for 4 different fractal dimensions. Spatially-resolved measurements performed on these samples indicated a close correspondence between measured optical coefficients (μ′_s,exp, γ_exp) and numerically computed ones (μ′_s,eff, γ_eff) and also demonstrated a linear correlation between γ_eff and α.

Keeping in mind that these experimental and theoretical investigations were conducted on fractal distributions of spherical scatterers that only partially feature the scattering complexity encountered in biological tissue samples, the obtained results nevertheless suggest that the analysis of γ_exp measurements conducted on biological turbid media may provide relevant information to the physician in term of cells and sub-cellular components’ distribution and their potential pathologic alterations.