Analytic function for predicting light fluence rate of circular fields on a semi-infinite turbid medium

Yi Hong Ong; Timothy C. Zhu

doi:10.1364/OE.24.026261

1. Introduction

Light propagation in tissue has been extensively studied for the past three decades to understand how light energy is distributed and absorbed within the target tissue. The in vivo light fluence distribution is influenced by numerous factors including the tissue optical properties, the tissue geometry and the irradiation geometry. Accurate determination of light fluence distribution in tissue is important to ensure successful photodynamic therapy (PDT) treatment as the rate of photochemical reactions depends on the light dose delivered to the tissue.

In vivo light fluence distribution can be directly measured by inserting an isotropic detector interstitially [1]. Due to the invasiveness, this technique is limited to measurement of light fluence rate at selected points or along a catheter rather than in a full three-dimensional volume. Model calculations such as Monte Carlo simulation [2] and diffusion approximation theory [3] are therefore more practical and have been widely employed to estimate the light distribution in tissue noninvasively, given that the 3D distribution of the tissue optical properties is known. Diffusion theory is attractive due to its simplicity in calculating light fluence rate quickly and accurately in scattering media. It has been found useful in various tissue optics applications as its validity condition of $μ_{a}$ << $μ_{s}^{'}$ is satisfied for most biological tissues at visible to near infrared wavelengths [4]. However, for regions that are close to sources or boundaries, diffusion theory is not appropriate as photons have a preferential direction of travel and cannot be approximated as a diffusion process. The light distribution described by diffusion theory is therefore less accurate in the first several hundred microns of tissue which can be an important region for PDT light-tissue interactions for surface irradiation [4]. Monte Carlo modeling, on the other hand, can be used to precisely simulate the light distribution in tissue with complex geometries and inhomogeneities and is often used as the gold standard method to model light transport in turbid media. However, the requirement of intensive computation load to achieve results with high accuracy has become the main limitation of MC method in application where speed is desirable [5]. Simple analytical or empirical expressions which allow for fast calculation of light distribution that mimic the accuracy of Monte Carlo simulations are therefore desirable in PDT treatment planning to optimize light fluence rate distribution during treatment.

Literature search reveals only analytical expression for broad beams, when the lateral dimension is much larger than the mean-free path [6, 7]. The light distribution of broad beams in tissue is only dependent on the two dimensional tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ) which can be simplified as a one-dimensional function of transport single scattering albedo, $a' = μ_{s}^{'} / (μ_{a} + μ_{s}^{'})$ , or diffuse reflectance, R_d, so long as all spatial dimensions are scaled by the mean-free path, $l = 1 / (μ_{a} + μ_{s}^{'})$ [6]. For small diameter circular light fields, the distribution of light fluence rate becomes dependent on the field size when the lateral dimensions are comparable to the effective mean-free-path of photons for optical properties relevant in PDT. Although this field size effect is not of great concern in most clinical PDT in which expanded or broad beams are usually used, it has to be taken into careful consideration in the treatment planning of experimental PDT using small animals where small diameter beams are often used to treat surface lesion as small as 1-cm in diameter [8, 9]. For clinical cases such as PDT of skin and recurrent breast cancers on the chest wall, small circular light fields of ≤ 2 cm diameter are sometimes used. It is common that the size of the treatment beam is determined based on the tumor size and the identical light fluence rate is used for the PDT treatment irrespective of the beam size [10].

MC simulation of light distribution clearly demonstrated that in vivo light fluence distribution changes with incident field size with significant reduction of light fluence rate at the central axis for small diameter circular beam from that for a broad beam [11–13]. Light scattering causes the light fluence distribution in tissue to be dependent on the incident beam diameter. Photons are scattered and diffuse radially away from the source. Closer to the edge of the beam, there are reduced scattering photons thus creating a large radial gradient in fluence rate at the beam perimeter. The edge effects overlap at the center of the beam when the incident beam diameter is narrow causing the central light fluence rate to be significantly reduced [12]. As the beam diameter increases, the central axis is less affected by the edge effects and the radius dependence of the in vivo light fluence distribution reduces. When the incident beam diameter is broader than 3 cm, the fluence rate at the center of beam approaches that for an infinitely wide beam with uniform incident irradiance [12]. It is therefore important for PDT treatments using small irradiance beam to be calibrated carefully for radius dependence in order to deliver optimal light dose which would improve the PDT treatment outcome.

In this paper, we will determine the central longitudinal light fluence rate distribution of circular fields with radii of 0.5 cm to 8 cm for a range of tissue optical properties which are relevant in PDT and many medical applications (absorption coefficients ( $μ_{a}$ ) between 0.01 and 1 cm⁻¹ and reduced scattering coefficients ( $μ_{s}^{'}$ ) between 2 and 40 cm⁻¹) based on a review on the in-vivo tissue optical properties [14]. These results are expressed as simple analytical functions in term of beam radius and tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ), and will be compared with Monte-Carlo simulations in term of maximum relative deviations as summarized in Tables 1 and 2. We used Monte-Carlo simulation because the diffusion theory is not valid when the lateral dimension of beam geometry becomes comparable to the mean-free-path of the photons or when the region of interest is near the air-tissue interface. These results will then be compared against the existing expressions in the literature for both accuracy and applicable range.

Table 1. Summary of maximum relative deviation (%) of the 6-parameter fitting using Eqs. (11)–(16) to the MC calculated $ϕ / ϕ_{a i r}$ for optical properties (μ_a, μ_s’) used in model training. Maximum relative deviation of 4 parameters fit for broad beam (radius = 8 cm) using Eqs. (6)–(9) are shown in bracket.

View Table | View all tables in this article

Table 2. Summary of maximum relative deviation (%) of the 6-parameter fitting using Eqs. (11)–(16) to the MC calculated $ϕ / ϕ_{a i r}$ for optical properties (μ_a, μ_s’) not used in model training. Maximum relative deviation of 4 parameters fit for broad beam (radius = 8 cm) using Eqs. (6)–(9) are shown in bracket.

View Table | View all tables in this article

2. Materials and Methods

2.1 MC Simulation

The setup geometry to be calculated by Monte-Carlo simulation for a semi-infinite medium with uniform optical properties, i.e. the absorption coefficient µ_a, the scattering coefficient µ_s, the scattering anisotropy g ( = 0.9), and the index of refraction (n₂ = 1.4 for tissue) is shown in Fig. 1(a). The outside medium is air (n₁ = 1). The light field is parallel and uniform inside a circle with radius R and incident toward the air-tissue interface. (R = 0.5, 0.75, 1, 1.5, 2, 3 and 8 cm.) Most tissues have a preferential scattering in the forward direction with a scattering anisotropy g = 0.90 [15]. In the diffusion approximation, the tissue scattering is converted to isotropic conditions (g = 0) and the tissue scattering property is described by a reduced scattering coefficient $μ_{s}^{'} = μ_{s} (1 - g)$ .

Fig. 1 (a) Setup geometry for Monte-Carlo simulation for a semi-infinite turbid medium. The light fluence rate is calculated along the central-axis on the center of the field. (b) Pencil beam setup is actually used for the Monte Carlo simulation. The ring scoring voxel has radius width dr and thickness dz at depth z and radius r in a cylindrical geometry. (c) Simplified setup geometry for Monte-Carlo simulation for a semi-infinite turbid medium in a water-tissue interface, where n₁ = 1.33 for water, n₂ = 1.4 for tissue and n₀ is always 1.

Download Full Size | PDF

To improve the scoring efficiency, Monte-Carlo simulation is performed for a pencil-beam incident on a semi-infinite tissue phantom as shown in Fig. 1(b). The photon density ( $ρ_{n}$ ) is scored in a cylindrical geometric relationship consisting of rings of thickness (dz) and width (dr). The light fluence rate (ϕ) per unit incident light power of the pencil beam (in unit of cm⁻²) is calculated as $ρ_{n} / N_{i n c} / μ_{a}$ , where $ρ_{n}$ is the photon density (in unit of cm⁻³), $N_{i n c}$ is the number of incident photons, and $μ_{a}$ is the absorption coefficient (in unit of cm⁻¹). Reciprocity theorem [16] is used to calculate ratio of the light fluence rate ϕ to the incident irradiance (or light fluence rate in air, $ϕ_{a i r}$ ) on the central axis of a circular field with radius R by an area integral [11]:

ϕ / ϕ_{a i r} = \int_{0}^{R} ρ_{n} / N_{i n c} / μ_{a} 2 π r d r,

We have also scored the fractional photons escaping the air-tissue interface into the air as a function of lateral radius (with a resolution of 0.05 cm) and emitting angle, θ (with a resolution of 1°). The reciprocity theorem is used to calculate the angular distribution of differential diffuse reflectance for circular fields as a function of radius and emitting angle. The light fluence rate above the tissue surface is calculated using the diffuse reflectance r_d for circular field with radius r, as previously described [11, 17, 18],

ϕ / ϕ_{a i r} = 1 + 2 r_{d},

Equation (2) states that the light fluence at a tissue surface is composed of 1 part from the incident photons and 2r_d part from the backscattered photons as described in the literature [19].

The Monte Carlo algorithm is similar to approaches used previously in the literature [2, 11, 20, 21]. The effect of reflection at the air-tissue interface, resulting from the refractive index mismatch, was modeled according to the Fresnel reflection coefficient for unpolarized light [20]. The initial weight of each launched photon is assigned to be unity. Photons are propagated through the medium step by step with random step size and scattering angle sampled based on the optical properties of the tissue model and their respective probability distributions. Each step size, $Δ s$ , is determined by $Δ s = l \cdot \ln ξ$ , where $ξ$ is a uniform distributed random number and l is the mean free path. The scattering angle is selected from the phase function determined by the Henyey-Greenstein phase function [21] with anisotropy g = 0.90. At the end of each step, the weight of the photon is reduced according to the absorption probability. At least 2 millions incident photons are launched normal to the interface along the z direction for each Monte Carlo simulation. We followed each of them until it escapes from the tissue model or when it is completely absorbed. The MC result provides an approximation to light distribution in tissue with a typical statistical uncertainty of less than 0.1%. The small standard variation for 2 millions photons is possible because a convolution of pencil-beam result is used to calculate the broad beam parameters.

2.2 Analytical expression and fitting algorithm

For small diameter circular fields, we found that Eq. (3) can fit the resulting MC data very well, typically to a depth of at least 0.5 cm or a few mean-free-path, l, whichever is larger. This fitting equation is termed 6-parameter model since it includes 6 fitting parameters: λ₁, λ₂, λ₃, b, C₂, C₃:

ϕ / ϕ_{a i r} = (1 - b \cdot e^{- λ_{1} d}) (C_{2} e^{- λ_{2} d} + C_{3} e^{- λ_{3} d}),

where

ϕ_{a i r} = P / A

,

A = π R^{2}

is called in-air fluence rate that can be calculated using the total power, P (in mW), and the treatment area, A, of the collimated beam. The 6-parameter model (

ϕ / ϕ_{a i r} = (1 - b \cdot e^{- λ_{1} d}) (C_{2} e^{- λ_{2} d} + C_{3} e^{- λ_{3} d})

) includes multiplication of two terms: the first term (

1 - b \cdot e^{- λ_{1} d}

) describes the buildup of scattered photons, the second term (

C_{2} e^{- λ_{2} d} + C_{3} e^{- λ_{3} d}

) describes the attenuation of the scattered photons. However, a single exponential decay is found to be insufficient to describe the attenuation of the scattered photons for small field sizes, thus we have added a second exponential decay term. λ₁ is the rate of buildup of scattered photons, 1-b is the ratio of fluence rate at surface to in-air fluence rate,

ϕ / ϕ_{a i r}

(d = 0), C₂ is the backscatter coefficients for scattered photons of broad beam and λ₂ is the effective attenuation coefficient of scattered photons,

C_{3} e^{- λ_{3} d}

term describes the reduction of the backscatter coefficients for scattered photons near the surface for small field sizes.

For a broad beam, existing literature shows that a 4-parameter model (λ₁, λ₂, b, C₂) is sufficient to fit the MC data [6, 7]:

ϕ / ϕ_{a i r} = (1 - b \cdot e^{- λ_{1} d}) \cdot C_{2} e^{- λ_{2} d},

We have modified the original form ( $ϕ / ϕ_{a i r} = C_{2}^{'} e^{- λ_{2} d} - C_{1}^{'} e^{- λ_{1} d}$ ) [6, 7] to be Eq. (4) because the term $(1 - b \cdot e^{- λ_{1} d})$ that fits the buildup region only can be taken out of the final results, and one will automatically get the analytical expression for diffusion approximation, i.e., C₂ is the backscattering coefficient and $λ_{2} = μ_{e f f}$ is the effective attenuation coefficient.

The objective function for minimization is defined as the least square difference between analytical fit of fluence rate, $ϕ_{f i t}$ , and the MC simulation, $ϕ_{M C}$ , for a range of tissue optical properties μ_a between 0.01 to 1 cm⁻¹ and μ_s’ between 2 and 40 cm⁻¹ and beam radius between 0.5 and 8 cm.

f = \sqrt{\frac{{\sum_{i = 1}^{N} (1 - \frac{ϕ_{f i t} (z_{i}, p a r a m s)}{ϕ_{M C} (z_{i})})}^{2}}{N (N - 1)}},

where

ϕ_{f i t}

is based on either Eq. (3) and/or Eq. (4) to result in either 6 or 4 parameters as the results of optimization.

ϕ_{M C}

is the MC simulation result as a function of discrete depth

z_{i}

, i = 1, 2, …, N. Typical value of

z_{i}

is between 0 and 2 cm in 0.01 cm steps, i.e., N = 200. Multi-variable optimization algorithm uses a differential evolution algorithm described elsewhere [22]. The non-linear optimization routine was performed using Matlab.

3. Results

3.1 Fitting results of MC simulation

Figure 2 shows the MC simulations of light fluence rate $ϕ / ϕ_{a i r}$ for air-tissue interface and their fits for four different optical properties ( $μ_{a}, μ_{s}^{'}$ ) = (a) (0.01, 2), (b) (0.1, 10), (c) (0.5, 10), and (d) (1, 40) cm⁻¹ and seven different beam radii, 0.5, 0.75, 1, 1.5, 2, 3, and 8 cm. For optical properties (a) and (b), the MC simulated light fluence rates inside tissue for small circular fields and broad beam were fitted using 6-parameters Eq. (3) and 4-parameters Eq. (4), respectively. For optical properties (c) and (d), only 4-parameters fit is necessary to fit the MC data for all radii. The light fluence rates on tissue surface were fitted using Eq. (17) as shown later in this paper.

Fig. 2 Light fluence rate $ϕ / ϕ_{a i r}$ for tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ) = (a) (0.01, 2), (b) (0.1, 10), (c) (0.5, 10), and (d) (1, 40) cm⁻¹ and ratio of the index of refraction, n₂/n₁ = 1.4. Solid lines are MC simulations. Dashed lines in (a) and (b) are 6-parameters fits using Eq. (3) for small circular fields and dotted lines are 4-parameters fits using Eq. (4) for broad beam. Dashed lines in (c) and (d) are for 4-parameters fits using Eq. (4). The curves in each plot are for seven different beam radii, from bottom to top, r = 0.5, 0.75, 1, 1.5, 2, and 3 cm for small fields, and 8 cm for broad beam.

Download Full Size | PDF

The fitting results are presented in Tables 3 to 7 in the appendix section. Tables 3, 4, and 5 show the summary results of 6-parameter fit to the MC simulation data for small circular fields with diameter of 1 cm, 1.5 cm, and 2 cm, respectively, for depth up to 3.5 cm deep using Eq. (3). For 1 cm and 1.5 cm diameter circular beams, C₃ approaches zero for a subset of data when μ_eff > 2 cm⁻¹ while for 2 cm diameter circular beam, C₃ approaches zero when μ_eff > 1 cm⁻¹ so that only a 4-parameter fit is necessary. In these cases, the value of λ₃ is unnecessary and its value is expressed as “—” in Tables 3, 4 and 5. Tables 6 and 7 show the summary results of 4-parameter fit using Eq. (4) to the broad beam MC simulation data (diameter of 8 cm) for an air-tissue interface (n₂/n₁ = 1.4) and water-tissue interface (n₂/n₁ = 1.053), respectively.

Table 3. Summary of the 6 fitting parameters (b, λ₁, λ₂, C₂, λ₃, C₃) using Eq. (3) to the MC calculated ϕ/ϕ_air for 1-cm diameter circular beam for optical properties (µ_a, µ_s’) and index of refraction, n₂/n₁ = 1.4. MC calculated diffuse reflectance for broad beam, R_d, is also listed.

View Table | View all tables in this article

Table 4. Summary of the 6 fitting parameters (b, λ₁, λ₂, C₂, λ₃, C₃) using Eq. (3) to the MC calculated ϕ/ϕ_air for 1.5-cm diameter circular beam for optical properties (µ_a, µ_s’) and index of refraction, n₂/n₁ = 1.4.

View Table | View all tables in this article

Table 5. Summary of the 6 fitting parameters (b, λ₁, λ₂, C₂, λ₃, C₃) using Eq. (3) to the MC calculated ϕ/ϕ_air for 2-cm diameter circular beam for optical properties (µ_a, µ_s’) and index of refraction, n₂/n₁ = 1.4.

View Table | View all tables in this article

Table 6. Summary of the 4 fitting parameters (b, λ₁, λ₂, C₂) using Eq. (4) to the MC calculated ϕ/ϕ_air for broad beam for optical properties (µ_a, µ_s’) and index of refraction, n₂/n₁ = 1.4.

View Table | View all tables in this article

Table 7. Summary of the 4 fitting parameters (b, λ₁, λ₂, C₂) using Eq. (4) to the MC calculated ϕ/ϕ_air for broad beam for optical properties (µ_a, µ_s’) and index of refraction, n₂/n₁ = 1.053.

View Table | View all tables in this article

3.2 Relationship between fitting parameters and the tissue optical properties for air-tissue interface

3.2.1 Broad beam

For light fluence rate of broad beams, we have found that all four fitting parameters in Eq. (4) are one-dimensional function of R_d as used by Jacques [6], where R_d is the diffuse reflectance of broad beam. The analytic expression of R_d can be found in the discussion section of this paper. The two fitting parameters (λ₁ and λ₂) should be proportional to μ_eff . We have found the following relationships:

Inside tissue:

b = 1.154 - 1.025 \cdot R_{d}^{0.1823},

\frac{λ_{1}}{μ_{e f f}} = 2.284 \cdot e^{3.19 R_{d}} - 1.81,

C_{2} = 5.996 R_{d} + 2.931,

\frac{λ_{2}}{μ_{e f f}} = 3.039 R_{d}^{5} - 8.023 R_{d}^{4} + 8.342 R_{d}^{3} - 4.306 R_{d}^{2} + 1.144 R_{d} + 0.865,

On tissue surface

ϕ / ϕ_{a i r} = 1 + 2 R_{d} .

3.2.2 Small circular field

For small circular field with beam diameter smaller or equal to 2 cm, we found that all the 6 parameters are as a function of radius, r. The parameters b, C₂ and λ₁ are two dimensional functions of R_d and radius, r and the parameters λ₂ is two dimensional function of μ_eff and radius, r. C₃ and λ₃ are three dimensional functions of $μ_{a}, μ_{s}^{'}$ and r. For air-tissue interface (n₂/n₁ = 1.4), we have found the following relationships for light fluence rate inside and on tissue surface:

Inside tissue:

b = ((e^{- 4.2 r} + 0.63) - (0.1 r^{- 1.35} + 0.28) e^{- 7 (1 - a')}) {(1 - a')}^{0.157 - 0.1 r^{- 0.35}},

\frac{λ_{1}}{μ_{e f f}} = (0.33 e^{- 0.96 r} + 2.16) e^{(3.26 - 3.99 e^{- 3.057 r}) R_{d}} - 1.8,

C_{2} = 6 R_{d} + 0.06 r^{- 1.8} + 2.945,

\frac{λ_{2}}{μ_{e f f}} = (1 - 0.6 e^{- 0.483 r}) μ_{e f f}^{0.236 e^{- 0.434 r}} + (4.37 e^{- 1.5 r}) μ_{e f f}^{- 1},

C_{3} = ((8.243 - 1.662 r) μ_{a} - 0.056 r^{- 3.395} - 1.251) \cdot ((0.802 r + 6.103) μ_{s}^{'}^{- (1.196 r^{8.315} + 0.212)} + 2.988 r^{4.152} - 3.131),

λ_{3} = ((30 - 0.005 r^{- 11.44}) μ_{a} + (3.39 - 2.404 r^{3.402})) \cdot {(0.973 - 0.0373 r^{2.345})}^{μ_{s}^{'}},

On tissue surface:

ϕ / ϕ_{a i r} = 1 + 2 R_{d} (1 - e^{- (0.56 (μ_{a}^{0.4} + 0.43) μ_{s}^{'}^{0.76} + 0.34) \cdot r}),

where

a' = \frac{μ_{s}^{'}}{μ_{s}^{'} + μ_{a}} .

4. Discussions

Figure 3 shows the MC simulated light fluence rate above tissue surface and their fits using Eq. (17) for five different optical properties, incident beam radii ranging from 0 to 8 cm and ratio of the index of refraction, n₂/n₁ = 1.4. The radius dependence of $ϕ / ϕ_{a i r}$ above tissue surface can be described as a function of tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ) and beam radius, r, as shown in Eq. (17). When the incident beam radius is very small, $1 - e^{- (0.56 (μ_{a}^{0.4} + 0.43) μ_{s}^{'}^{0.76} + 0.34) \cdot r}$ term approaches zero and the light fluence rate on tissue surface decreases to 1 (1 part from the incident photons). For broad beam where the r is large, $1 - e^{- (0.56 (μ_{a}^{0.4} + 0.43) μ_{s}^{'}^{0.76} + 0.34) \cdot r}$ approaches 1 and the light fluence rate above tissue surface becomes 1 + 2R_d (1 part from the incident photons and 2 part from the backscattered photons), as described by Marijnissen and Star [19].

Fig. 3 Radius dependence of light fluence rate $ϕ / ϕ_{a i r}$ above tissue surface (n₂/n₁ = 1.4)

Download Full Size | PDF

Figure 4 shows the fitting results for the 4-parameters Eq. (4) for broad beam, and Eqs. (6) – (9) and Eqs. (20)–(23) for their respective parameters at different refractive index mismatch. The fitting is performed using the curve fitting tool from Matlab global optimization toolbox (cftool.m). R² of the fitting results are all very good (R² > 0.96). For comparison, the fitting results for broad beam using the original form of Eq. (4) reported by Jacques [6] are included in Fig. 4. For broad beam, the diffuse reflectance can be used in place of the pairs of optical properties ( $μ_{a}, μ_{s}^{'}$ ) [6]. The diffuse reflectance values for broad beams R_d are obtained from MC simulation [11, 23] and are listed in Tables 3 to 7 and Fig. 5. Analytical expression based on the same form as a diffusion approximation [3, 24] can be determined between R_d and the albedo, a’ as follows:

R_{d} = 0.4843 a' \cdot (1 + e^{- 4.428 \sqrt{1 - a'}}) e^{- 2.65 \sqrt{1 - a'}},

where the albedo, a’ is defined in Eq. (18).

Fig. 4 Fitting results between the 4 parameters (b, C₂, λ₁/µ_eff, and λ₂/µ_eff) and the diffuse reflectance (R_d) for broad beams with air-tissue interface (n₂/n₁ = 1.4) and water-tissue interface (n₂/n₁ = 1.053). Solid lines are Eqs. (6)-(9) and dotted lines are Eqs. (20)-(23) for their respective parameters. Black circles and diamonds represent fitting results using 4 parameters from Table 6 and Table 7, respectively. Red triangles represent parameters using formulas from Ref [6]. R² refers to the goodness of fit between the fitting results of each parameter to their respective equation.

Download Full Size | PDF

Fig. 5 Relationship between R_d and a´/(1-a´). R² refers to the goodness of fit between the MC calculated R_d and Eqs. (19) and (24). Solid and dotted lines are Eqs. (19) and (24) respectively. Circles represent MC simulation results for air-tissue interface (n₂/n₁ = 1.4) while diamonds represent MC simulation results for water-tissue interface (n₂/n₁ = 1.053).

Download Full Size | PDF

Figure 5 shows the goodness of the fit between R_d and $a' / (1 - a')$ using Eqs. (19) and (24), for air-tissue and water-tissue interface, respectively.

Figure 6 shows the fitting results for the broad beam using 4-parameters, Eqs. (6)–(9) and small circular beams with diameter of 1 cm, 1.5 cm and 2 cm using 6-parameters, Eqs. (11)–(16). Notice that for small circular field, the scaling theorem is no longer in existence that most of the parameters have to be expressed as at least two dimensional function of diffuse reflectance, R_d or tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ) and beam radius, r. R² of the fitting results are also very good (R² > 0.95, not shown in Fig. 6). To minimize the potential of overfitting in the multi-variable optimization, we started the fitting using 4 parameters (Eq. (4)) to as deep a depth as possible that still fit the data, e.g., 0.5 cm for μ_a = 0.1 cm⁻¹ and μ_s’ = 10 cm⁻¹, of the MC data to get the best fit for parameters b and λ₁, which account for the buildup region of the light fluence rate, and C₂ and λ₂. Then, by making parameters b and λ₁ constants with values obtained from the first fit and C₂ value to be the same as that for broad beam, parameters λ₂, C₃ and λ₃ are refitted using Eq. (3) to MC data to a greater depth, at least 3.5 cm. By using this method, we are able to obtain consistent fitting results. A differential evolution algorithm [22] is used to perform the optimization and ideally find a global minimum for either 6 parameters or 4 parameters fitting results to the MC data. The differential evolution algorithm is a method that optimizes a problem by iteratively trying to improve a candidate solution with regard to a given measure of quality. This method is also known as metaheuristic as it makes few or no assumptions about the problem being optimized and can search very large spaces of candidate solutions. We found that 6 parameters are necessary to fit the data for small circular field data, particularly those with beam diameters smaller than 3 times light penetration depth, δ. Light penetration depth is the reciprocal of effective attenuation coefficient, µ_eff, where $μ_{e f f} = \sqrt{3 μ_{a} (μ_{a} + μ_{s}^{'})}$ . C₃ is found to be smaller than zero and the value decreases for smaller beam diameter and tissue optical properties ( $μ_{a}, μ_{s}^{'}$ ). C₃ term is necessary to improve the fitting result of small circular field to the MC results as it accounts for field size effect which reduces the light fluence rate at the central cylinder of the beam. C₃ increases with both beam diameter and tissue optical properties due to the increase in photon scattering by tissue and approaches zero when beam diameter equal to or larger than 3/µ_eff. Similarly, λ₃ can be described as a function of beam diameter and tissue optical properties. λ₃ is larger when beam diameter and tissue scattering decrease. For beam with diameter larger than 3/µ_eff, only 4 parameters are needed to fit the data to MC results as C₃ term approaches zero and λ₃ became irrelevant.

Fig. 6 Fitting results between the 4 and 6 parameters (b, C₂, C₃, λ₁/µ_eff, λ₂/µ_eff, and λ₃) and the optical properties (µ_a, µ_s’) for circular beams with radii = 0.5 cm, 0.75 cm, 1 cm and 8 cm and refractive index mismatch, n₂/n₁ = 1.4. Black solid lines are Eqs. (6)-(9) while the red, blue and green lines are Eqs. (11)-(16) for their respective parameters. Circles with their respective color represent MC fitting results using 4 or 6 parameters from Tables 3 to 6.

Download Full Size | PDF

Figure 7 shows the dependence of parameters b, C₂, λ₁ and λ₂ to beam radius using Eqs. (11)-(14). b accounts for escaping reflectance which causes a depletion of fluence rate near the surface. From Fig. 7 we can see that b increases when beam radius is smaller than ~1.5 cm for albedo smaller or equal to 0.8 and when beam radius is smaller than 1 cm for albedo > 0.8. b shows higher radius dependence for smaller albedo. b increases by 23% when beam radius decreases from 1.5 cm to 0.5 cm for albedo = 0.6, and the increase in b is about 10% for the similar decrease in beam radius when albedo is 0.95. λ₁/µ_eff increases rapidly with beam radius from 0.5 cm to 2 cm when R_d is high. The radius dependence of λ₁/µ_eff is smaller when R_d is small and λ₁/µ_eff does not change when beam radius is larger than 2 cm. Combining both parameters in the buildup term $(1 - b \cdot e^{- λ_{1} d})$ in Eqs. (3) and (4), we can see a higher drop in the fluence rate near the surface when beam radius is small and a faster build up rate when beam radius and R_d are large. C₂ shows low dependence on beam radius and Eq. (13) approaches Eq. (8) for all beam radii, where C₂ can be described as a one-dimensional function of R_d. λ₂/µ_eff increases with decreasing beam radius and approaches to 1 when beam radius is larger than 2.5 cm. For narrow field with beam radius ≤ 2.5 cm, radius dependence of λ₂/µ_eff increases with decreasing µ_eff. For beam radius = 0.5 cm, λ₂ is 4.5 times of µ_eff when µ_eff is 0.5 cm⁻¹, 2.5 times when µ_eff is 1 cm⁻¹, 1.5 times when µ_eff is 2.5 cm⁻¹ and the ratio approaches to 1 when µ_eff is larger than 7.5 cm⁻¹.

Fig. 7 Relationship between parameters b, C₂, λ₁/µ_eff, and λ₂/µ_eff using Eqs. (11)–(14) and beam radius.

Download Full Size | PDF

Figure 8(a) shows the comparison between MC simulation results of $ϕ / ϕ_{a i r}$ broad beam using the fitting parameters in Table 6, as well as Eqs. (6)–(9). For comparison, we have also included the fitting results in the literature [6]. We can see that the fitting results from Ref [6]. are good for μ_a ≥ 0.1 cm⁻¹ and it is unable to fit the buildup region accurately for μ_a = 0.01 cm⁻¹. The fitting results using 4 parameters and Eqs. (6)–(9) are very good for all three optical properties presented in Fig. 8(a) for depth up to 3.5 cm. Equation (4), which is modified from the original form from Ref [6]. is able to improve the fitting crosstalks between parameters b and λ₁, which correspond to the buildup region, as shown in Fig. 4. Figure 8(b) shows the comparison between MC simulation results of $ϕ / ϕ_{a i r}$ of small diameter circular beam with radius of 0.5 cm, 0.75 cm and 1 cm using the fitting parameters in Tables 3, 4, and 5, respectively, as well as Eqs. (11)–(16). We can see greater radius dependence of $ϕ / ϕ_{a i r}$ in small circular diameter beams when μ_a and μ_s’ decreases. In fact, Zhu et. al. [11] has reported that the radius dependence of $ϕ / ϕ_{a i r}$ is a function of μ_a and μ_s’. The relative deviations of both broad beam and small diameter circular beams for these three tissue optical properties over a depth of 3.5 cm are presented in Fig. 9. Relative deviation is defined as the difference between the fit and the MC result relative to the maximum fluence rate.

Fig. 8 $ϕ / ϕ_{a i r}$ for (a) broad beam and (b) small diameter circular beams for three tissue optical properties (µ_a, µ_s’) = (i) (1, 40), (ii) (0.1, 10), and (iii) (0.01, 2) cm⁻¹. Solid lines are MC simulations. Dotted lines are fittings using Eqs. (6)–(9) for (a) and Eqs. (11)–(16) for (b). Dashed line are 4 parameters fits for (a) and 6 parameters fits for (b). Black lines in (b) correspond to beam with radius = 0.5 cm, blue lines correspond to beam with radius = 0.75 cm, and red lines correspond to beam with radius = 1 cm.

Download Full Size | PDF

Fig. 9 Relative deviation of the fitting results to the MC simulation data for (a) broad beam and (b) small diameter circular beam for three tissue optical properties (μ_a, μ_s’) = (i) (1, 40), (ii) (0.1, 10), and (iii) (0.01, 2) cm⁻¹. Data for (i) and (ii) are vertically shifted for clarity: + 25% for (ii) and + 30% for (i) for (a) and + 10% for (ii) and + 16% for (i) for (b), respectively. Black lines in (b) correspond to beam with radius = 0.5 cm, blue lines correspond to beam with radius = 0.75 cm, and red lines correspond to beam with radius = 1 cm.

Download Full Size | PDF

From Fig. 9(a), we can see that Eqs. (6)–(9) fits all broad beam MC data to a relative deviation of ± 3% for three optical properties up to 3.5 cm depths, while parameters from Table 6 fits the MC data to within 1% relative deviation. Parameters from Ref [6]. fit most of the data except for = 0.01 cm⁻¹, where 20% error is observed. Among the difference between our MC data and Jacque’s fit, 3% can be attributed to difference in the value of index of refractions used for the MC simulation: 1.4 and 1.38 between us and Jacques, as well as other subtle differences in MC simulation. Figure 9(b) shows the relative deviation of small diameter circular beam fits to the MC simulation results. For optical properties = 1 cm⁻¹ and = 40 cm⁻¹, both 6 parameters fitting and Eqs. (11)–(16) are able to fit the MC calculated $ϕ / ϕ_{a i r}$ accurately to within 1.5% relative deviation. It should be noted that C₃ and λ₃ term can be ignored for this optical properties as they approach to zero for all the beam radii studied. A 4-parameters fitting using Eq. (4) should give similar accuracy to the results presented. For optical properties μ_a = 0.1 cm⁻¹ and μ_s’ = 10 cm⁻¹, the maximum relative deviation of 6 parameters fitting and Eqs. (11)–(16) to the MC calculated $ϕ / ϕ_{a i r}$ is 4.5%. For optical properties μ_a = 0.01 cm⁻¹ and μ_s’ = 2 cm⁻¹, 6 parameters fitting and Eqs. (11)–(16) are able to fit the MC calculated $ϕ / ϕ_{a i r}$ to within 4% and 9% relative deviation. The maximum relative deviations of the fittings using Eqs. (6)–(9) and Eqs. (11)–(16) to broad beam and small diameter circular beams for all optical properties studied, up to a depth of 3.5 cm, are summarized in Tables 1 and 2.

Table 1 shows the maximum relative deviation of the 4 and 6 parameters fitting to the MC simulation results for optical properties used in model training (µ_a = 0.01, 0.1, 0.5, 1 cm⁻¹ and µ_s’ = 2, 10, 23.5, 34.5, 40 cm⁻¹) and beam radii = 0.5, 0.75, 1, 3, 8 cm. Independent MC simulations using different combinations of optical properties (µ_a = 0.05, 0.3 and 0.8 cm⁻¹; µ_s’ = 5, 15, 25 cm⁻¹) and same beam radii were run to compare the performance of the 6-parameter model. The maximum relative deviations of the fitting at these additional optical properties to the MC simulation data are shown in Table 2. It should be noted that 4 parameters fitting using Eqs. (11)–(14) are performed for small circular fields with beam diameter larger than 3 times light penetration depth, δ. The maximum relative deviation is defined as the maximum difference between the fit and the MC result relative to the maximum fluence rate. 4-parameters fitting using Eqs. (6)-(9) is sufficient to fit for broad beam for all optical properties to within 4.5% relative deviation of the MC results. For small diameter circular beam (radius = 0.5 cm, 0.75 cm, 1 cm and 3 cm), 6-parameters fitting using Eqs. (11)–(16) is able to fit for higher scattering coefficient (μ_s’ ≥ 23.5 cm⁻¹) to within 6% relative deviation and for lower scattering coefficient (μ_s’ ≤ 10 cm⁻¹) to within 10.6% relative deviation of MC results. For small diameter beam, when beam diameter is greater than 3 times of light penetration depth, C₃ approaches zero and both C₃ and λ₃ terms has minimal impact on the overall fitting to the MC calculated $ϕ / ϕ_{a i r}$ . When beam diameter is smaller than 3 times of light penetration depth (or 3/µ_eff), C₃ and λ₃ terms are significant and has greater impact on the overall fitting to the MC calculated $ϕ / ϕ_{a i r}$ as beam diameter decreases. For instance, the light penetration depth at optical properties μ_a = 0.01 cm⁻¹ and μ_s’ = 2 cm⁻¹ is 4.07 cm. The diameter of a 1-cm beam at this optical properties is smaller than 3/µ_eff and the maximum relative deviation for 6-parameters fit to MC results is 3.5%. The maximum relative deviation increases to 140% (result not shown in Table 1) when the C₃ and λ₃ terms are excluded from the fitting.

Another advantage of our expressions, Eqs. (3) and (4), is that the fitting based on diffusion approximation can be obtained by dividing them by $(1 - b \cdot e^{- λ_{1} d})$ . The fitting results without the buildup term are shown in Fig. 10 for (a) broad beam and (b) small diameter circular beam, respectively. For optical properties (i) and (ii), Eqs. (6)–(9) and Eqs. (11)–(16) without the buildup region fit the broad beam and small diameter circular beam data very well. For low absorption coefficient (μ_a = 0.01 cm⁻¹), the buildup region term has a larger impact on the light fluence rate distribution on the central axis for both the broad beam and small diameter circular beam, extending up to 0.5 cm depth.

Fig. 10 $ϕ / ϕ_{a i r}$ for (a) broad beam with MC simulation and 4-parameter fitting using Eqs. (6)–(9) without the buildup term and (b) small diameter circular beam with MC simulation and 6-parameters fitting using Eqs. (11)–(16) without the buildup region for three tissue optical properties (μ_a, μ_s’) = (i) (1, 40), (ii) (0.1, 10), and (iii) (0.01, 2) cm⁻¹. Black lines in (b) correspond to beam with radius = 0.5 cm, blue lines correspond to beam with radius = 0.75 cm, and red lines correspond to beam with radius = 1 cm.

Download Full Size | PDF

Figure 11 shows the relative deviation between the fit without the buildup term and the MC simulation, normalized to the maximum MC calculated fluence rate. Beyond the buildup region, the maximum error is 1% for the broad beam and around 5% for the small diameter circular beam. Within the buildup region, the maximum error is 18% for the broad beam and 25% for the small diameter circular beam. The buildup term has minimal impact on the light fluence distribution on the central axis for both the broad beam and small diameter circular beam when the absorption and scattering coefficients are high (μ_a ≥ 0.1 cm⁻¹, μ_s’ ≥ 10 cm⁻¹) but has a larger impact on the central light fluence rate distribution for low absorption and scattering coefficients (μ_a = 0.01 cm⁻¹, μ_s’ = 2 cm⁻¹).

Fig. 11 (a) Relative deviation of the 4-parameters fitting using Eqs. (6)–(9) without the buildup term to the MC simulation data for broad beam. Relative standard deviation of 6-parameters fitting using Eqs. (11)–(16) without the buildup term to the MC simulation data for (b) 1-cm beam, (c) 1.5-cm beam and (d) 2-cm beam.

Download Full Size | PDF

4.1 Effect of refraction index mismatch

Intracavitory clinical PDT cases such as intraperitoneal and pleural PDT often involve water-tissue interface with a different refractive index mismatch than air-tissue interface [17, 25–27]. At water tissue interface the refractive index mismatch is about 1.053 with a critical angle of 71.8° compared to air-tissue interface with refractive index mismatch of 1.4 and critical angle of 48.8°. Fewer diffusely backscattered light that reaches the surface of the water-tissue boundary is reflected back into the medium. The reduction in total internal reflection and increment in escaping reflectance augments the behavior of light fluence rate near the surface. The expressions for the 4 parameters in Eqs. (6)-(9) with refractive index mismatch of 1.4 has been adjusted to mimic accurate Monte Carlo simulations of water-tissue interface as shown in Eqs. (20)–(23).

b = 0.943 - 0.531 R_{d}^{0.1946},

\frac{λ_{1}}{μ_{e f f}} = 1.572 e^{3.19 R_{d}} + 1.135 \times 10^{- 9} e^{26 R_{d}} - 1.8,

C_{2} = 2.388 R_{d} + 2.992,

\frac{λ_{2}}{μ_{e f f}} = 2.626 R_{d}^{5} - 7.451 R_{d}^{4} + 8.436 R_{d}^{3} - 4.81 R_{d}^{2} + 1.426 R_{d} + 0.812,

The light fluence rate on the tissue surface can be calculated using Eq. (10) where it is composed of 1 part from the incident photons and 2 R_d part from the backscattered photons. The value of R_d depends on both the optical properties and also the refractive mismatch at the tissue surface boundary. At water-tissue interface, the analytic expression of R_d can be determined as follows.

R_{d} = 0.5013 a' \cdot (1 + e^{- 2.789 \sqrt{1 - a'}}) e^{- 1.732 \sqrt{1 - a'}},

Good agreement is observed between the empirical solution of R_d and the MC simulation. The goodness of the fit between R_d using Eq. (24) and MC data as a function of a’/(1- a’) can be found in Fig. 5. It should be noted that the expression of diffuse reflectance R_d for Eqs. (19) and (24) depend only upon the transport albedo a’ and the ratio of the index of refraction, n₂/n₁, and one can calculated the total diffuse reflectance for any given optical properties and refractive indices using the analytical expression based on the same form as a diffusion approximation given by Farrell [3]. The radius dependence of R_d can be calculated using Eq. (17).

5. Conclusion

We have performed MC simulation for air-tissue interface and provided a 4-parameters and 6-parameters expressions for the longitudinal light fluence rate distribution on the central axis with an accuracy of better than a relative standard deviation of 4.5% for broad beams and 10.6% for the small diameter circular beams for the full range of the in-vivo tissue optical properties studied. Simple analytical expressions are provided to determine the parameters as a one dimensional function of diffuse reflectance, R_d, for broad beam and as at least two dimensional functions of tissue optical properties, μ_a and μ_s’, or diffuse reflectance and beam radius, r. We have also independently validated the analytical fitting result for broad beam and small diameter circular beams. We found that 4-parameters expression is able to fit accurately for circular beam with beam diameters larger than 3 times optical penetration depth for all optical properties studied. 6-parameters expression is necessary to fit for small circular diameter beams with beam diameters smaller than 3 times optical penetration depth. The additional two terms in 6-parameters expression, C₃ and λ₃, account for the field size effect for small diameter circular beams which reduces the light fluence rate at the central axis of the beam near the surface. Both 4 and 6 parameters expressions are simple to use and can accurately mimic the Monte Carlo simulations of light transport including at small depths close to air-tissue interface where diffusion theory fails. The accuracy and simplicity of these radial dependent analytic expressions would be extremely helpful to PDT in optimizing treatment plan involving various irradiance beam sizes that would improve the treatment outcome.

6 Appendix: Fitting parameters

Funding

National Institute of Health (NIH) (R01 CA154562, P01 CA87971).

References and links

1. A. Dimofte, J. C. Finlay, and T. C. Zhu, “A method for determination of the absorption and scattering properties interstitially in turbid media,” Phys. Med. Biol. 50(10), 2291–2311 (2005). [CrossRef] [PubMed]

2. L. Wang, S. L. Jacques, and L. Zheng, “MCML - Monte Carlo modeling of light transport in multi-layered tissues,” Comput. Methods Programs Biomed. 47(2), 131–146 (1995). [CrossRef] [PubMed]

3. T. J. Farrell, M. S. Patterson, and B. Wilson, “A diffusion theory model of spatially resolved, steady-state diffuse reflectance for the noninvasive determination of tissue optical properties in vivo,” Med. Phys. 19(4), 879–888 (1992). [CrossRef] [PubMed]

4. S. L. Jacques and B. W. Pogue, “Tutorial on diffuse light transport,” J. Biomed. Opt. 13(4), 041302 (2008). [CrossRef] [PubMed]

5. C. Zhu and Q. Liu, “Review of Monte Carlo modeling of light transport in tissues,” J. Biomed. Opt. 18(5), 050902 (2013). [CrossRef] [PubMed]

6. S. L. Jacques, “Light Distributions from Point, Line and Plane Sources for Photochemical Reactions and Fluorescence in Turbid Biological Tissues,” Photochem. Photobiol. 67(1), 23–32 (1998). [CrossRef] [PubMed]

7. C. M. Gardner, S. L. Jacques, and A. J. Welch, “Light transport in tissue: Accurate expressions for one-dimensional fluence rate and escape function based upon Monte Carlo simulation,” Lasers Surg. Med. 18(2), 129–138 (1996). [CrossRef] [PubMed]

8. H. Qiu, M. M. Kim, R. Penjweini, and T. C. Zhu, “Macroscopic singlet oxygen modeling for dosimetry of Photofrin-mediated photodynamic therapy: an in-vivo study,” J. Biomed. Opt. 21(8), 088002 (2016). [CrossRef] [PubMed]

9. R. Penjweini, M. M. Kim, B. Liu, and T. C. Zhu, “Evaluation of the 2-(1-Hexyloxyethyl)-2-devinyl pyropheophorbide (HPPH) mediated photodynamic therapy by macroscopic singlet oxygen modeling,” J. Biophotonics, published online Sept 22 (2016).

10. Q. Chen, B. C. Wilson, M. O. Dereski, M. S. Patterson, M. Chopp, and F. W. Hetzel, “Effects of light beam size on fluence distribution and depth of necrosis in superficially applied photodynamic therapy of normal rat brain,” Photochem. Photobiol. 56(3), 379–384 (1992). [CrossRef] [PubMed]

11. T. C. Zhu, A. Dimofte, S. M. Hahn, and R. A. Lustig, “Light dosimetry at tissue surfaces for small circular fields,” Proc SPIE Int Soc Opt Eng 4952, 56–67 (2003). [CrossRef] [PubMed]

12. M. Keijzer, S. L. Jacques, S. A. Prahl, and A. J. Welch, “Light distributions in artery tissue: Monte Carlo simulations for finite-diameter laser beams,” Lasers Surg. Med. 9(2), 148–154 (1989). [CrossRef] [PubMed]

13. W. M. Star, “Light dosimetry in vivo,” Phys. Med. Biol. 42(5), 763–787 (1997). [CrossRef] [PubMed]

14. J. L. Sandell and T. C. Zhu, “A review of in-vivo optical properties of human tissues and its impact on PDT,” J. Biophotonics 4(11-12), 773–787 (2011). [CrossRef] [PubMed]

15. W. F. Cheong, S. A. Prahl, and A. J. Welch, “A review of the optical properties of biological tissues,” IEEE J. Quantum Electron. 26(12), 2166–2185 (1990). [CrossRef]

16. F. H. Attix, Introduction to Radiological Physics and Radiation Dosimetry (1st. John Wiley & Sons, 1983), Chap. 3.

17. T. G. Vulcan, T. C. Zhu, C. E. Rodriguez, A. Hsi, D. L. Fraker, P. Baas, L. H. Murrer, W. M. Star, E. Glatstein, A. G. Yodh, and S. M. Hahn, “Comparison between isotropic and nonisotropic dosimetry systems during intraperitoneal photodynamic therapy,” Lasers Surg. Med. 26(3), 292–301 (2000). [CrossRef] [PubMed]

18. A. Dimofte, T. C. Zhu, S. M. Hahn, and R. A. Lustig, “In vivo light dosimetry for motexafin lutetium-mediated PDT of recurrent breast cancer,” Lasers Surg. Med. 31(5), 305–312 (2002). [CrossRef] [PubMed]

19. J. P. Marijnissen and W. M. Star, “Calibration of isotropic light dosimetry probes based on scattering bulbs in clear media,” Phys. Med. Biol. 41(7), 1191–1208 (1996). [CrossRef] [PubMed]

20. B. C. Wilson and G. Adam, “A Monte Carlo model for the absorption and flux distributions of light in tissue,” Med. Phys. 10(6), 824–830 (1983). [CrossRef] [PubMed]

21. S. A. Prahl, M. Keijzer, S. L. Jacques, and A. J. Welch, “A Monte Carlo model of light propagation in tissue,” Proc. SPIE 5, 102–111 (1989).

22. T. C. Zhu, B. E. Bjärngard, Y. Xiao, and M. Bieda, “Output ratio in air for MLC shaped irregular fields,” Med. Phys. 31(9), 2480–2490 (2004). [CrossRef] [PubMed]

23. T. C. Zhu, J. C. Finlay, A. Dimofte, and S. M. Hahn, “Light Dosimetry at Tissue Surfaces for Oblique Incident Circular Fields,” Proc SPIE Int Soc Opt Eng 5315, 113–124 (2004). [CrossRef] [PubMed]

24. A. Dimofte, T. C. Zhu, J. Finlay, M. Culligan, C. E. Edmonds, J. S. Friedberg, K. Cengel, and S. M. Hahn, “In-vivo Light Dosimetry for Pleural PDT,” Proc. SPIE 7164, 71640A (2009). [CrossRef]

25. C. B. Simone 2nd and K. A. Cengel, “Photodynamic Therapy for Lung Cancer and Malignant Pleural Mesothelioma,” Semin. Oncol. 41(6), 820–830 (2014). [CrossRef] [PubMed]

26. T. C. Zhu, M. M. Kim, X. Liang, B. Liu, J. L. Meo, J. C. Finlay, A. Dimofte, C. Rodriguez, C. B. Simone 2nd, K. Cengel, and J. S. Friedberg, “Real-time treatment feedback guidance of Pleural PDT,” Proc SPIE 8568, 85680O (2013). [CrossRef] [PubMed]

27. S. M. Hahn, D. L. Fraker, R. Mick, J. Metz, T. M. Busch, D. Smith, T. Zhu, C. Rodriguez, A. Dimofte, F. Spitz, M. Putt, S. C. Rubin, C. Menon, H. W. Wang, D. Shin, A. Yodh, and E. Glatstein, “A phase II trial of intraperitoneal photodynamic therapy for patients with peritoneal carcinomatosis and sarcomatosis,” Clin. Cancer Res. 12(8), 2517–2525 (2006). [CrossRef] [PubMed]

µ_a	µ_s’	R_d	µ_eff	b	λ₁	λ₂	λ₃	C₂	C₃
0.01	2	0.6907	0.246	0.22	2.5	2.18	3.03	7.3	−4.6
0.01	10	0.8317	0.548	0.19	8	2.32	2.12	8.1	−1.91
0.01	23.5	0.8777	0.840	0.17	14.5	2.55	1.6	8.35	−0.5
0.01	34.5	0.8936	1.017	0.17	19	2.66	1.06	8.45	−0.215
0.01	40	0.8991	1.096	0.17	21	2.71	0.8	8.5	0.015
0.1	2	0.3589	0.794	0.33	3	2.5	4.6	5.35	−2.95
0.1	10	0.6042	1.741	0.23	14.5	3	3.15	6.8	−1
0.1	23.5	0.7084	2.661	0.21	28	3.86	2	7.4	−0.135
0.1	34.5	0.7475	3.222	0.1853	38.6	4.27	—	7.578	—
0.1	40	0.7611	3.468	0.1763	43.7	4.43	—	7.665	—
0.5	2	0.1334	1.936	0.63	3	3	—	4	—
0.5	10	0.3589	3.969	0.3579	17.1	4.88	—	5.207	—
0.5	23.5	0.4947	6.000	0.2618	34	6.57	—	6.175	—
0.5	34.5	0.5518	7.246	0.2354	48.2	7.66	—	6.525	—
0.5	40	0.573	7.794	0.2252	55.2	8.15	—	6.64	—
1	2	0.071	3.000	0.66	3.1	3.7	—	3.8	—
1	10	0.252	5.745	0.4081	16.5	6.33	—	4.699	—
1	23.5	0.3848	8.573	0.3068	40	8.84	—	5.472	—
1	34.5	0.4463	10.320	0.2812	58.1	10.5	—	5.845	—
1	40	0.4694	11.091	0.2513	69.8	11.2	—	5.887	—

µ_a
0.01	2	0.6907	0.246	0.22	3.5	1.67	2.6	7.2	−4.1
0.01	10	0.8317	0.548	0.18	12	1.78	1.69	8.05	−1.04
0.01	23.5	0.8777	0.840	0.17	21.5	1.96	1	8.35	−0.271
0.01	34.5	0.8936	1.017	0.17	28.5	2	0.5	8.417	0.027
0.01	40	0.8991	1.096	0.1658	31.4	2.06	0.2	8.461	0.078
0.1	2	0.3589	0.794	0.32	4	2	5.09	5.2	−2.34
0.1	10	0.6042	1.741	0.22	19.5	2.46	3.1	6.75	−0.698
0.1	23.5	0.7084	2.661	0.2076	35.3	3.19	1.96	7.322	−0.119
0.1	34.5	0.7475	3.222	0.1763	55.7	3.61	—	7.484	—
0.1	40	0.7611	3.468	0.1725	63.8	3.81	—	7.564	—
0.5	2	0.1334	1.936	0.57	3.2	2.54	—	3.85	—
0.5	10	0.3589	3.969	0.3217	18.3	4.3	—	5.224	—
0.5	23.5	0.4947	6.000	0.2491	44.8	6.13	—	6.015	—
0.5	34.5	0.5518	7.246	0.2395	70.7	7.31	—	6.345	—
0.5	40	0.573	7.794	0.2381	83.8	7.84	—	6.468	—
1	2	0.071	3.000	0.59	3.2	3.25	—	3.5	—
1	10	0.252	5.745	0.3739	18.7	5.85	—	4.56	—
1	23.5	0.3848	8.573	0.3002	49.9	8.55	—	5.325	—
1	34.5	0.4463	10.320	0.2871	73.8	10.3	—	5.723	—
1	40	0.4694	11.091	0.2517	87.1	11.1	—	5.781	—

µ_a
0.01	2	0.6907	0.246	0.22	4	1.07	1.13	7.15	−3.05
0.01	10	0.8317	0.548	0.18	15	1.45	0.8	8	−0.2
0.01	23.5	0.8777	0.840	0.18	26.5	1.61	0.387	8.407	0.0753
0.01	34.5	0.8936	1.017	0.17	33.7	1.68	—	8.362	—
0.01	40	0.8991	1.096	0.1659	37	1.71	—	8.442	—
0.1	2	0.3589	0.794	0.33	4	1.6	3.41	5.15	−1.6
0.1	10	0.6042	1.741	0.22	21.5	2.24	2	6.638	0.0009
0.1	23.5	0.7084	2.661	0.2028	43.3	2.89	—	7.259	—
0.1	34.5	0.7475	3.222	0.1747	68	3.38	—	7.421	—
0.1	40	0.7611	3.468	0.174	78.8	3.6	—	7.499	—
0.5	2	0.1334	1.936	0.51	3.2	2.3	—	3.78	—
0.5	10	0.3589	3.969	0.3047	19.7	4.08	—	5.171	—
0.5	23.5	0.4947	6.000	0.2466	49.7	6.02	—	5.958	—
0.5	34.5	0.5518	7.246	0.241	77.1	7.24	—	6.304	—
0.5	40	0.573	7.794	0.241	90.9	7.79	—	6.433	—
1	2	0.071	3.000	0.57	3.3	3	—	3.48	—
1	10	0.252	5.745	0.3628	19.9	5.71	—	4.494	—
1	23.5	0.3848	8.573	0.299	52.2	8.51	—	5.297	—
1	34.5	0.4463	10.320	0.288	76	10.3	—	5.708	—
1	40	0.4694	11.091	0.2517	89.1	11	—	5.77	—

µ_a	µ_s’	R_d	b	λ₁	λ₂	C₂
0.01	2	0.6907	0.23	4.82	0.25	7.10
0.01	10	0.8317	0.16	16	0.55	7.89
0.01	23.5	0.8777	0.16	30	0.85	8.2
0.01	34.5	0.8936	0.15	38	1.03	8.24
0.01	40	0.8991	0.14	42	1.12	8.27
0.1	2	0.3589	0.29	4.21	0.79	5.1
0.1	10	0.6042	0.2	22.3	1.74	6.54
0.1	23.5	0.7084	0.2	52.6	2.67	7.19
0.1	34.5	0.7475	0.18	80.6	3.23	7.36
0.1	40	0.7611	0.18	89.7	3.48	7.46
0.5	2	0.1334	0.45	3.81	1.86	3.69
0.5	10	0.3589	0.29	22.1	3.93	5.08
0.5	23.5	0.4947	0.25	52.9	5.98	5.93
0.5	34.5	0.5518	0.24	77.4	7.23	6.30
0.5	40	0.5730	0.24	90.62	7.78	6.43
1	2	0.0710	0.52	3.74	2.78	3.25
1	10	0.2520	0.36	21.66	5.63	4.42
1	23.5	0.3848	0.3	52.95	8.5	5.29
1	34.5	0.4463	0.29	76.26	10.28	5.71
1	40	0.4694	0.25	89.21	11.04	5.77

µ_a
0.01	2	0.8047	0.44	4.65	0.29	4.97
0.01	10	0.9073	0.42	25.5	0.56	5.16
0.01	23.5	0.9386	0.42	61.4	0.85	5.21
0.01	34.5	0.9493	0.43	92.6	1.04	5.24
0.01	40	0.9527	0.43	105	1.12	5.25
0.1	2	0.5060	0.47	4.49	0.79	4.21
0.1	10	0.7354	0.45	24.7	1.74	4.76
0.1	23.5	0.8183	0.44	60.4	2.67	4.97
0.1	34.5	0.8478	0.43	89.4	3.23	4.99
0.1	40	0.8576	0.42	102	3.49	5.03
0.5	2	0.2291	0.55	3.77	1.88	3.57
0.5	10	0.5060	0.47	22.1	3.96	4.23
0.5	23.5	0.6398	0.44	57.1	5.98	4.51
0.5	34.5	0.6908	0.45	83.6	7.24	4.66
0.5	40	0.7095	0.46	100	7.79	4.69
1	2	0.1333	0.58	3.63	2.80	3.28
1	10	0.3853	0.50	22.7	5.64	3.87
1	23.5	0.5336	0.47	54.8	8.53	4.26
1	34.5	0.5940	0.47	85.7	10.25	4.40
1	40	0.6160	0.47	96.4	11.04	4.48

Analytic function for predicting light fluence rate of circular fields on a semi-infinite turbid medium

Abstract

Corrections

1. Introduction

2. Materials and Methods

2.1 MC Simulation

2.2 Analytical expression and fitting algorithm

3. Results

3.1 Fitting results of MC simulation

3.2 Relationship between fitting parameters and the tissue optical properties for air-tissue interface

3.2.1 Broad beam

3.2.2 Small circular field

4. Discussions

4.1 Effect of refraction index mismatch

5. Conclusion

6 Appendix: Fitting parameters

Funding

References and links

Cited By

Figures (11)

Tables (7)

Equations (24)

Optics Express

μ_a (cm⁻¹)	μ_s’ (cm⁻¹)
μ_a (cm⁻¹)	μ_s’ (cm⁻¹)	0.5	0.75	1	3	8
0.01	2	3.5	8.0	8.2	4.5	1.9 (0.9)
0.01	10	4.4	3.7	6.1	2.0	1.5(0.8)
0.01	23.5	2.7	1.6	2.1	1.0	1.3(1.2)
0.01	34.5	4.7	2.8	2.0	1.1	1.5(1.4)
0.01	40	5.2	3.9	2.2	1.5	1.7(1.9)
0.1	2	5.8	8.6	5.3	2.7	2.7(2.8)
0.1	10	3.1	4.5	3.3	1.7	1.0(0.8)
0.1	23.5	3.0	2.0	1.9	1.9	1.3(1.1)
0.1	34.5	1.2	1.3	1.3	1.2	0.8(0.6)
0.1	40	1.9	2.2	1.4	1.3	0.9(0.6)
0.5	2	6.8	6.8	3.5	5.4	4.8(1.1)
0.5	10	3.2	0.7	1.7	3.0	2.3(2.3)
0.5	23.5	1.1	0.6	0.9	0.4	0.6(0.7)
0.5	34.5	1.2	0.9	1.2	0.9	1.0(0.9)
0.5	40	1.2	1.0	1.2	1.1	1.1(0.9)
1	2	9.5	6.7	5.4	3.5	5.8(3.2)
1	10	5.9	4.4	3.7	3.9	3.0(2.0)
1	23.5	2.9	1.3	1.0	2.6	1.7(1.7)
1	34.5	2.2	1.0	1.0	2.9	1.5(1.7)
1	40	1.5	0.8	0.6	2.4	0.5(0.7)