Statistics of the depth probed by cw measurements of photons in a turbid medium

Photon migration in a turbid medium has been modeled in many different ways. The motivation for such modeling is based on technology that can be used to probe potentially diagnostic optical properties of biological tissue. Surprisingly, one of the more effective models is also one of the simplest. It is based on statistical properties of a nearest-neighbor lattice random walk. Here we develop a theory allowing one to calculate the number of visits by a photon to a given depth, if it is eventually detected at an absorbing surface. This mimics cw measurements made on biological tissue and is directed towards characterizing the depth reached by photons injected at the surface. Our development of the theory uses formalism based on the theory of a continuous-time random walk~CTRW!. Formally exact results are given in the Fourier-Laplace domain, which, in turn, are used to generate approximations for parameters of physical interest. @S1063-651X~98!15411-9#


I. INTRODUCTION
Many forms of optical technology are now being explored for potential use in medical diagnosis, because they do not involve the use of potentially hazardous ionizing radiation, ͓1,2͔.As an example, Bonner and Nossal proved the feasibility of measuring the number density and rms speed of moving red cells in vivo by optical methods ͓3͔.A somewhat similar technology has also been investigated with engineering applications in mind ͓4͔.Since light is multiply scattered and can be absorbed in a tissue, it is necessary to have a theoretical basis with which to translate experimental data into optical parameters potentially capable of distinguishing between healthy and diseased states.The most frequently used of these are s Ј , the transport-corrected scattering co- efficient and a , the absorption coefficient.Many types of theories have been applied to reducing data from optical measurements on turbid media.These include detailed multiple-scattering formalism ͓5͔, a variety of forms of transport and diffusion theory ͓6͔, and most recently the theory of lattice random walks ͓7,8͔.
All biomedical applications that make use of optical methods require the capability of characterizing and controlling photon trajectories, insofar as is possible.This is obviously desirable in optical imaging applications in which one looks for regions that may have anomalous scattering or absorption properties ͓9-11͔, but it is also true in applications involving homogeneous media.In many instances it is important to characterize the depth to which photons can penetrate into a medium so as to associate an optical response to underlying physiology and to remove possible artifacts due to uninteresting inhomogeneities.One descriptive parameter that relates to the statistical properties of a photon trajectory is the maximum depth to which it penetrates the medium of interest.This was first discussed in ͓7͔.A related, but somewhat more involved, calculation is that of the average depth sampled by a photon in a cw measurement ͓12͔.Aspects of this theory were experimentally confirmed in ͓13͔.The calculations in both ͓7͔ and ͓12͔ were based on random walk theory but they could as well have been found using diffusion theory.
The theory developed here deals with experiments in which photons are injected into tissue at a point, some of which are later reemitted at the tissue interface at another point.Three basic types of experiments are envisioned here: ͑1͒ The cw experiment, in which a continuous beam of photons enters the tissue, and the intensity of reemitted photons is measured as a function of distance from the entrance point.͑2͒ Time-gated experiments, in which injecting and detecting optodes are kept at a fixed distance apart and the intensity of reemitted photons is measured as a function of time.͑3͒ Frequency modulated experiments, in which a response is measured as a function of frequency to periodically injected photons.In each of these the intensity of reemitted light is measured, from which one can derive estimates of the optical parameters s Ј and a .Each of the three types of experi- ments is potentially capable of providing information related to optical parameters hidden in bulk tissue, and each has both practical advantages and disadvantages.This paper deals only with the simplest possible case of tissues with homogeneous optical properties.More interesting problems are posed by tissues having different forms of heterogeneity, particularly in light of possible imaging applications.
The average depth calculated in ͓12͔ for the cw experiment is equivalent to the expected value of a random variable known as the occupancy of a random walk ͓14͔.This, in turn, is related to the mathematical notion of local time ͓15͔, which has been used in the context of optical imaging by a number of other investigators ͓16,17͔.The average value of a random variable can only crudely characterize its properties.This suggests the utility of developing a more detailed calculation able to more accurately quantitate information related to photon trajectories.
In this paper we derive a formally exact theory for the distribution of the number of visits to a given depth in a semi-infinite medium, rather than just its first moment, potentially allowing us to calculate higher moments for such systems.A similar calculation for transillumination experi-ments ͑i.e., transmission through slabs͒ is described elsewhere ͓18͔.Here we consider only the case of cw measurements made on a semi-infinite medium, which is also the model treated in ͓7͔.The extension that allows us to deal with time-gated experiments will be left for future discussion, but is essentially based on the mathematical analysis developed in the present paper.
In contrast to the analysis in ͓12͔, which assumed that steps were made at equally spaced intervals in time, here we make use of the continuous-time random walk model ͑CTRW͒ ͓19,14͔, which has been shown to simplify a number of calculations required to elucidate the theory of photon migration ͓20͔.In the CTRW model the number of steps taken at any time is random, as will be seen from the detailed analysis to follow.The present analysis will be restricted to the case of the cw experiment as in ͓12͔ but can also be extended to analyze the other experiments mentioned earlier.Optical techniques as a whole are widely used to measure meteorological parameters, and similar problems arise in such applications ͓21͔.

II. DESCRIPTION OF THE MODEL
A schematic picture of the model to be analyzed is shown in Fig. 1.Photons are injected into tissue by a laser beam, idealized as being a line.Since typical penetration depths in tissue are of the order of millimeters, the medium being probed is modeled as being a semi-infinite bulk and the skin, or interface, is assumed to be planar.In the course of their migration photons are randomly scattered by inhomogeneities within the medium, exemplified by cytoplasmic organelles and cell surfaces.Some of the photons may be absorbed internally, the remaining ones reaching the interface where they appear as emitted light.The measurable light intensity considered as a function of distance along the surface can be shown to contain information about optical properties of the underlying tissue which is potentially useful for diagnostic purposes ͓7͔.
The analysis that follows is based on the following assumptions.
͑1͒ The internal structure of the medium is replaced by a simple cubic lattice.A site on this lattice will be denoted by the vector rϭ(x,y,z), where the components of the vector are integers.While these coordinates are convenient for the mathematical analysis, they can later be converted to dimensional units by r ¯ϭr&/ s Ј where s Ј is the transport- corrected scattering coefficient.Positive values of z correspond to points interior to the medium, the surface of the medium being zϭ0.The coordinates transverse to the z axis will be labelled x and y.The ranges of the three coordinates are 0рzрϱ and Ϫϱрx, yрϱ.
͑2͒ The isotropic random walk that is a model for a photon trajectory allows for motion to nearest-neighboring sites only, i.e., on a given step the random walker moves from ͑x,y,z͒ to one of the sites (xϮ1,yϮ1,zϮ1), the particular site being chosen with probability 1/6.This picture cannot account for preferential forward scattering, which is important at very short times, but the isotropic model has been shown to be adequate for many experimental applications ͓8͔.Corrections to account for anisotropic scattering can be accounted for either rigorously using a full transport theory or possibly by using a more phenomenological theory based on the telegrapher's equation rather than the diffusion equation.This was originally suggested by Ishimaru ͓22͔, and investigated more recently and in greater detail by Durian and Rudnick ͓23͔.
͑3͒ The time between successive steps is a random variable.The probability density that describes this time, designated by (t), is a negative exponential: where the rate constant k is related to the ͑assumed͒ constant speed of light in the medium, c, and s Ј by kϭc s Ј .In the ensuing analysis we will make use of the dimensionless time ϭkt, which is equivalent to setting kϭ1.͑4͒ Absorption in the medium is governed by Beer's law.That is to say, the probability that a random walker survives inside the lattice for a dimensionless time without being absorbed is equal to exp(Ϫ) where is expressible in terms of the experimentally measurable absorption coefficient a as the ratio ϭ a / s Ј .

͑5͒
The planar interface zϭ0 is comprised only of absorbing points.A photon reaching the surface from within the tissue is immediately removed from the system, and at the arrival instant, is detectable as reemitted light.Because of this property the entry point, designated by r 0 , is at ͑0,0,1͒.This is the equivalent to the assumption customarily made in diffusion analyses, that considers the initial positions of photons to be approximately a single scattering length below the surface in the medium ͓24͔.
We will be interested in finding the distribution function for the total number of times the random walker has visited level zϭZ conditional on its reaching the surface at a point (x,y,0)ϭ(,0) at time where ϭ(x,y).Our aim is to calculate the probability that the depth Z has been visited exactly k times conditional on the random walker reaching the surface at ͑,0͒ at the ͑dimensionless͒ time .This conditional probability will be denoted by v k (Z͉,).The joint probability for visiting zϭZ exactly k times and reaching the surface at ͑,0͒ at time will be denoted by v k (Z,,).The conditional probability v k (Z͉,) can be converted to information about how much time has been spent by the photon trajectory at a given depth below the surface.

A. Representation of the propagators
As background to the forthcoming calculations we need to introduce several propagators, where by a propagator we mean the probability that the random walker is at r at time starting from a specified point.The propagator for a nearestneighbor random walk on a translationally invariant, i.e., unbounded, simple cubic lattice will be denoted by g 0 (r,).This is the probability that a random walker initially at 0 is at r at time .Taking account the possibility of internal absorption this propagator can be shown to have the form ͓20͔ By using the asymptotic expansion of I x (/3) for →ϱ it can be shown that g 0 (r,) can be approximated by the Gaussian function where r 2 ϭx 2 ϩy 2 ϩz 2 .This equation will be valid for distances such that r is of the order of 1/2 or less.By keeping higher-order terms in expansions of the Bessel functions one can establish that the Gaussian form in this last equation breaks down at times for which rϭO().The breakdown of the Gaussian approximation is to be expected because the photons whose distance increases linearly with time are ballistic photons.Since the Gaussian structure of Eq. ͑3.2͒ is a consequence of the central-limit theorem we expect this sort of a breakdown to occur in the tails of the curve.Because of the lattice structure, Eq. ͑3.1͒ is only symmetric with respect to interchanges in x and y or to changes in sign of these variables.The limiting behavior in Eq. ͑3.2͒ shows that radial symmetry only appears at times long enough so that the number of ballistic photons is negligible.Numerical calculations that suggest the validity of the assumption of radial symmetry are to be found in ͓20͔.
When the absorbing boundary at zϭ0 and the initial position of the random walker r 0 ϭ(0,0,1) are taken into account Eq. ͑3.1͒ is replaced by Since the aim of our analysis is to derive information related to visits to the level zϭZ, we will also need an expression for the propagator when there are absorbing boundaries both at zϭ0 and zϭZ.This has been shown to be equal to ͓20͔.The propagators in Eqs.͑3.1͒, ͑3.3͒, and ͑3.4͒ are exact for the CTRW model, in contrast with the Gaussian approximation used in much of the earlier work described in ͓8͔.
The subscripts that are used to distinguish between the propagators in Eqs.͑3.1͒, ͑3.3͒, and ͑3.4͒ refer to the number of boundaries.Because the product I x (/3)I y (/3) will appear so frequently in later calculations we henceforth replace the product by a function to be denoted by I(;) defined by

͑3.5͒
This function can be shown to be approximately Gaussian at long times by using the same argument that leads to Eq. ͑3.2͒.

B. Probability of never reaching the depth Z
The simplest problem asks for the probability that the maximum depth is less than Z conditional on reaching ͑,0͒ at time .We will, in fact, calculate the probability that the photon never reaches the depth Z, conditional on reaching ͑,0͒ at time .This will be denoted by v 0 (Z͉,).This probability can be found from the joint probabilityprobability density that the random walker does not visit z ϭZ and is later absorbed at ͑,0͒ during the time interval (,ϩd).It is assumed in all the following calculations that the random walker's initial position is ͑0,1͒, which allows us to drop the initial position as an argument.The required relation between the conditional and joint probabilityprobability density is since v 0 (ϱ,,)d is just the probability that the random walker is absorbed at ͑,0͒ between and ϩd in a semiinfinite medium.To exclude the possibility that depth Z is never reached by photons that are later detected, we put an absorbing boundary at that depth and use g 1 (r,) ͓Eq.͑3.3͔͒ as the propagator to calculate the denominator and g 2 (r,;Z) ͓Eq.͑3.4͔͒ to calculate the numerator of Eq. ͑3.6͒.
The probability density v 0 (ϱ,,) is just which follows by reasoning that in order for the random walker to reach ͑,0͒ during (,ϩd) it must have reached ͑,1͒ at time р, paused for a time Ϫ, then made a single step to ͑,0͒ with a probability equal to 1/6.During that time it must not have been absorbed internally.The numerator in Eq. ͑3.6͒ takes a slightly more complicated form because the propagator g 2 (r,) replaces the function I 1 (/3)/ in the integral in Eq. ͑3.7͒.One finds that which, together with Eq. ͑3.7͒, provides an explicit expression for the conditional probability v 0 (Z͉,) in the time domain.The integral can be evaluated and leads to an expression for v 0 (Z͉,) that has the form

͑3.9͒
Notice that this is independent of the absorption parameter , which follows from our conditioning the random walk to be absorbed at a point on the surface rather than internally.An analysis of cw experiments requires finding not v 0 (Z,,) but rather the photon intensity integrated over all time, ͐ 0 ϱ v 0 (Z,,)d.
While the evaluation of Eq. ͑3.9͒ requires that the term in the denominator be integrated numerically, it is possible to derive a relatively simple expression for the total intensity of light emerging from the surface at time that never reached the depth Z.This is This can be simplified by making use of the identity which, together with Eqs.͑3.9͒ and ͑3.10͒, yields 2ϩcos͑ j/Z ͒ .

͑3.12͒
The associated conditional probability that the random walker has not reached the depth Z conditional on reaching the surface at time is An interesting consequence of the representation in Eq. ͑3.9͒ is that after a short time v 0 (Z͉,) becomes almost independent of so that for all practical purposes v 0 (Z͉,) Ϸv 0 (Z͉).To see why this should be so we observe that for fixed and large one can approximate the product I(;/3) by which is seen to be independent of .Since a Bessel function of the form I x (/3) increases nearly exponentially in at large enough values of this parameter we can, for example, write the denominator of Eq. ͑3.9͒ as where can be any fixed parameter less than .Since both the numerator and denominator of Eq. ͑3.9͒ contain similar integrals we conclude that v 0 (Z͉,) is almost independent of .Detailed numerical calculations indicate that for ZϾ3 the error due to the assumption of isotropy is generally less than 5%.Figure 2 shows curves of v 0 (Z͉) plotted as a function of for Zϭ5 and Zϭ7, which reinforces the obvious point that as the observation time increases, the probability that a random walker has not reached a given level before reaching the surface must decrease.

The probability of reaching Z k times
kϭ1.In the remainder of this section we derive an expression for a joint Fourier-Laplace transform of v k (Z,,), where v k (Z,,)d is the probability that the depth Z has been reached by the random walker exactly k times before it is absorbed at ͑,0͒ during the time interval (,ϩd).We will show that the joint Fourier-Laplace transform of v k (Z,,) can be related to that of the expected number of visits to level Z.It is instructive, although not strictly necessary, to first analyze the case of kϭ1 because one form of the expression for the joint Fourier-Laplace transform of v k (Z,,) for kϾ1 is proportional to the joint transform of v 1 (Z,,).
An expression for v 1 (Z,,) is obtained by accounting for a sequence of five events, which we write symbolically as ͑ 0,1͒→͑Ј,ZϪ1͒→͑Ј,Z͒→͑Ј,ZϪ1͒→͑,1͒→͑,0͒, ͑3.16͒ where Ј can be any point in the ͑x,y͒ plane.It should be understood that once the photon reaches zϭZ it never penetrates to any greater depth nor does it move around in the plane zϭZ.The second, third, and last transition in Eq. ͑3.16͒ each occur with a probability of 1/6, the total time consumed in making these transitions being described by a probability density 3 (), which is calculated in terms of the probability density for the time of a single step, ͑͒.The relationship is most succinctly written in terms of the Laplace transform of (t), which will be denoted by ˆ(s) ϭ1/(1ϩs).The required relation is ˆ3(s)ϭ ˆ3(s) so that the transitions (Ј,ZϪ1)→(Ј,Z)→(Ј,ZϪ1) and (,1) →(,0) are accounted for in the Laplace transform domain by multiplying by a term ˆ3(s)/216.
These considerations yield the following expression for v 1 (Z,,): Observe that the lattice is translationally invariant in any plane perpendicular to z and that the random walk is isotropic in any plane parallel to the surface.This allows us to invoke the transformation which can be used to transform Eq. ͑3.18͒ into a spatial convolution integral.The delta function constraint is equivalent to a convolution in the domain.
In order to further simplify Eq. ͑3.18͒ we introduce the joint Fourier-Laplace transform of v 1 (Z,,) defined by v ˆ1͑Z,,s ͒ϭ ͚ By the transform g ˆ2(,ZϪ1,s) we mean F L ͓g 2 (,Z Ϫ1,͉0,1)͔ where F ͕•͖ denotes the two-dimensional Fou- rier transform with the transform parameter vector and L ͕•͖ is a Laplace transform with a transform parameter s.
Notice that our use of a Beer's law approximation is readily incorporated into the Laplace transform formalism by replacing s on the right-hand side by sϩ.This means that calculations can be carried for the case of no internal absorption without loss of generality and the factor e Ϫ inserted at the conclusion of any calculation.

͑3.21͒
This transform can be evaluated exactly since it is possible to find an exact expression for g ˆ2 (,ZϪ1,s).The calculation of the two-dimensional Fourier transform makes use of the identity and the subsequent Laplace transform is readily found.Following this set of substitutions we find that g ˆ2͑,ZϪ1, sϩ ͒ϭ 2 in which cϭcos 1 ϩcos 2 .If this relation is substituted into Eq.͑3.21͒ the resulting set of transforms can be inverted exactly, leading to an expression for v 1 (Z,,) as where

͑3.25͒
The conditional probability v 1 (Z͉,) is found by dividing Eq. ͑3.24͒ by Eq. ͑3.7͒.Because the factor appears only because of the factors I(;/3) in the numerator and denominator of the resulting expression we might expect that the dependence on is a relatively weak one.This turns out to be correct provided that у20.To test this conjecture we compared the cases ϭ(3,4) with ͑5,0͒ finding that at ϭ10 the relative discrepancy amounted to 6.2%, at ϭ20 it amounted to 1.5%, and at ϭ30 it amounted to 0.4%.A similar comparison for ͑6,8͒ and ͑10,0͒ showed that ϭ10 the relative discrepancy was 17.2%, at ϭ20 it amounted to 2.3% and at ϭ30 it amounted to 0.3%.This suggests that for у30 it is reasonable to approximate v 1 (Z͉,) by v 1 (Z͉,)Ϸv 1 (Z͉) with an insignificant error except at very large values of the distance .kϾ1.To treat the case in which zϭZ is visited exactly k(Ͼ1) times it will be necessary to work in terms of transforms.The details of doing so are presented in the subsequent analysis.In the case considered here the cycle in Eq. ͑3.16͒ is replaced by where the j are arbitrary points in the ͑x,y͒ plane.Notice that any transition of the form ( l ,Z)→( lϩ1 ,Z) allows the random walk to either travel at levels above or below Z except for the first and last visits to that plane.We must therefore account for components of the cycle of the form ( i ,Z)→( iϩ1 ,Z).Since sites in any ͑x,y͒ plane are translationally invariant we can, without loss of generality, choose i ϭ0 and, for convenience of notation set iϩ1 ϭЈ.(Ј,Z͉0,Z), jϭ1,2,...,kϪ1͖ where F ( j) (Ј,Z͉0,Z)d is the probability that, starting from (0, Z), a random walker reaches zϭZ for the jth time at a point Ј in the x,y plane at some time during the interval (,ϩd).The jth passage probabilities can be found by recursion from where we have invoked the property of translational invariance in the plane to pass from the first line of this equation to the second.Since the recursion relation is a convolution in both space and time, it assumes a simpler form in terms of joint Fourier-Laplace transforms defined analogously to Eq. ͑3.20͒.Let F ˆs ( j) (,Z͉0,Z) denote the transform of F ( j) (Ј,Z͉0,Z).The recursion step in Eq. ͑3.27͒ is then equivalent to or which is required for calculating the transform of v k (Z,,).The cycle scheme in Eq. ͑3.26͒ can be decomposed into a sum of three contributions allowing us to abbreviate the sequence as C 1 →C 2 →C 3 .In this symbolic scheme C 1 takes the photon from its initial position to the plane zϭZ for the first time.The propagator corresponding to C 2 is so that, according to Eq. ͑3.29͒, the equivalent in the transform domain is just ͓F ˆsϩ (1) (,Z͉0, Z)͔ kϪ1 .The Fourier-Laplace transform of the propagator for C 1 is just g ˆ2(,Z Ϫ1,sϩ͉0,1) ˆ(sϩ)/6 and the one corresponding to C 3 is g ˆ2(,ZϪ1,sϩ͉0,1) 2 (sϩ)/36.The product of these two factors is equal to v ˆ1(Z,,s) as can be seen from Eq. ͑3.21͒ so that v ˆ1(Z,,s) can be expressed as It is awkward to work with the function v ˆ1 but this expression can be reformulated in terms of a more convenient relation in terms of the joint transform of the expected number of visits to the level Z.This quantity will be designated E ˆ(Z,,s) and can be written as so that Eq. ͑3.31͒ can be rewritten as Thus, establishing the dependence of v ˆk(Z,,s) on k for general k requires being able to find expressions for, and suitable approximations to, the transforms F ˆsϩ (1) (,Z͉0,Z) and E ˆ(Z,,sϩ).A further consequence of Eq. ͑3.33͒ is that we can express the joint transform of the second moment, S(Z,,sϩ), in terms of E ˆ (Z,,sϩ).By summing the geometric series that comes from Eq. ͑3.32͒ we find

͑3.34͒
A derivation in detail of the representation for E ˆ(Z,,s ϩ) is given in ͓25͔.Expressions for the transforms of higher-order moments can also be expressed just in terms of the functions that appear in Eq. ͑3.34͒ because of the geometric form of v ˆk(Z,,sϩ) in Eq. ͑3.31͒.
C. An evaluation of F ˆs "1… ",Zͦ0,Z… The Fourier-Laplace transform F ˆs (1) (,Z͉0,Z) can be calculated in terms of the transform of the propagator g 1 (,Z,͉0,Z) by invoking the relation where we have made use of the translational invariance that holds in any ͑x,y͒ plane.On taking the joint transform of the identity in Eq. ͑3.35͒ and solving for F ˆs (1) (,Z͉0,Z) we find which is the final building block required for the evaluation of Eq. ͑3.31͒.An appeal to the explicit formula for g 1 (,Z,͉0,Z) in Eq. ͑3.3͒ as well as the identity in Eq. ͑3.22͒ yields an expression for g ˆ1(,Z,s͉0,Z), which is in which the parameter a sϩ is where, as before, cϭcos 1 ϩcos 2 .With these exact results in hand, we can return to the problem of evaluating v ˆk(Z,,sϩ) for the cw experiment.If we define the functions which provides a basis for deriving both exact results and approximations to these results.In the final step numerical results can be found by numerically integrating the two-dimensional integral ͑3.41͒ These were evaluated by using an FFT routine with a grid of 64ϫ64 points.In the following subsection we present some of the results obtained by a numerically integrating the last equation, and finally comparing them to approximations derived from our exact results.

D. Approximations
All of the analysis to this point has been aimed at deriving exact results.All of these have been given in the jointtransform domain.It is possible, in principle, to derive results in the space-time domain by numerically inverting Eq. ͑3.41͒ with respect to both and s.However, enough is known about the order of magnitude of parameters in the context of biological media to introduce approximations allowing one to derive some qualitative implications from the analytical results.
To begin with we derive a usable approximation for the k-dependent term ͓F ˆs (1) (,Z͉0,Z)͔ kϪ1 that appears in the exact expression for the joint transform in Eq. ͑3.31͒.This will be done by restricting the values of Z and to be large and by considering the cw experiment, which allows us to set s ϭ0 in that equation thereby reducing the problem to that of inverting just a two-dimensional Fourier transform.
The first restriction is that Z be large enough so that the results are practically unaffected by the presence of a boundary.A lattice spacing for biological media is approximately equal to 1 mm.In practical applications involving biological media we can therefore safely assume that Zу5, restricting the penetrations to 5 mm or more since otherwise boundary effects can be quite significant.This restriction justifies our ignoring the second term in brackets in Eq. ͑3.37͒ since it is dominated by the factor 1 that appears in the brackets.
As a second approximations we restrict ourselves to the consideration of distances to the detector, which satisfies the condition 2 ϭ•ӷ1.This restriction is needed to satisfy the physically motivated requirement that details of lattice structure should not appear in the final result.This restriction in the space domain is equivalent to in the Fourier domain, allowing us to rewrite Eq. ͑3.38͒ as to lowest order in .The restriction to large Z allows us to approximate to the parameters ␣ and ␤ in Eq. ͑3.39͒ as follows: where, for many biological tissues of interest the dimensionless absorption parameter is quite small; a value of around 0.01 is a typical order of magnitude ͓26͔.Taking all of these approximations into account we can approximate to v ˆk(Z,,) in terms of a sum of three terms: In order to invert Eq. ͑3.45͒ it is necessary to invert the Fourier transform of terms of the form U ˆm.
Because our assumptions produce what is essentially a continuum limit, the Fourier series may be replaced by a Fourier transform.Further, since the transforms depend only on the magnitude and not the full vector the function v k (Z,,) will be a function of the magnitude .This allows us to write the inverse transform of U ˆm() as Equation ͑3.45͒ indicates that v k (Z,,) can be written in terms of these functions as from which our numerical results will be calculated.Notice that this approximation depends only on the magnitude rather than on the full vector , which indicates that we are in the largeregime.This has been shown to be a reasonably good approximation in ͓20͔.
In Fig. 3 we compare values of log 10 ͓v k (5͉)͔ to those obtained from Eq. ͑3.49͒ by plotting them as a function of k for different values of with ϭ0.02.The figure suggests that the accuracy of the approximation improves with increasing values of as has already been anticipated.Although there is a systematic deviation between the approximation and numerically calculated results at ϭ3 the error in using Eq.͑3.49͒ is still less than 15% at the lowest value of k and less than that at the intermediate values in the range shown.In Fig. 4 we compare numerical values of the expected value of the number of visits to depth Z as found numerically, against the formula ͗k͉͘ϭ given in ͓25͔.The agreement between results produced by this approximation and the numerical results is seen from the data shown in the figure to be quite good at lower values of at all of the values of Z.The accuracy of the approximation is seen to decrease as is increased.Again, the accuracy improves as the target depth is increased.Figure 5 shows plots of the standard deviation of the number of visits to level Z as a function of the distance .Further exploration of the data seems to suggest that the coefficient of variation, (k͉)/͗k͉͘ decreases as increases.This is to be expected on the consideration that as the radial distance increases, the photon trajectory tends to spend more time at points close to the surface thereby minimizing the possibility of internal absorption.

E. Some concluding remarks
In the present paper we have produced results on the statistical properties of the depth probed by a photon in a cw measurement.The results are given in the Fourier-Laplace domain, and potentially allows us to examine these statistical properties as a function of both space and time.We have only discussed results for the cw experiments, which can be obtained by simply setting the Laplace transform parameter equal to zero.However, the results can equally well be used  to investigate analogous properties for time-gated experiments by a further inversion, most likely numerical, of the Laplace transform.We hope to present this extension at some later date.It is also easy to analyze frequency-domain spectroscopy ͓27,28͔ using the present formalism by replacing the parameter s by i, where is the dimensionless frequency.
While the model that forms the basis of this work is somewhat specialized we would expect it to furnish reliable results except possibly at the very earliest times.The use of diffusion-based models or random walk models is suspect at extremely short times due to the neglect of anisotropic scattering effects.

FIG. 1 .
FIG. 1. Schematic diagram of the lattice random model of a semi-infinite medium with laser injection at a point.The parameter z measures the distance in integer units into the medium, and is the transverse distance along the surface.

FIG. 2 .
FIG.2.Curves of the probability that level Z is never visited conditional on reaching some point of the trapping surface at a dimensionless time .The curves are plotted as a function of the trapping time .

FIG. 3 .
FIG.3.Plot of numerically generated values of log 10 ͓v k (5͉)͔ ͑circles͒ as a function of k, compared to the approximation in Eq. ͑3.48͒ ͑lines͒.The function v k (5͉) is the probability that zϭ5 has been visited k times by a photon later absorbed at a distance from the point at which photons enter the tissue.In generating the curves the dimensionless Beer's law parameter was taken to be ϭ0.02.

FIG. 4 .
FIG. 4. Numerically generated values of the average number of visits to depth Z in a cw experiment plotted as a function of the distance from the input point at which the photon emerges at the trapping surface.The circles were obtained numerically while the dashed lines indicate approximations obtained in ͓12͔ and ͓25͔.The dimensionless Beer's law parameter is equal to 0.02.