Monte Carlo model of the depolarization of backscattered linearly polarized light in the sub-diffusion regime

We present a predictive model of the depolarization ratio of backscattered linearly polarized light from spatially continuous refractive index media that is applicable to the sub-diffusion regime of light scattering. Using Monte Carlo simulations, we derived a simple relationship between the depolarization ratio and both the sample optical properties and illumination-collection geometry. Our model was validated on tissue simulating phantoms and found to be in good agreement. We further show the utility of this model by demonstrating its use for measuring the depolarization length from biological tissue in vivo. We expect our results to aid in the interpretation of the depolarization ratio from sub-diffusive reflectance measurements. ©2014 Optical Society of America OCIS codes: (170.6935) Tissue characterization; (290.5855) Scattering, polarization. References and links 1. N. Ghosh and I. A. Vitkin, “Tissue polarimetry: concepts, challenges, applications, and outlook,” J. Biomed. Opt. 16(11), 110801 (2011). 2. E. Akkermans, P. E. Wolf, R. Maynard, and G. Maret, “Theoretical study of the coherent backscattering of light by disordered media,” J. Phys. France 49, 77–98 (1988). 3. L. F. Rojas-Ochoa, D. Lacoste, R. Lenke, P. Schurtenberger, and F. Scheffold, “Depolarization of backscattered linearly polarized light,” J. Opt. Soc. Am. A 21(9), 1799–1804 (2004). 4. M. Xu and R. R. Alfano, “Light depolarization by tissue and phantoms,” Proc. SPIE 60840, 60840T (2006). 5. J. D. Rogers, I. R. Capoğlu, and V. Backman, “Nonscalar elastic light scattering from continuous random media in the Born approximation,” Opt. Lett. 34(12), 1891–1893 (2009). 6. A. J. Gomes, S. Ruderman, M. DelaCruz, R. K. Wali, H. K. Roy, and V. Backman, “In vivo measurement of the shape of the tissue-refractive-index correlation function and its application to detection of colorectal field carcinogenesis,” J. Biomed. Opt. 17(4), 047005 (2012). 7. A. Radosevich, J. Rogers, V. Turzhitsky, N. Mutyal, J. Yi, H. Roy, and V. Backman, “Polarized enhanced backscattering spectroscopy for characterization of biological tissues at subdiffusion length-scales,” IEEE J. Sel. Top. Quantum Electron. 18(4), 1313–1325 (2011). 8. M. Moscoso, J. B. Keller, and G. Papanicolaou, “Depolarization and blurring of optical images by biological tissue,” J. Opt. Soc. Am. A 18(4), 948–960 (2001). 9. J. M. Schmitt and G. Kumar, “Turbulent nature of refractive-index variations in biological tissue,” Opt. Lett. 21(16), 1310–1312 (1996). 10. M. Xu and R. R. Alfano, “Fractal mechanisms of light scattering in biological tissue and cells,” Opt. Lett. 30(22), 3051–3053 (2005). 11. V. Turzhitsky, A. Radosevich, J. D. Rogers, A. Taflove, and V. Backman, “A predictive model of backscattering at subdiffusion length scales,” Biomed. Opt. Express 1(3), 1034–1046 (2010). 12. M. Hunter, V. Backman, G. Popescu, M. Kalashnikov, C. W. Boone, A. Wax, V. Gopal, K. Badizadegan, G. D. Stoner, and M. S. Feld, “Tissue self-affinity and polarized light scattering in the born approximation: a new model for precancer detection,” Phys. Rev. Lett. 97(13), 138102 (2006). 13. C. J. Sheppard, “Fractal model of light scattering in biological tissue and cells,” Opt. Lett. 32(2), 142–144 (2007). 14. O. Nadiarnykh, R. B. LaComb, M. A. Brewer, and P. J. Campagnola, “Alterations of the extracellular matrix in ovarian cancer studied by Second Harmonic Generation imaging microscopy,” BMC Cancer 10(1), 94 (2010). #176090 $15.00 USD Received 12 Sep 2012; revised 27 Oct 2012; accepted 15 Nov 2012; published 28 Feb 2014 (C) 2014 OSA 10 March 2014 | Vol. 22, No. 5 | DOI:10.1364/OE.22.005325 | OPTICS EXPRESS 5325 15. V. M. Turzhitsky, A. J. Gomes, Y. L. Kim, Y. Liu, A. Kromine, J. D. Rogers, M. Jameel, H. K. Roy, and V. Backman, “Measuring mucosal blood supply in vivo with a polarization-gating probe,” Appl. Opt. 47(32), 6046– 6057 (2008). 16. D. Bicout, C. Brosseau, A. S. Martinez, and J. M. Schmitt, “Depolarization of multiply scattered waves by spherical diffusers: Influence of the size parameter,” Phys. Rev. E Stat. Phys. Plasmas Fluids Relat. Interdiscip. Topics 49(2), 1767–1770 (1994). 17. M. Xu and R. R. Alfano, “Random walk of polarized light in turbid media,” Phys. Rev. Lett. 95(21), 213901 (2005). 18. Y. Liu, Y. Kim, X. Li, and V. Backman, “Investigation of depth selectivity of polarization gating for tissue characterization,” Opt. Express 13(2), 601–611 (2005). 19. X. Guo, M. F. G. Wood, N. Ghosh, and I. A. Vitkin, “Depolarization of light in turbid media: a scattering event resolved Monte Carlo study,” Appl. Opt. 49(2), 153–162 (2010). 20. M. Ahmad, S. Alali, A. Kim, M. F. Wood, M. Ikram, and I. A. Vitkin, “Do different turbid media with matched bulk optical properties also exhibit similar polarization properties?” Biomed. Opt. Express 2(12), 3248–3258 (2011). 21. N. Ghosh, H. Patel, and P. Gupta, “Depolarization of light in tissue phantoms effect of a distribution in the size of scatterers,” Opt. Express 11(18), 2198–2205 (2003). 22. V. Sankaran, M. J. Everett, D. J. Maitland, and J. T. Walsh, Jr., “Comparison of polarized-light propagation in biological tissue and phantoms,” Opt. Lett. 24(15), 1044–1046 (1999). 23. S. A. Prahl, M. J. van Gemert, and A. J. Welch, “Determining the optical properties of turbid media by using the adding-doubling method,” Appl. Opt. 32(4), 559–568 (1993). 24. R. Michels, F. Foschum, and A. Kienle, “Optical properties of fat emulsions,” Opt. Express 16(8), 5907–5925 (2008). 25. A. J. Gomes and V. Backman, “Analytical light reflectance models for overlapping illumination and collection area geometries,” Appl. Opt. 51(33), 8013–8021 (2012). 26. I. Georgakoudi, B. C. Jacobson, J. Van Dam, V. Backman, M. B. Wallace, M. G. Müller, Q. Zhang, K. Badizadegan, D. Sun, G. A. Thomas, L. T. Perelman, and M. S. Feld, “Fluorescence, reflectance, and lightscattering spectroscopy for evaluating dysplasia in patients with Barrett’s esophagus,” Gastroenterology 120(7), 1620–1629 (2001). 27. A. N. Bashkatov, E. A. Genina, V. I. Kochubey, A. A. Gavrilova, S. V. Kapralov, V. A. Grishaev, and V. V. Tuchin, “Optical properties of human stomach mucosa in the spectral range from 400 to 2000 nm art. no. 673401,” Proc. Soc. Photo-Opt. Ins. 6734, 73401 (2007). 28. C. Holmer, K. S. Lehmann, J. Wanken, C. Reissfelder, A. Roggan, G. Mueller, H. J. Buhr, and J. P. Ritz, “Optical properties of adenocarcinoma and squamous cell carcinoma of the gastroesophageal junction,” J. Biomed. Opt. 12(1), 014025 (2007). 29. A. Amelink, H. J. Sterenborg, M. P. Bard, and S. A. Burgers, “In vivo measurement of the local optical properties of tissue by use of differential path-length spectroscopy,” Opt. Lett. 29(10), 1087–1089 (2004). 30. A. Kim, M. Roy, F. Dadani, and B. C. Wilson, “A fiber optic reflectance probe with multiple source-collector separations to increase the dynamic range of derived tissue optical absorption and scattering coefficients,” Opt. Express 18(6), 5580–5594 (2010). 31. R. Reif, M. S. Amorosino, K. W. Calabro, O. A’Amar, S. K. Singh, and I. J. Bigio, “Analysis of changes in reflectance measurements on biological tissues subjected to different probe pressures,” J. Biomed. Opt. 13(1), 010502 (2008). 32. A. J. Gomes, V. Turzhitsky, S. Ruderman, and V. Backman, “Monte Carlo model of the penetration depth for polarization gating spectroscopy: influence of illumination-collection geometry and sample optical properties,” Appl. Opt. 51(20), 4627–4637 (2012). 33. M. J. Everett, K. Schoenenberger, B. W. Colston, Jr., and L. B. Da Silva, “Birefringence characterization of biological tissue by use of optical coherence tomography,” Opt. Lett. 23(3), 228–230 (1998). 34. A. J. Radosevich, J. D. Rogers, I. R. Capoğlu, N. N. Mutyal, P. Pradhan, and V. Backman, “Open source software for electric field Monte Carlo simulation of coherent backscattering in biological media containing birefringence,” J. Biomed. Opt. 17(11), 115001 (2012).


Introduction
Polarized light scattering methods have found numerous biomedical applications ranging from early cancer detection to glucose sensing owing to their ability to both suppress multiple scattering and obtain intrinsic tissue polarization metrics [1].Relating the depolarization of light to sample properties is thus an important area of research.Strides in this area have been made in the case of the transmission and reflection of diffuse polarized light.Beginning with Akermans et al [2], an expression for the depolarization ratio was derived for the case of completely isotropic scattering (anisotropic factor g = 0).The validity of the expression was further extended across the range of anisotropy values by the work of Rojas-Ochoa et al [3].
Further refinements led to the understanding of diffuse polarized transmission through a slab of fractal continuous media [4].
Less research has focused on the sub-diffusive regime of light scattering which is an important area in biophotonics.For many biomedical applications, it is desirable to target superficial tissue structures.Many precancerous dysplastic lesions arise from the epithelium and mucosae of several organs (breast, lung, colon, etc.).The thickness of the mucosal layer is typically on the order of several hundred microns.Targeting these shallow layers often necessitates the design of fiber-optic probes or systems that have small illumination/collection areas, minimal separation between the source and collection areas, and/or angled illumination/collection beams.Clinical considerations, such as the desire for small probes to fit through hospital instrument channels, also dictate probe design.For these types of system geometries, the diffusion approximation does not apply.Consequently, much previous work that studied depolarization using the diffusion approximation are not applicable to probes or systems that specifically target small tissue volumes dominated by the collection of subdiffusive photons.In order to use depolarization to characterize turbid media from subdiffusive light scattering measurements, it is then necessary to formulate an alternative framework to model the effect of the optical properties and system geometry on the reflectance depolarization ratio.
In this paper, we empirically develop a simple analytical model that relates the reflectance depolarization ratio to the sample transport mean free path and the depolarization length scale, as well as to aspects of the illumination and collection geometry such as the collection angle and size of the illumination/collection area.Our model was first created using polarization-sensitive Monte Carlo simulations that were validated against experimental phantoms.An important aspect of these simulations is that they employed a phase function based on the Whittle-Matérn model [5] of light scattering from media with spatially continuous refractive index fluctuations such as those found in biological tissue.We tested our model on a tissue simulating Intralipid phantom and found that our model could accurately measure the depolarization length of the phantom.We further applied our model to in vivo polarization-gated fiber-optic probe measurements from human esophageal tissue and demonstrate that the depolarization ratio, the transport mean free path, and the depolarization length are all altered with dysplasia.The work presented in this paper should aid in the use of the depolarization ratio of sub-diffusive reflectance to characterize the intrinsic optical properties of turbid media.

Polarization-sensitive Monte Carlo simulations
The basis for our Monte Carlo (MC) simulations has been described in detail elsewhere [6,7].In brief, the simulations employ a Whittle-Matérn phase function for continuous refractive index media.Models of continuous refractive index media characterized by weak index fluctuations have proven useful in explaining the power spectra, phase functions, and depolarization properties both biological cells and tissues [8][9][10].The shape of the phase function depends on the shape parameter (m) of the refractive index correlation function, and the correlation length l c of the correlation function [11].The correlation function shape can approach a Gaussian (m → ∞), a decaying exponential (m = 2), a stretched exponential (1.5 < m < 2), or a power law (m < 1.5) in which case the medium has a mass fractal distribution with mass fractal dimension D mf = 2m.Expected values of m from biological tissue range from 1 to 2 [5].Previous light scattering measurements have characterized some tissues as having fractal structural organization implying that m < 1.5 in these cases [9,10,12,13].Ovarian tissue has been found to have a measured m between 1.5 and 2 implying a stretched exponential correlation function [14].The value of kl c , where k is the wavenumber, is expected to be larger than one in tissue [5].The anisotropy factor (g), which is defined as the average cosine of the scattering angle, can be expressed as a function of m and kl c .Detailed expressions for g have been given in [5].The illumination and collection of geometry of our MC simulations were configured to mimic a polarization-gated fiber-optic probe we have used in previous studies [15].The probe delivers collimated light within a circular area to the medium and collects backscattered light that emerges from within the illumination area.For this type of probe, the geometry can be characterized by two parameters: the radius (R) of the illumination-collection area and the angle (θ c ) at which collected photons exit the medium relative to the surface normal.The Stokes vectors of the collected photons were summed and Mueller matrix multiplication was used to calculate the intensity of reflectance with polarization parallel ( || I ) to the incident beam and the intensity of reflectance with polarization perpendicular ( I ⊥ ) to the incident beam.We tracked the full intensity depolarization ratio (d) previously defined as [3] || || ,

Depolarization ratio model
The starting point for our model development was based on previous research concerning the diffuse transmission of polarized light through a slab.In these studies, the depolarization ratio is a function of a length scale ξ which itself is a function of the depolarization length (l p ) and the transport mean free path (l t ): ξ = 0.5 [4,16].We assume that d for the sub-diffusive reflectance case should also be a function of l p and l t .Increased l p would be expected to increase d due to increased polarization preservation while increased l t would be expected to increase d due to less depolarizing scattering events.In the Whittle-Matérn model for continuous refractive index, the properties of the refractive index correlation function summarized by the variance of the refractive index dn 2 , the correlation length l c , and the shape parameter m determine the scattering properties of the medium.Explicit relationships between the index correlation parameters and scattering properties of the medium such as l t are given in [5].To determine the dependence of l p on the scattering properties of the medium under the Whittle-Matérn model, we first use a result previously derived from a random walk of vector photons [17]: In these expressions S 1 and S 2 are elements of the amplitude scattering matrix which for soft tissue satisfy S 1 = S and S 2 = S cos θ.The amplitude scattering matrix is related to the spectral density Φ via [10]: where k is the wavenumber and θ is the scattering angle.The relationship between l p and S under an alternative fractal-continuous model has been studied by Xu and Alfano [4].In the Whittle-Matérn model of light scattering, the spectral density at spatial frequency κ is defined as [5]: where N is an normalization constant.Examination of Eqs. ( 2)-( 4) reveals that l p /l s will be a function of kl c and m.The quantity l p /l s , since it is normalized by l s , depends on the phase function for a single scattering event.The Whittle-Matérn phase function can also be written as a function of kl c and m [11].In addition to the optical properties, d should also be a function of the illumination/collection radius R and collection angle θ c .It is known that the fraction of co-polarized photons collected will be higher as R decreases [18] and as θ c goes closer to exact backscattering [19].With these considerations in mind, we developed the following expression for the depolarization ratio in the sub-diffusion regime using leastsquares surface fitting: where f 1 , f 2 , and f 3 are linear functions of θ c and their values are given in Table 1.The functions f 1 -f 3 envelope the θ c dependence of d.The R/l t term in Eq. ( 5) captures the effect of both illumination-collection area and the effective scattering power of the sample on d while the R/l p term describes the effect of both illumination-collection area and the phase function (embodied in the l p /l s dependence) on d.It is important to stress that Eq. ( 5) assumes that the illumination and collection areas completely overlap on the sample surface.
The intrinsic properties of the sample that compose Eq. ( 5) are the transport mean free path l t and the characteristic length scale of depolarization l p .These properties in turn depend on the properties of the refractive index correlation function as summarized by dn 2 , l c , and m.Explicit relationships between the index correlation parameters and scattering properties of the medium under the Whittle-Matérn model are given by Rogers et al. [5].To determine how the changes in the correlation function parameters influence the terms in Eq. ( 5) and the resulting light depolarization two things need to be calculated: 1.)The effect of the index correlation on the scattering parameters which can be determined from [5] and 2.) The effect of the index correlation function on the depolarization length scale which can be determined from Eqs. ( 2)-(4).As an example, for a fixed kl c = 10, an increase in m from 1.6 to 1.9 is expected to decrease l t (decrease d due to increased scattering) according to Fig. 5 of [5] while l p /l t is decreased (decrease d) according to Eqs. ( 2)-( 4).The decrease in the depolarization length combines with the decrease in l t such that d is decreased with increasing m.We performed MC simulations as described in Section 2.1 and calculated d from each simulation according to Eq. ( 1).We then fit Eq. ( 5) to the MC data using a least squares approach to determine the coefficients f 1 , f 2 , and f 3 .
There are two important aspects of Eq. ( 5).The first is that it is applicable to systems having continuous refractive index profiles and also potentially applicable to systems having discrete refractive index profiles.The difference will involve in how the theoretical value l p is calculated for a given sample.For example, Mie theory can be used in the case of discrete spherical particles to calculate the amplitude scattering matrix used in Eqs. ( 2)-( 4) while Whittle-Matérn theory is used for continuous media.The second important consequence of Eq. ( 5) is that it implicitly predicts that samples with matched bulk optical properties (l t ) but different phase functions (and hence l p 's) will have different depolarization properties.This is in agreement with experimental evidence showing that tissue and phantoms depolarize light differently and that the size distribution of scatterers affects light depolarization [20][21][22].

Probe instrumentation
Depolarization ratio measurements from phantoms and human esophageal tissue were taken with a polarization-gated spectroscopy probe that has been characterized and described in detail previously [15].In brief, the probe consists of three 200 μm diameter fibers (NA = 0.22): an illumination fiber and two collection fibers.Linear polarizer placement allows the simultaneous collection of co-polarized ( I  ) and cross-polarized ( I ⊥ ) reflected light relative to the incident linearly polarized beam.A gradient refractive index (GRIN) lens collimates the incident light (half-angle divergence of 3 degrees) and focuses light backscattered from the sample onto the collections fibers.A white light-emitting diode (WT&T) is used for illumination while fiber optic spectrometers (Ocean Optics) processes light from the collection fibers as a function of wavelength.The radius (R) of the illumination/collection areas on the sample surface for the current probe is 475 μm and the angle (θ c ) between illumination and collection beams is ~14°.

Experimental phantoms
We carried out several experiments on microsphere phantoms with the polarization-gated probe and compared the depolarization ratios observed with the probe with results from Mie based MC simulations of the microsphere phantom and probe geometry.The purpose of these experiments was to determine if the MC simulations could accurately predict the depolarization ratio observed by the probe.We prepared phantoms consisting of polystyrene microspheres (Thermo-Scientific) of different sizes and concentrations in deionized water.The optical properties of the prepared phantoms were calculated using Mie theory and the phantom properties are summarized in Table 2.A series of 20 individual spectral measurements were taken with the fiber-optic probe immersed in the microsphere solution.The solution was illuminated with light from a white light-emitting diode (WT&T).Reflected light was collected by two Ocean Optics fiber-optic spectrometers.These measurements were background subtracted and normalized to a diffuse reflectance standard.The signal intensity at a wavelength of 589 nm was used to calculate the depolarization from each microsphere solution.These results were then directly compared with Mie-based MC simulations of the exact phantom optical properties and probe geometry.From these experiments, we calculated a correction coefficient to apply to Eq. ( 5) that takes into account incomplete overlap of illumination-collection areas of the probe as well as the below 100% contrast polarizers used by the probe.The correction coefficient is described in Section 3.4.
Next, we took spectral measurements with the polarization-gated probe on a 20% Intralipid (Sigma-Aldrich) solution.The purpose of this experiment was to test if our depolarization model in Eq. ( 5) could be utilized to predict the depolarization length of the sample.Equation ( 5) was derived from MC simulations using a Whittle-Matérn model of light scattering from continuous refractive index media.The refractive index variations from Intralipid cannot be considered continuous.However, Eq. ( 5) suggests that the depolarization ratio will only depend on l t and l p .These values can also be calculated for discrete systems [4].We determined the l t spectrum of the Intralipid phantom using integrating sphere measurements coupled with the adding-doubling algorithm [23,25].From Eqs. ( 2) and (3), l p depends on the elements of the amplitude scattering matrix and l s .The parameter l s can be determined from the l t measurements using the anisotropy factor (g) of 20% Intralipid which has been previously measured [24].The elements of the amplitude scattering matrix for 20% Intralipid can be computed by superimposing the Mie theory results for each particle size in the Intralipid particle size distribution which has also been provided in [24].The measured d value from the probe and the l t values of the phantom were used to calculate experimentally measured l p values using Eq. ( 5) which were then compared with the theoretical l p values for 20% Intralipid.

Measurement of depolarization length with polarization-gating
From Eq. ( 5), the depolarization length l p can be extracted if both d and l t can be measured.The polarization-gated probe directly measures d while l t can be quantified with the crosspolarized reflectance signal.The cross-polarized reflectance signal is advantageous because its signal intensity has no added relationship to the phase function than what is incorporated into l t [6].The physical basis for this observation is that since it takes several scattering events to depolarize the incident light, the cross-polarized signal originates from multiple scattering even when the illumination-collection radius is small.Using this fact, we have previously shown that spectral measurement of the cross-polarized intensity can be used to measure both l t and m using a least-squares fitting method [25].In brief, the reflectance intensity of the cross-polarized signal can be written as where c is a calibration constant determined experimentally on a phantom with known optical properties, f 1 and f 2 are functions of m and θ c , a is a fitting a parameter, a μ is the absorption coefficient, and L ⊥ is the mean effective path length whose functional dependence on optical properties and illumination-collection geometry has been given previously [25].In the Whittle-Matérn model, the quantity 2 4 m aλ − determines the value of ( ) kl >> which is satisfied for most biological tissue [5].By fitting Eq. ( 6) to the collected spectral intensity of the cross-polarized signal using a and m as free parameters, both m and l t can be measured.

Human studies
The study was conducted at Mayo Clinic Florida in, Jacksonville FL and approved by the Institutional Review Board.Data were collected between 2010 and 2012 from a total of 63 patients undergoing upper GI endoscopy.Barrett's esophagus (BE) surveillance was the primary reason for endoscopy though some patients were evaluated for anemia, dyspepsia, or heartburn.After informed consent, patients underwent upper endoscopy with conscious sedation using midazolam, fentanyl and meperidine or propofol.After the initial diagnostic evaluation but prior to any tissue biopsy sampling, approximately 10 readings, each taking 15 ms, were taken with a polarization-gated fiber optic probe passed down the endoscope accessory channel.The probe was placed in gentle contact with the gut mucosa surface in BE regions within the distal esophageal mucosa.If a dysplastic nodule was seen, probe measurements were also taken from the nodule itself.The co-polarized and cross-polarized reflectance intensities were analyzed at a wavelength of 630 nm to mitigate the effect of hemoglobin absorption.The depolarization ratios of tissue sites from the same region were averaged together to form a single d for that region and patient.Of the 63 patients analyzed, 6 had dysplastic nodules from which probe measurements were taken.

Dependence of depolarization ratio on sample optical properties
The main optical properties that compose the d model in Eq. ( 5) are l p and l t .We examine the effect of m and g on l p /l s as shown in Fig. 1.In Fig. 1, we plot l p /l s as a function of m for different g values.When g is fixed and m is allowed to vary or vice versa, the value of kl c is implicitly also changed since m, kl c , and g are functionally related [5].For a fixed m value, a higher g leads to an increased l p /l s value.This means that it takes more scattering events to depolarize the incident linearly polarized light in a high g medium when m is fixed.This is consistent with the observation that small-angle scattering events (high g) depolarize the incident light less than small-scattering events (low g) because there is less randomization of the electric field vector direction.In addition, Fig. 1 demonstrates that m also modulates l p , particularly for higher g values.The increased l p at lower m values can be explained by the fact that for a fixed g, a phase function with lower m has higher phase function intensity in the backward direction [11].This characteristic has been shown to maintain light polarization in the case of discrete Mie particles [19].Figure 2(a) illustrates the relationship between d and m for different g values and a fixed R/l s = 8.It is observed that d in Fig. 2(a) follows the same trend as l p /l s in Fig. 1.This is expected as d of the collected light will depend on the ability of the light to maintain polarization in the medium as quantified by the l p /l s parameter.Finally, in Fig. 2(b) we show the dependence of d on the dimensionless parameter R/l s .As l s increases for a fixed R, there is less scattering to depolarize the incident light leading to the expected increase in d.

Dependence of depolarization ratio on illumination-collection geometry
Two elements of the illumination-collection geometry that will be considered here are the size of the illumination-collection radius (R) on the sample surface and the collection angle θ c relative to the exact backscattering direction.Due to the scaling property of Monte Carlo, the d can be plotted versus the dimensionless parameter R/l s as was done in Fig. 2(b).Thus, the effect of R on d for a fixed l s can also be inferred from Fig. 2 is maximal as R → 0 and approaches 0 as R → ∞ [18].Next, we consider the effect of θ c on d in Fig. 3.The parameter d decreases as θ c increases, a finding in accordance with the observation that large angle scattering events are more depolarizing than small angle scattering events [19].

Monte Carlo model for the depolarization ratio
Equation ( 5) was fit to the data generated by MC simulations for a variety of optical properties and for several R and θ c .An example fit is shown in Fig. 4(a) for θ c = 0-18°.The individual black data points are MC outcomes while the shaded surface is the best fit to the data points (r 2 > 0.99).The functional relationship between [f 1 f 2 f 3 ] and θ c was found to be linear of the form x 1 + x 2 θ c and these values are shown in Table 1.Using the values of these parameters, we computed d generated by Eq. ( 5

Experimental verification of the depolarization ratio model
We measured d from microsphere phantoms with the polarization-gated probe and compared it to d calculated from Mie-based MC simulations of the corresponding phantom.The purpose of this was two-fold.The first objective was to demonstrate that our MC simulations were predictive of the results observed experimentally with the probe.The second objective was to quantify how the non-ideal polarization behavior of a realistic probe would cause deviations from Eq. ( 5).Equation ( 5) was derived under the assumptions of ideal polarizers and completely overlapping illumination-collection areas.In our probe, we have observed that the polarizer and GRIN contrast is on the order of 94% while there is a small 30 µm center-tocenter separation between the illumination-collection areas.To examine these effects, we conducted two Mie-based MC simulations: one assuming ideal polarization and complete overlap and the other simulating polarizers with 94% contrast and 30 µm center-to-center separations.The experimental comparison with these two MC simulations is summarized in Table 2.Each row corresponds to a separate microsphere phantom and the columns give the microsphere diameter of the phantom, the g and ' experimentally with the probe, and the d's calculated using ideal and non-ideal MC simulations.Table 2 demonstrates that there is an appreciable difference between the d computed by the ideal and non-ideal MC simulations.The probe measurement of d closely matches the non-ideal MC simulation with an average percent error of <5%.This demonstrates that the MC model of the probe geometry accurately predicts the depolarization observed experimentally by the probe.
The percent error between the probe experiment d value and the ideal MC d value was ~12%.These results suggest that an error term needs to be introduced into Eq.( 5) to account for the non-ideal behavior of realistic fiber-optic probes.To determine this error term we calculated the d from Whittle-Matérn based MC simulations for ideal and non-ideal cases across a wide range of optical properties described in Section 2.1.The probe geometry was held constant with R = 475 µm and θ c = 0-18°.The result of this comparison is shown in Fig. 5 where it seen that there is a consistent proportional relationship (r 2 > 0.99) between the d calculated from the ideal and non-ideal simulations.For this specific probe geometry, Eq. ( 5) can then be modified as For probes with other R, θ c , or center-to-center separation the scaling factor may be different than 0.89.The scaling factor should be calculated independently for each probe design.
Next we sought to determine whether our model for the depolarization ratio in Eq. ( 7) could be used to measure the depolarization length of the Intralipid phantom.Equation ( 7) can be inverted to find l p using both the measured depolarization ratio from the probe and the l t values of the phantom.The l t values can be determined using an integrating sphere measurement coupled with the adding-doubling algorithm or by analysis of the crosspolarized reflectance signal as described in Section 2.5.In Fig. 6 we present the l p /l t spectrum of the Intralipid sample.The l t spectrum was measured with an integrating sphere and the theoretical l p spectrum of the Intralipid was calculated from Eqs. ( 2)-(4) as well as the known particle size distribution and g spectrum of Intralipid as described in Section 2.4.Using the known l t spectrum and the measured depolarization ratio as a function of wavelength, we computed an experimentally measured l p at each wavelength.These experimentally measured l p values are compared with the theoretical l p values in Fig. 7(a).There is a strong linear correlation between the experimentally derived and theoretical l p values with a linear correlation coefficient greater than 0.99.We further quantified the percent error between the measured and theoretical l p values as illustrated in Fig. 7 the theoretical l p by an average of 7%.The error magnitude decreases for higher l p which is consistent with Fig. 4(c) showing a lower error for higher d (and hence higher l p ). Sources of discrepancy include noise in the polarization-gated measurement, possible error in the measurement of l t , as well as uncertainty regarding the precise value of the scaling factor in Eq. (7).Despite these sources of error, Fig. 7 demonstrates good correlation between the measured and theoretical l p values.

Application of the depolarization model to biological tissue
Having shown in the previous section that our model of the depolarization ratio can be used to accurately measure the depolarization length of a tissue simulating phantom, we next took polarization-gated measurements from biological tissue in vivo as an illustrative example of applying the model.Polarization-gated probe measurements were taken from regions of the esophageal mucosa with and without dysplasia.We found that the depolarization ratio was significantly altered in dysplastic sites as shown in Fig. 8(a).The depolarization ratio increases from a mean value of 0.56 in non-dysplastic regions to a mean value of 0.7 in dysplastic ones with a statistically significant p-value of 0.005.Figure 8(b) illustrates that this elevation of d with dysplasia is driven by both an increase in l t (from 0.08 to 0.12 cm, p-value = 0.02) and l p (from 0.12 to 0.19 cm, p-value = 0.02).The value of l t was determined from the cross-polarized signal as described in Section 2.5, while l p was determined from Eq. ( 7).Our observations of an increase in l t with dysplasia in Barrett's esophagus are consistent with previous studies using diffuse reflectance spectroscopy [26].Encoded in the ratio p t l l is information about the scattering phase function.This is apparent from Eqs. ( 2)-( 4) that show p s l l depends on the variables m and kl c which are the parameters that govern the shape of the Whittle-Matérn phase function [11].Furthermore, we arrive at p t l l by multiplying p s l l by (1-g) where g is average cosine of the scattering angle.The value g can be expressed in terms of m and kl c such that p t l l is entirely a function of m and kl c [5].Since we measure m from the cross-polarized reflectance signal (from Eq. ( 6)), the value of kl c can be determined from the measurement p t l l .This would allow us to reconstruct the phase function shape and determine the anisotropy factor from each tissue site.However, the measured g is very sensitive to the precise magnitude of p t l l .
For example, for a fixed m = 1.6, a shift of 10% in p t l l from 1.42 to 1.29 can alter the measured g from 0.8 to 0.95.Thus very precise and accurate measurements of p t l l are required and due to inherent limitations on signal-to-noise ratio, accuracy of the l p , l t , and m measurements, accuracy of the depolarization ratio model, as well as tissue inhomogeneity and variability, we do not believe that we can obtain a robust measurement of g from each tissue site at this time.In an attempt to minimize these sources of error, we took the average values of all the depolarization ratios and m values measured from each tissue site in our dataset.These average ± standard deviation values were m = 1.7 ± 0.09 and d = 0.57 ± 0.06.Using these values in conjunction with of the value of l t for esophagus, the parameter p t l l may be determined.The value of l t for esophagus at ~ 630 nm varies in the literature from ~ 0.06-0.13cm with our own measured value 0.08 cm falling in between this range [26,28].The resulting range of p t l l would be from 1.19 -1.67.These values correspond to klc = 1.26 -30.9 and g = 0.58 -0.99.It is encouraging that these values are within the expected range for the upper gastrointestinal tract [27,28].The p t l l ratio could potentially be used to measure the tissue anisotropy factor from each tissue site but the wide range of g values obtainable from the l t range suggests that the ratio must be estimated precisely.
In our depolarization model, for the sake of simplicity, we have not considered tissue inhomogeneity or other polarization effects such as birefringence.Our model is based on the assumption, common in the biomedical optics literature, that the scattering properties are distributed homogeneously within the medium [29][30][31].Tissue, however, is a complex layered and inhomogeneous structure.The depolarization signal from our probe, for example, will sample both the epithelium and stroma of the esophageal tissue which are expected to have different optical properties.The optical properties we measure will have contributions from both layers.This can be mitigated by designing a probe whose penetration depth for the cross-polarized signal (which has the highest depth) does not exceed the depth of the top layer.This can be achieved by decreasing R and/or increasing θ c [32].Finally, tissue can be expected to have linear birefringence due to structural alignment of collagen fibers, lipid bilayers, muscle fibers, etc.To test whether this property could influence the depolarization ratio, we performed and compared two MC simulations: one where the sample had no birefringence and another where the sample possessed random birefringence (the case applicable to esophageal tissue where collagen fibers are oriented randomly) having a difference in ordinary and extraordinary refractive indices of 0.001 which is the maximum difference expected in tissue [33].The method for the MC simulation with birefringence followed that given in [34].We found that the depolarization ratio varied less than 4.2% between these two simulations indicating that the birefringence will have a minimal effect on the observed depolarization ratio.This is consistent with the finding in [34] that birefringence had minimal impact on the spatial reflectance profiles for linearly polarized light.

Conclusion
The methods outlined in this paper provide a model to predict the depolarization of backscattered linearly polarized light for different sample optical properties as well as different illumination-collection geometries.The model can be applied to probes and systems that collect sub-diffusive photons.Our model was validated on tissue simulating phantoms and we applied our depolarization model to calculate the depolarization length from biological tissue in vivo using a simple polarization-gated measurement.We expect that our depolarization model will help in the assessment of tissue optical properties using polarization-gating spectroscopy.Future work will further explore the effect of tissue inhomogeneity on the observed depolarization ratio as well as obtaining site-specific measurements of the anisotropy factor.The full range of applicability of the depolarization model to discrete systems is also of interest.

Fig. 1 .
Fig. 1.The scaled depolarization length as a function of the shape parameter (m) of the refractive index correlation function for different anisotropy (g) values.

Fig. 2 .
Fig. 2. The relationship between the depolarization ratio (d) and the optical properties of the medium.(a) d as a function of the parameters of the phase function m and g for R/l s = 8 and θ c = 0-18°.(b) The relationship between d and the dimensionless parameter R/l s for m = 1.3, g = 0.9 and θ c = 0-18°.

Fig. 3 .
Fig. 3. Relationship between the depolarization ratio from Monte Carlo (d MC ) and the collection angle around the exact backscattering direction.The other constants were R/l s = 8, m = 1.3, and g = 0.9.

Figure 2 (
b) shows that d decreases as R increases.This can be understood by studying the radial probability distribution (P(r)) of light exiting the sample surface.Previous studies have shown that ) (d Model ) and compared it with d from the actual MC calculations (d MC ).In Fig. 4(b), d MC is plotted versus d Model with the result showing a strong linear correlation (R 2 > 0.99).The percentage error between d MC and d Model is shown in Fig. 4(c) as a function of d MC .The mean percentage error was 8 percent across the entire range of d MC while it was 3 percent in the more biologically relevant regime of d MC > 0.5.

Fig. 4 .
Fig. 4. (a) The depolarization ratio as a function of R/l t and R/l p from MC simulations (filled circles).The displayed surface is the best fit of Eq. (5) to the MC data for a θ c of 14°.(b) The depolarization ratio from MC simulations (d MC ) compared with the depolarization ratio as computed by the model in Eq. (5) (d Model ).(c) The percent error between d MC and d Model as a function of d MC .

Figure 7
(b).The percent error is defined as (b) shows that the experimental l p values are generally lower than

Fig. 5 .
Fig. 5. Depolarization ratio (d) calculated from MC simulations incorporating ideal polarizers and complete illumination-collection area overlap (d ideal ) versus the depolarization ratio as calculated from MC simulations incorporating polarizers with 94% contrast and a 30 µm center-to-center separation between illumination and collection areas (d non-ideal ).

Fig. 6 .
Fig. 6.Calculated ratio of the Intralipid depolarization length to the transport mean free path as a function of wavelength.

Fig. 8 .
Fig. 8. (a) The depolarization ratio (d) is significantly elevated in dysplastic regions of the esophagus.(b) The elevation in d is driven by increases in both l t and l p with dysplasia.Intervals are means ± 95% confidence intervals.